Index: soft/giet_vm/applications/rosenfeld/src-par/mca.c
===================================================================
--- soft/giet_vm/applications/rosenfeld/src-par/mca.c	(revision 805)
+++ soft/giet_vm/applications/rosenfeld/src-par/mca.c	(revision 821)
@@ -15,12 +15,15 @@
 #include <malloc.h>
 
-#ifdef CLI
 #include "nrc_os_config.h"
+#include "config.h"
 #include "nrc.h"
-#endif
-
+
+#if TARGET_OS == GIETVM
+    #include <giet_config.h>
+#endif
 
 #include "util.h"
 #include "ecc_common.h"
+#include "ecc_features.h"
 #include "mca_matrix_dist.h"
 
@@ -34,18 +37,11 @@
 
 
-// ----------------------
+// -----------------------
 void MCA_Error(char * msg)
-// ----------------------
+// -----------------------
 {
     printf("MCA ERROR: %s\n", msg);
     printf("now exiting to system...\n");
-    exit(1);   
-}
-
-
-// ---------------------
-void MCA_Zero(MCA * mca)
-// ---------------------
-{
+    exit(1);
 }
 
@@ -58,8 +54,5 @@
     mca = (MCA *) malloc(sizeof(MCA));
 
-    if (mca) {
-        MCA_Zero(mca);
-    }
-    else {
+    if (mca == NULL) {
         MCA_Error("allocation failed in MCA_pConstructor_Empty");
     }
@@ -91,5 +84,5 @@
 
     mca->j0 = 0;
-    mca->j1 = width-1;
+    mca->j1 = width - 1;
 }
 
@@ -102,5 +95,5 @@
   
     mca->i0 = 0;
-    mca->i1 = height-1;
+    mca->i1 = height - 1;
 }
 
@@ -161,19 +154,17 @@
     int j0_par, j1_par;
     int height_par, height_mod;
-    
+
     int pw2;
     int32 ne_par;          // quantite par bande
     uint32 nemax_par;      // la puissance de 2 >=
     uint32 e0_par, e1_par; // indice par bande [start..end]
+    int nb_level;
     
     MCA ** mcas;
     MCA *  mca_par;
     
-    MCA_VERBOSE1(printf("====================\n"));
-    MCA_VERBOSE1(printf("== MCA_Initialize ==\n"));
-    MCA_VERBOSE1(printf("====================\n"));
-    
-    MCA_VERBOSE1(printf("height = %d\n", height));
-    MCA_VERBOSE1(printf("width  = %d\n", width));
+    printf("*** %s ***\n", __func__);
+    MCA_VERBOSE1(printf("   height = %d\n", height));
+    MCA_VERBOSE1(printf("   width  = %d\n", width));
     
     // array of pointers to mca workers
@@ -187,16 +178,57 @@
     height_par = height / np;
     height_mod = height % np;
-    //if(height % np) height_par++;
-    MCA_VERBOSE1(printf("height_par = %d x %d + %d\n", height_par, np, height_mod));
-    
-    MCA_VERBOSE1(printf("========================\n"));
+
+    MCA_VERBOSE1(printf("   height_par = %d x %d + %d\n", height_par, np, height_mod));
+    MCA_VERBOSE1(printf("   ========================\n"));
     
     i1_par_previous = 0;
-    
-    for (int p = 0; p < np; p++) {
+
+    // puissance de 2 de chaque bande
+    ne_par = height_par * width + 1;
+    MCA_VERBOSE1(printf("   ne_par    = %d\n", ne_par));
+    pw2 = i32log2(ne_par);
+    if (ne_par > (1 << pw2)) {
+        pw2++;
+    }
+    nemax_par = 1 << pw2;
+
+    MCA_VERBOSE1(printf("   nemax_par = %d\n", nemax_par));
+
+    nb_level = i32log2(np);
+    if ((1 << nb_level) < np) {
+        nb_level++;
+    }
+
+#if PYR_BARRIERS
+    // ------------------------------------------
+    // -- Allocation des barriÃšres pyramidales --
+    // ------------------------------------------
+
+    pthread_barrier_t * barriers = NULL;
+    if (nb_level > 0) {
+        barriers = malloc(sizeof(pthread_barrier_t) * nb_level);
+
+        // Initially all threads are active except thread 0
+        int nb_active = np - 1;
+        pthread_barrier_init(&barriers[0], NULL, nb_active);
+        for (int i = 1; i < nb_level; i++) {
+            // thread 0 never does any merge
+            for (int p = 1; p < np; p++) {
+                if ((p + (1 << (i - 1))) % (1 << i) == 0) {
+                    // thread inactive at level i
+                    nb_active -= 1;
+                }
+            }
+            pthread_barrier_init(&barriers[i], NULL, nb_active);
+        }
+    }
+#endif
+
+    for (int p = 0; p < np; p++) {
+
         // ----------------- //
         // -- constructor -- //
         // ----------------- //
-        MCA_VERBOSE1(printf("-- p = %d ----------------\n", p));
+        MCA_VERBOSE2(printf("-- p = %d ----------------\n", p));
     
         // alloc of mca workers into array of pointers
@@ -208,4 +240,11 @@
         mca_par->p   = p;
         mca_par->mca = mca; // pointer to master
+#if TARGET_OS == GIETVM
+        int x, y; // cluster coordinates
+        // We have p == 4 => x = 0; y = 1
+        x = (p / NB_PROCS_MAX) / Y_SIZE;
+        y = (p / NB_PROCS_MAX) % Y_SIZE;
+        MCA_VERBOSE2(printf("p = %d (x = %d, y = %d)\n", p, x, y));
+#endif
         
         // ------------------------------------- //
@@ -214,7 +253,4 @@
         
         // hauteur de chaque bande
-        
-        //printf("i1_par_previous = %d\n", i1_par_previous);
-        
         if (p == 0) {
             i0_par = 0;
@@ -232,19 +268,6 @@
         i1_par_previous = i1_par;
         
-        MCA_VERBOSE1(printf("i0_par = %d\n", i0_par));
-        MCA_VERBOSE1(printf("i1_par = %d\n", i1_par));
-        
-        // puissance de 2 de chaque bande
-        ne_par = height_par * width + 1;
-        //if (p == 0) {
-        //    ne_par++;
-        //}
-        MCA_VERBOSE1(printf("ne_par    = %d\n", ne_par));
-        pw2 = i32log2(ne_par);
-        if (ne_par > (1 << pw2)) {
-            pw2++;
-        }
-        nemax_par = 1 << pw2;
-        MCA_VERBOSE1(printf("nemax_par = %d\n", nemax_par));
+        MCA_VERBOSE2(printf("i0_par = %d\n", i0_par));
+        MCA_VERBOSE2(printf("i1_par = %d\n", i1_par));
         
         // etiquettes
@@ -260,22 +283,50 @@
         MCA_VERBOSE2(printf("e0_par = %d\n", e0_par));
         MCA_VERBOSE2(printf("e1_par = %d\n", e1_par));
-        
+
         mca_par->width  = width;
         mca_par->height = height_par;
-    
         mca_par->i0 = i0_par;
         mca_par->i1 = i1_par;
         mca_par->j0 = 0;
         mca_par->j1 = width - 1;
-    
         mca_par->e0 = e0_par;
         mca_par->e1 = e1_par;
-        
         mca_par->alpha = pw2;
+        mca_par->np = np;
+        // Pour les barriÃšres pyramidales
+        mca_par->nb_level = nb_level;
+#if PYR_BARRIERS
+        mca_par->barriers = barriers;
+#else
+        mca_par->barriers = NULL;
+#endif
+        mca_par->F = NULL; // default init
+        mca_par->stats = NULL; // default init
  
         // ---------------- //
         // -- allocation -- //
         // ---------------- //
-    
+#if TARGET_OS == GIETVM
+        mca_par->X = remote_ui8matrix(i0_par, i1_par, 0, width - 1, x, y);
+        mca_par->E = remote_dist_ui32matrix(i0_par, i1_par, 0, width - 1, x, y); // distributed matrix with border
+        
+        if (p == 0) {
+            mca_par->T = remote_ui32vector(e0_par - 1, e1_par, x, y); // car e0 = 1, on a besoin que T[0] = 0 pour FindRoot
+#if FEATURES
+            mca_par->stats = remote_RegionStatsVector(e0_par - 1, e1_par, x, y);
+#endif
+        }
+        else {
+            mca_par->T = remote_ui32vector(e0_par, e1_par, x, y);
+#if FEATURES
+            mca_par->stats = remote_RegionStatsVector(e0_par, e1_par, x, y);
+#endif
+        }
+        
+        mca_par->D = (uint32 **) remote_vvector(0, np - 1, x, y);
+#if FEATURES
+        mca_par->F = (RegionStats **) remote_vvector(0, np - 1, x, y);
+#endif
+#else // !GIETVM
         mca_par->X = ui8matrix (i0_par, i1_par, 0, width - 1);
         mca_par->E = dist_ui32matrix(i0_par, i1_par, 0, width - 1); // distributed matrix with border
@@ -283,22 +334,29 @@
         if (p == 0) {
             mca_par->T = ui32vector(e0_par - 1, e1_par); // car e0 = 1, on a besoin que T[0] = 0 pour FindRoot
+#if FEATURES
+            mca_par->stats = RegionStatsVector(e0_par - 1, e1_par);
+
+#endif
         }
         else {
             mca_par->T = ui32vector(e0_par, e1_par);
+#if FEATURES
+            mca_par->stats = RegionStatsVector(e0_par, e1_par);
+#endif
         }
         
         mca_par->D = (uint32 **) vvector(0, np - 1);
-        
+#if FEATURES
+        mca_par->F = (RegionStats **) vvector(0, np - 1);
+#endif
+#endif   
         MCA_VERBOSE2(printf("X = %p\n", mca_par->X));
         MCA_VERBOSE2(printf("E = %p\n", mca_par->E));
         MCA_VERBOSE2(printf("T = %p\n", mca_par->T));
         MCA_VERBOSE2(printf("D = %p\n", mca_par->D));
-    
     } // p
     
-    // pour debug
-    MCA_VERBOSE2(printf("init des tables d'EQ a l'identite\n"));
-    for (int p = 0; p < np; p++) {
-    
+
+    for (int p = 0; p < np; p++) {
         MCA * mca_par = mcas[p];
         
@@ -319,5 +377,4 @@
     MCA_VERBOSE2(printf("display des tables d'EQ\n"));
     for (int p = 0; p < np; p++) {
-        
         MCA * mca_par = mcas[p];
         
@@ -326,10 +383,10 @@
         uint32 e1 = mca_par->e1;
         
-        MCA_VERBOSE1(printf("p = %d T[%d..%d]\n", p, e0, e1));
+        MCA_VERBOSE2(printf("p = %d T[%d..%d]\n", p, e0, e1));
         if (p == 0) {
-            MCA_VERBOSE1(display_ui32vector_number(T, e0 - 1, e0 + 10, "%5d", "T"));
-        }
-        else {
-            MCA_VERBOSE1(display_ui32vector_number(T, e0, e0 + 10, "%5d", "T"));
+            MCA_VERBOSE2(display_ui32vector_number(T, e0 - 1, e0 + 10, "%5d", "T"));
+        }
+        else {
+            MCA_VERBOSE2(display_ui32vector_number(T, e0, e0 + 10, "%5d", "T"));
         }
         MCA_VERBOSE2(printf("\n"));
@@ -343,17 +400,17 @@
     // table d'indirection distribuee D
     MCA_VERBOSE2(printf("nemax_par = %d\n", nemax_par));
-    
-    for (int p = 0; p < np; p++) {
-    
+    for (int p = 0; p < np; p++) {
         MCA * mca_p = mcas[p];
         uint32 ** D = mca_p->D;
+        RegionStats ** F  = mca_p->F;
         
         for (int k = 0; k < np; k++) {
-            //mcas[p]->D[k] = (void*) (&(mcas[k]->T[mcas[k]->e0]))  - mcas[k]->e0; //k * nemax_par;
-            
             MCA * mca_k = mcas[k];
             uint32 * T = mca_k->T;
-            
             D[k] = T + k * nemax_par; // il faut soustraire le "MSB"
+#if FEATURES
+            RegionStats * stat = mca_k->stats;
+            F[k] = stat + k * nemax_par; // il faut soustraire le "MSB"
+#endif
         } // k
     } // p
@@ -376,5 +433,5 @@
             MCA_VERBOSE2(display_ui32vector(D[k], 0, 9, "%5d", "D\n"));
         }
-        printf("\n");
+        MCA_VERBOSE2(printf("\n"));
     }
     
@@ -413,25 +470,24 @@
 // -----------------------------------
 {
-    int p, np = mca->np;
+    int np = mca->np;
     
     MCA ** mcas = mca->mcas;
     MCA *  mca_par;
-    
-    printf("============================\n");
-    printf("== MCA_Display_Parameters ==\n");
-    printf("============================\n");
-    
-    printf("height = %d\n", mca->height);
-    printf("width  = %d\n", mca->width);
-    printf("np     = %d\n", mca->np);
-    
-    for (p = 0; p < np; p++) {
+    (void) mca_par;
+    
+    printf("*** MCA_Display_Parameters ***\n");
+    
+    MCA_VERBOSE1(printf("   height = %d\n", mca->height));
+    MCA_VERBOSE1(printf("   width  = %d\n", mca->width));
+    MCA_VERBOSE1(printf("   np     = %d\n", mca->np));
+    
+    for (int p = 0; p < np; p++) {
         mca_par = mcas[p];
         
-        printf("Display MCA[%d]\n", p);
-        printf("p = %d\n", mca_par->p);
-        printf("i0 = %8d  i1 = %8d\n", mca_par->i0, mca_par->i1);
-        printf("j0 = %8d  j1 = %8d\n", mca_par->j0, mca_par->j1);
-        printf("e0 = %8d  e1 = %8d\n", mca_par->e0, mca_par->e1);
+        MCA_VERBOSE2(printf("Display MCA[%d]\n", p));
+        MCA_VERBOSE2(printf("p = %d\n", mca_par->p));
+        MCA_VERBOSE2(printf("i0 = %8d  i1 = %8d\n", mca_par->i0, mca_par->i1));
+        MCA_VERBOSE2(printf("j0 = %8d  j1 = %8d\n", mca_par->j0, mca_par->j1));
+        MCA_VERBOSE2(printf("e0 = %8d  e1 = %8d\n", mca_par->e0, mca_par->e1));
     }
 }
@@ -442,5 +498,5 @@
 // -------------------------
 {
-    int p, np = mca->np;
+    int np = mca->np;
     
     MCA ** mcas = mca->mcas;
@@ -451,10 +507,11 @@
     uint32 e0, e1;
     
-    printf("==================\n");
-    printf("== MCA_Finalize ==\n");
-    printf("==================\n");
-    
-    for (p = 0; p < np; p++) {
-    
+    printf("*** MCA_Finalize ***\n");
+    
+#if PYR_BARRIERS
+    free(mcas[0]->barriers);
+#endif
+
+    for (int p = 0; p < np; p++) {
         mca_par = mcas[p];
     
@@ -475,21 +532,26 @@
         if (p == 0) {
             free_ui32vector(mca_par->T, e0 - 1, e1); // car e0 = 1, on a besoin que T[0] = 0 pour FindRoot
+#if FEATURES
+            free_RegionStatsVector(mca_par->stats, e0 - 1, e1);
+#endif
         }
         else {
             free_ui32vector(mca_par->T, e0, e1);
+#if FEATURES
+            free_RegionStatsVector(mca_par->stats, e0, e1);
+#endif
         }
         
         free_vvector((void **) mca_par->D, 0, np - 1);
-        
-    }
-    printf("[MCA_Finalize]: fin boucle\n");
+#if FEATURES
+        free_vvector((void **) mca_par->F, 0, np - 1);
+#endif
+        free(mca_par);
+    }
     free(mcas);
-    printf("[MCA_Finalize]: mcas freed\n");
-    free(mca); // plante si free XET, ne plante pas si pas de free ...
-    printf("[MCA_Finalize]: mca freed\n");
-}
-
-
-// @QM check my modifs...
+    free(mca);
+}
+
+
 // -------------------------------
 void MCA_Scatter_ImageX(MCA * mca)
@@ -498,11 +560,9 @@
     // diffusion de l'image binaire source
     
-    int      np = mca->np;
+    int np = mca->np;
     uint8 ** X  = mca->mca->X;
     
     if (mca->p == 0) { 
-        MCA_VERBOSE1(printf("------------------------\n"));
-        MCA_VERBOSE1(printf("-- MCA_Scatter_ImageX --\n"));
-        MCA_VERBOSE1(printf("------------------------\n"));
+        printf("*** MCA_Scatter_ImageX ***\n");
     }
     
@@ -524,5 +584,4 @@
 
 
-// @QM check my modifs...
 // ------------------------------
 void MCA_Gather_ImageL(MCA * mca)
@@ -532,4 +591,8 @@
     int np = mca->np;
     uint32 ** E = mca->mca->E;
+
+    if (mca->p == 0) { 
+        printf("*** MCA_Gather_ImageL ***\n");
+    }
 
     int i0 = mca->i0;
