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

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