source: soft/giet_vm/applications/rosenfeld/src-par/mca.c @ 823

Last change on this file since 823 was 823, checked in by meunier, 8 years ago
  • Improved scripts for simulations and graphes
  • Continued to clean up the lib nrc2 (from nrio2x.x to nrmem1.c)
  • Added a version (Fast - Parmerge - No stats)
File size: 15.8 KB
Line 
1/* ------------- */
2/* --- mca.c --- */
3/* ------------- */
4
5/*
6 * Copyright (c) 2016 Lionel Lacassagne, LIP6, UPMC, CNRS
7 * Init  : 2016/03/03
8 * modif : 2016/03/07
9 */
10
11#include <stdio.h>
12#include <stdlib.h>
13#include <string.h>
14#include <math.h>
15#include <malloc.h>
16
17#include "nrc_os_config.h"
18#include "config.h"
19#include "nrc.h"
20
21#if TARGET_OS == GIETVM
22    #include <giet_config.h>
23#endif
24
25#include "util.h"
26#include "ecc_common.h"
27#include "ecc_features.h"
28#include "mca_matrix_dist.h"
29
30#include "palette.h"
31#include "bmpNR.h"
32#include "str_ext.h"
33
34
35// -- local --
36#include "mca.h"
37
38
39// -----------------------
40void MCA_Error(char * msg)
41// -----------------------
42{
43    printf("MCA ERROR: %s\n", msg);
44    exit(1);
45}
46
47
48// -------------------------------
49MCA * MCA_pConstructor_Empty(void)
50// -------------------------------
51{
52    MCA * mca;
53    mca = (MCA *) malloc(sizeof(MCA));
54
55    if (mca == NULL) {
56        MCA_Error("allocation failed in MCA_pConstructor_Empty");
57    }
58    return mca;
59}
60
61
62// ---------------------------------------
63void MCA_Set_ImageX(MCA * mca, uint8 ** X)
64// ---------------------------------------
65{
66    mca->X = X;
67}
68
69
70// --------------------------------------
71void MCA_Set_ImageL(MCA * mca, uint32 ** E)
72// --------------------------------------
73{
74    mca->E = E;
75}
76
77
78// ------------------------------------
79void MCA_Set_Width(MCA * mca, int width)
80// ------------------------------------
81{
82    mca->width = width;
83
84    mca->j0 = 0;
85    mca->j1 = width - 1;
86}
87
88
89// --------------------------------------
90void MCA_Set_Height(MCA * mca, int height)
91// --------------------------------------
92{
93    mca->height = height;
94 
95    mca->i0 = 0;
96    mca->i1 = height - 1;
97}
98
99
100// ------------------------------------------------
101void MCA_Set_Size(MCA * mca, int width, int height)
102// ------------------------------------------------
103{
104    MCA_Set_Width(mca, width);
105    MCA_Set_Height(mca, height);
106}
107
108
109// -----------------------------------------------------
110void MCA_Set_Dimension(MCA * mca, int width, int height)
111// -----------------------------------------------------
112{
113    MCA_Set_Width(mca, width);
114    MCA_Set_Height(mca, height);
115}
116
117
118// -------------------------------
119void MCA_Set_NP(MCA * mca, int np)
120// -------------------------------
121{
122    mca->np = np;
123}
124
125
126// -------------------------------
127void MCA_Set_NR(MCA * mca, int nr)
128// -------------------------------
129{
130    mca->nr = nr;
131}
132
133
134
135// ------------------------------------------------------------------
136uint32 MCA_CalcMaxLabels(int connection, uint32 height, uint32 width)
137// ------------------------------------------------------------------
138{
139    uint32 nemax = 0;
140   
141    if (connection == 4) {
142        nemax =  ((height + 1) / 2) * ((width + 1) / 2) + (height / 2) * (width / 2); // 4C
143    }
144    if (connection == 8) {
145        nemax = ((height + 1) / 2) * ((width + 1) / 2); // 8C
146    }
147    return nemax;
148}
149
150
151// ---------------------------
152void MCA_Initialize(MCA * mca)
153// ---------------------------
154{
155    // input
156    int np     = mca->np;
157    int nr     = mca->nr;
158    int width  = mca->width;
159    int height = mca->height;
160   
161    // variables
162    int i0_par, i1_par, i1_par_previous;
163    int j0_par, j1_par;
164    int height_par, height_mod;
165
166    int pw2;
167    int32 ne_par;          // quantite par bande
168    uint32 nemax_par;      // la puissance de 2 >=
169    uint32 e0_par, e1_par; // indice par bande [start..end]
170    int nb_level;
171   
172    MCA ** mcas;
173    MCA *  mca_par;
174   
175    MCA_VERBOSE1(printf("*** %s ***\n", __func__));
176    MCA_VERBOSE2(printf("   height = %d\n", height));
177    MCA_VERBOSE2(printf("   width  = %d\n", width));
178   
179    // array of pointers to mca workers
180    mcas = (MCA **) malloc(np * sizeof(MCA *));
181    if (mcas == NULL) {
182        MCA_Error("MCA_Initialize1");
183    }
184    mca->mcas = mcas;
185
186    // hauteur de chaque bande
187    height_par = height / np;
188    height_mod = height % np;
189
190    MCA_VERBOSE2(printf("   height_par = %d x %d + %d\n", height_par, np, height_mod));
191    MCA_VERBOSE2(printf("   ========================\n"));
192   
193    i1_par_previous = 0;
194
195    // puissance de 2 de chaque bande
196    ne_par = height_par * width + 1;
197    MCA_VERBOSE2(printf("   ne_par    = %d\n", ne_par));
198    pw2 = i32log2(ne_par);
199    if (ne_par > (1 << pw2)) {
200        pw2++;
201    }
202    nemax_par = 1 << pw2;
203    mca->alpha = pw2;
204
205    MCA_VERBOSE2(printf("   nemax_par = %d\n", nemax_par));
206
207    nb_level = i32log2(np);
208    if ((1 << nb_level) < np) {
209        nb_level++;
210    }
211
212#if PYR_BARRIERS
213    // ------------------------------------------
214    // -- Allocation des barriÚres pyramidales --
215    // ------------------------------------------
216
217    pthread_barrier_t * barriers = NULL;
218    if (nb_level > 0) {
219        barriers = malloc(sizeof(pthread_barrier_t) * nb_level);
220
221        // Initially all threads are active except thread 0
222        int nb_active = np - 1;
223        pthread_barrier_init(&barriers[0], NULL, nb_active);
224        for (int i = 1; i < nb_level; i++) {
225            // thread 0 never does any merge
226            for (int p = 1; p < np; p++) {
227                if ((p + (1 << (i - 1))) % (1 << i) == 0) {
228                    // thread inactive at level i
229                    nb_active -= 1;
230                }
231            }
232            pthread_barrier_init(&barriers[i], NULL, nb_active);
233        }
234    }
235#endif
236
237    for (int p = 0; p < np; p++) {
238
239        // ----------------- //
240        // -- constructor -- //
241        // ----------------- //
242        MCA_VERBOSE3(printf("-- p = %d ----------------\n", p));
243   
244        // alloc of mca workers into array of pointers
245        mca_par = MCA_pConstructor_Empty();
246        if (mca_par == NULL) {
247            MCA_Error("MCA_Initialize2\n");
248        }
249        mcas[p]      = mca_par;
250        mca_par->p   = p;
251        mca_par->mca = mca; // pointer to master
252#if TARGET_OS == GIETVM
253        int x, y; // cluster coordinates
254        // We have p == 4 => x = 0; y = 1
255        x = (p / NB_PROCS_MAX) / Y_SIZE;
256        y = (p / NB_PROCS_MAX) % Y_SIZE;
257        MCA_VERBOSE3(printf("p = %d (x = %d, y = %d)\n", p, x, y));
258#endif
259       
260        // ------------------------------------- //
261        // -- calcul des parametres: passe #1 -- //
262        // ------------------------------------- //
263       
264        // hauteur de chaque bande
265        if (p == 0) {
266            i0_par = 0;
267        }
268        else {
269            i0_par = i1_par_previous + 1;
270        }
271        if (height_mod) {
272            i1_par = i0_par + height_par;
273            height_mod = height_mod - 1;
274        }
275        else {
276            i1_par = i0_par + height_par - 1;
277        }
278        i1_par_previous = i1_par;
279       
280        MCA_VERBOSE3(printf("i0_par = %d\n", i0_par));
281        MCA_VERBOSE3(printf("i1_par = %d\n", i1_par));
282       
283        // etiquettes
284        if (p == 0) {
285            e0_par = 1;
286            e1_par = nemax_par - 1;
287        }
288        else {
289            e0_par = p * nemax_par;
290            e1_par = e0_par + nemax_par - 1;
291        }
292   
293        MCA_VERBOSE3(printf("e0_par = %d\n", e0_par));
294        MCA_VERBOSE3(printf("e1_par = %d\n", e1_par));
295
296        mca_par->width  = width;
297        mca_par->height = height_par;
298        mca_par->i0 = i0_par;
299        mca_par->i1 = i1_par;
300        mca_par->j0 = 0;
301        mca_par->j1 = width - 1;
302        mca_par->e0 = e0_par;
303        mca_par->e1 = e1_par;
304        // à la premiÚre itération, on remet à 0 toute la table T
305        mca_par->ne_prev = e1_par;
306        mca_par->alpha = pw2;
307        mca_par->np = np;
308        mca_par->nr = nr;
309        // Pour les barriÚres pyramidales
310        mca_par->nb_level = nb_level;
311#if PYR_BARRIERS
312        mca_par->barriers = barriers;
313#else
314        mca_par->barriers = NULL;
315#endif
316        mca_par->F = NULL; // default init
317        mca_par->stats = NULL; // default init
318 
319        // ---------------- //
320        // -- allocation -- //
321        // ---------------- //
322#if TARGET_OS == GIETVM
323        mca_par->X = remote_ui8matrix(i0_par, i1_par, 0, width - 1, x, y);
324        mca_par->E = remote_dist_ui32matrix(i0_par, i1_par, 0, width - 1, x, y); // distributed matrix with border
325       
326        mca_par->T = remote_ui32vector(e0_par, e1_par, x, y);
327        mca_par->stats = remote_RegionStatsVector(e0_par, e1_par, x, y);
328       
329        mca_par->D = (uint32 **) remote_vvector(0, np - 1, x, y);
330        mca_par->F = (RegionStats **) remote_vvector(0, np - 1, x, y);
331#else // !GIETVM
332        mca_par->X = ui8matrix (i0_par, i1_par, 0, width - 1);
333        mca_par->E = dist_ui32matrix(i0_par, i1_par, 0, width - 1); // distributed matrix with border
334       
335        mca_par->T = ui32vector(e0_par, e1_par);
336        mca_par->stats = RegionStatsVector(e0_par, e1_par);
337       
338        mca_par->D = (uint32 **) vvector(0, np - 1);
339        mca_par->F = (RegionStats **) vvector(0, np - 1);
340#endif   
341        MCA_VERBOSE3(printf("X = %p\n", mca_par->X));
342        MCA_VERBOSE3(printf("E = %p\n", mca_par->E));
343        MCA_VERBOSE3(printf("T = %p\n", mca_par->T));
344        MCA_VERBOSE3(printf("D = %p\n", mca_par->D));
345    } // p
346   
347
348    for (int p = 0; p < np; p++) {
349        MCA * mca_par = mcas[p];
350       
351        uint32 * T = mca_par->T;
352        uint32 e0 = mca_par->e0;
353        uint32 e1 = mca_par->e1;
354   
355        MCA_VERBOSE3(printf("p = %d T[%d..%d]\n", p, e0, e1));
356        set_ui32vector_j(T, e0, e1);
357    }
358   
359    MCA_VERBOSE3(printf("display des tables d'EQ\n"));
360    for (int p = 0; p < np; p++) {
361        MCA * mca_par = mcas[p];
362       
363        uint32 * T = mca_par->T;
364        uint32 e0 = mca_par->e0;
365        uint32 e1 = mca_par->e1;
366       
367        MCA_VERBOSE3(printf("p = %d T[%d..%d]\n", p, e0, e1));
368        MCA_VERBOSE3(display_ui32vector_number(T, e0, e0 + 10, "%5d", "T"));
369        MCA_VERBOSE3(printf("\n"));
370    }
371    //exit(-1);
372   
373    // ------------------------------------------------------------- //
374    // -- calcul des parametres: passe #2 = parametres distribues -- //
375    // ------------------------------------------------------------- //
376   
377    // table d'indirection distribuee D
378    MCA_VERBOSE3(printf("nemax_par = %d\n", nemax_par));
379    for (int p = 0; p < np; p++) {
380        MCA * mca_p = mcas[p];
381        uint32 ** D = mca_p->D;
382        RegionStats ** F  = mca_p->F;
383       
384        for (int k = 0; k < np; k++) {
385            MCA * mca_k = mcas[k];
386            uint32 * T = mca_k->T;
387            D[k] = T + k * nemax_par; // il faut soustraire le "MSB"
388            RegionStats * stat = mca_k->stats;
389            F[k] = stat + k * nemax_par; // il faut soustraire le "MSB"
390        } // k
391    } // p
392   
393    MCA_VERBOSE3(printf("table d'indirection distribuee D\n"));
394   
395    for (int p = 0; p < np; p++) {
396        MCA_VERBOSE3(printf("== p = %d ==========\n", p));
397       
398        MCA * mca_p = mcas[p];
399        uint32 ** D = mca_p->D;
400       
401        for (int k = 0; k < np; k++) {
402            MCA * mca_k = mcas[k];
403            uint32 * T = mca_k->T;
404           
405            uint32 e0 = mca_k->e0;
406            uint32 e1 = mca_k->e1;
407            MCA_VERBOSE3(display_ui32vector_number(T, e0, e0 + 9, "%5d", "T"));
408            MCA_VERBOSE3(display_ui32vector(D[k], 0, 9, "%5d", "D\n"));
409        }
410        MCA_VERBOSE3(printf("\n"));
411    }
412   
413    for (int p = 0; p < np; p++) {
414        if (p > 0) {
415            //printf("i0_(%d) = %d i1_{%d} = %d\n", p, mcas[p]->i0, p-1, mcas[p-1]->i1);
416            mcas[p]->E[mcas[p]->i0 - 1] = mcas[p - 1]->E[mcas[p - 1]->i1];
417           
418            /*printf("E[%2d] = E[%2d] = %p\n", mcas[p    ]->i0 - 1,
419                                             mcas[p - 1]->i1,
420                                             mcas[p - 1]->E[mcas[p - 1]->i1]);*/
421        }
422        if (p < np - 1) {
423            //printf("i1_(%d) = %d i0_{%d} = %d\n", p, mcas[p]->i1, p+1, mcas[p-1]->i0);
424            mcas[p]->E[mcas[p]->i1 + 1] = mcas[p + 1]->E[mcas[p + 1]->i0];
425           
426            /*printf("E[%2d] = E[%2d] = %p\n", mcas[p    ]->i1 + 1,
427                                             mcas[p + 1]->i0,
428                                             mcas[p + 1]->E[mcas[p + 1]->i1]);*/
429        }
430    }
431}
432
433
434// -----------------------------------
435void MCA_Display_Parameters(MCA * mca)
436// -----------------------------------
437{
438    int np = mca->np;
439   
440    MCA ** mcas = mca->mcas;
441    MCA *  mca_par;
442    (void) mca_par;
443   
444    MCA_VERBOSE1(printf("*** MCA_Display_Parameters ***\n"));
445   
446    MCA_VERBOSE2(printf("   height = %d\n", mca->height));
447    MCA_VERBOSE2(printf("   width  = %d\n", mca->width));
448    MCA_VERBOSE2(printf("   np     = %d\n", mca->np));
449   
450    for (int p = 0; p < np; p++) {
451        mca_par = mcas[p];
452       
453        MCA_VERBOSE3(printf("Display MCA[%d]\n", p));
454        MCA_VERBOSE3(printf("p = %d\n", mca_par->p));
455        MCA_VERBOSE3(printf("i0 = %8d  i1 = %8d\n", mca_par->i0, mca_par->i1));
456        MCA_VERBOSE3(printf("j0 = %8d  j1 = %8d\n", mca_par->j0, mca_par->j1));
457        MCA_VERBOSE3(printf("e0 = %8d  e1 = %8d\n", mca_par->e0, mca_par->e1));
458    }
459}
460
461
462// -------------------------
463void MCA_Finalize(MCA * mca)
464// -------------------------
465{
466    int np = mca->np;
467   
468    MCA ** mcas = mca->mcas;
469    MCA *  mca_par;
470   
471    int i0, i1;
472    int j0, j1;
473    uint32 e0, e1;
474   
475    MCA_VERBOSE1(printf("*** MCA_Finalize ***\n"));
476   
477#if PYR_BARRIERS
478    free(mcas[0]->barriers);
479#endif
480
481    for (int p = 0; p < np; p++) {
482        mca_par = mcas[p];
483   
484        i0 = mca_par->i0;
485        i1 = mca_par->i1;
486        j0 = mca_par->j0;
487        j1 = mca_par->j1;
488        e0 = mca_par->e0;
489        e1 = mca_par->e1;
490   
491        // ---------- //
492        // -- free -- //
493        // ---------- //
494       
495        free_ui8matrix (mca_par->X, i0, i1, j0, j1);
496        free_dist_ui32matrix(mca_par->E, i0, i1, j0, j1);
497       
498        free_ui32vector(mca_par->T, e0, e1);
499        free_RegionStatsVector(mca_par->stats, e0, e1);
500       
501        free_vvector((void **) mca_par->D, 0, np - 1);
502        free_vvector((void **) mca_par->F, 0, np - 1);
503        free(mca_par);
504    }
505    free(mcas);
506    free(mca);
507}
508
509
510// -------------------------------
511void MCA_Scatter_ImageX(MCA * mca)
512// -------------------------------
513{
514    // diffusion de l'image binaire source
515   
516    int np = mca->np;
517    uint8 ** X  = mca->mca->X;
518   
519    if (mca->p == 0) { 
520        MCA_VERBOSE1(printf("*** MCA_Scatter_ImageX ***\n"));
521    }
522   
523    int i0    = mca->i0;
524    int i1    = mca->i1;
525    int width = mca->width;
526
527    uint8 **  X_par = mca->X;
528    uint32 ** E_par = mca->E;
529
530    //printf("copie de [%d..%d]x[%d..%d]\n", i0, i1, 0, width - 1);
531    for (int i = i0; i <= i1; i++) {
532        for (int j = 0; j <= width - 1; j++) {
533            X_par[i][j] = X[i][j];
534            E_par[i][j] = 0; // inutile normalement car ecriture de 0
535        }
536    }
537}
538
539
540// ------------------------------
541void MCA_Gather_ImageL(MCA * mca)
542// ------------------------------
543{
544    // recuperation de l'image d'etiquettes
545    int np = mca->np;
546    uint32 ** E = mca->mca->E;
547
548    if (mca->p == 0) { 
549        MCA_VERBOSE1(printf("*** MCA_Gather_ImageL ***\n"));
550    }
551
552    int i0 = mca->i0;
553    int i1 = mca->i1;
554    int width = mca->width;
555
556    uint32 ** E_par = mca->E;
557
558    for (int i = i0; i <= i1; i++) {
559        for (int j = 0; j <= width - 1; j++) {
560            E[i][j] = E_par[i][j];
561        }
562    }
563}
564
565
566// Local Variables:
567// tab-width: 4
568// c-basic-offset: 4
569// c-file-offsets:((innamespace . 0)(inline-open . 0))
570// indent-tabs-mode: nil
571// End:
572
573// vim: filetype=cpp:expandtab:shiftwidth=4:tabstop=4:softtabstop=4
574
575
Note: See TracBrowser for help on using the repository browser.