calcul_stats_image<<<1, 1>>>( *d_img_x, *d_img_x2, H, L, (uint64*)*d_stats_snake);
- cudaThreadSynchronize() ;
- toc(chrono, "\tTemps GPU");
- if(DEBUG_IMG_CUMUL)
+ cudaThreadSynchronize() ;
+ toc(chrono, "\tTemps GPU");
+ if(DEBUG_IMG_CUMUL)
{
//allocation memoire CPU
* generation snake en mem GPU
*/
int dist = 140 ;
-
- tic(&chrono, NULL);
- int MAX_DIAGOS = 1024*65536 ;
- int ret, Q = 120 ;
- int * d_n_diagos, h_n_diagos;
- uint4 * d_diagos_snake ;
- uint4 * h_diagos_snake = new uint4[MAX_DIAGOS];
-
-
- cudaMalloc( (void**) &d_n_diagos, sizeof(int)) ;
- cudaMalloc( (void**) &d_diagos_snake, MAX_DIAGOS*sizeof(uint4)) ;
-
- genere_diagos_rectangle<<<1,1>>>(d_diagos_snake, H,L,Q, d_n_diagos);
+ /* Test de determination du snake rectangle initial optimal*/
+ int div = 100;//nb de divisions de l'image : cela définit le pas. La valeur max découle du nb max de threads possible ds une grille
+ int Nperm = div*div*bs;//nb total de rectangles a tester. La distribution est ainsi irrégulière, mais plus simple.
+ double best_crit ;
+ int ind_best_crit ;
+
+ t_rectangle_snake * d_all_crit, d_best_crit;//tableaux pour les résultats des différents rectangles / le meilleur
+ t_rectangle_snake * h_all_crit = new t_rectangle_snake[Nperm];//correspondant CPU
- cudaMemcpy( &h_n_diagos, d_n_diagos, sizeof(int), cudaMemcpyDeviceToHost) ;
- ret = cudaMemcpy( h_diagos_snake, d_diagos_snake, MAX_DIAGOS*sizeof(uint4), cudaMemcpyDeviceToHost) ;
+ //allocations
+ cudaMalloc((void**) &d_all_crit, Nperm*sizeof(t_rectangle_snake));
+ cudaMalloc((void**) &d_best_crit, sizeof(t_rectangle_snake));
+
+ tic(&chrono, NULL);
- toc(chrono, "\tCalcul diagos");
+ //execution kernel
+ dim3 grid = dim3(H/div, L/div, 1);
+ calcul_contribs_snake4<<<grid, bs, CFI(bs)*sizeof(tcontribs) >>>(*d_snake, *d_img_x, *d_img_x2, H, L, *d_stats_snake, d_all_crit) ;
+ cudaThreadSynchronize();
+ toc(chrono, "\nCALCULS RECTANGLES");
- printf("COPY : %s, DIAGOS = %d / %d\n", (ret==0)?"OK":"ERREUR", h_n_diagos,MAX_DIAGOS);
- for (int d=0; d<200;d++){
- printf("| (%d,%d)-(%d,%d) ", h_diagos_snake[d].x, h_diagos_snake[d].y, h_diagos_snake[d].z, h_diagos_snake[d].w );
+ //recup data rectangles
+ int ret;
+ ret = cudaMemcpy( h_all_crit, d_all_crit, Nperm*sizeof(t_rectangle_snake), cudaMemcpyDeviceToHost) ;
+ printf("COPIE DATA = %s\n",(ret==0)?"OK":"ERR");
+
+ //optimum sur CPU
+ best_crit = h_all_crit[0].crit ;
+ ind_best_crit = 0 ;
+ for (int k=1; k<100; k++){
+ if ((h_all_crit[k].crit > 0) && (h_all_crit[k].crit < best_crit)) {
+ best_crit = h_all_crit[k].crit ;
+ ind_best_crit = k ;
+ }
+ printf("%d -> ( %d, %d )--( %d, %d) CRITERE = %f\n", k, h_all_crit[k].bpi, h_all_crit[k].bpj,
+ h_all_crit[k].opi, h_all_crit[k].opj, h_all_crit[k].crit );
}
+ printf("BEST RECTANGLE/%d tests : %d -> ( %d, %d )--( %d, %d) CRITERE = %f\n", Nperm, ind_best_crit, h_all_crit[ind_best_crit].bpi, h_all_crit[ind_best_crit].bpj,
+ h_all_crit[ind_best_crit].opi, h_all_crit[ind_best_crit].opj, best_crit );
+
exit(0);
- genere_snake_rectangle_4nodes_gpu<<< 1, 1>>>(*d_snake, 140, H, L) ;
+ /*fin test snake rectangle initial optimal*/
+
+ //genere_snake_rectangle_4nodes_gpu<<< 1, 1>>>(*d_snake, 140, H, L) ;
int nnodes = nb_nodes ;
d_stats_snake[3], d_stats_snake[4], d_stats_snake[5]);
}
+// kernel d'init rectangulaire au niveau snake
+// un block par point de base et un point opposé par thread du bloc
-__global__ void calcul_contribs_snake4(t_cumul_x * cumul_x, t_cumul_x2 * cumul_x2, int h, int l, tcontribs * gcontribs,
- uint64 SUM_1, uint64 SUM_X, uint64 SUM_X2)
+__global__ void calcul_contribs_snake4(snake_node_gpu * snake, t_cumul_x * cumul_x, t_cumul_x2 * cumul_x2, int h, int l,
+ int64 * sums, t_rectangle_snake * critere)
{
// nb de diagonales testees par bloc (ie. par point de base NO)
- int blockSize = blockDim.x ;
+ int blockSize = blockDim.x ;
// indice du second point de chaque diagonale (=Opposite Point, = point SE)
int OPib = threadIdx.x ;
// coordonnees de chaque point de base (NO)
int BPi = blockIdx.x ;
int BPj = blockIdx.y ;
//coordonnees de chaque Opposite Point (SE)
- int OPi = OPib / (l - BPj) ;
- int OPj = OPib - (l - BPj)*OPi ;
+ double incThread = ((h-BPi)*(l-BPj) + blockSize -1)/(double)blockSize ;
+ int OPi = (double)(OPib*incThread)/(l - BPj) ;
+ int OPj = OPib*incThread - OPi*(l-BPj) ;
OPi += BPi ;
OPj += BPj ;
+ if (OPi >= h) OPi = h-1 ;
+ if (OPj >= l) OPj = l-1 ;
//indices des pixels dans les images cumulees
int posG, posD;
//contrib 1 du snake
int C1 = (OPi - BPi)*(OPj - BPj) ;
-
-
+
//pour stocker contribs de chaque snake d'un block
- //TODO on peut utiliser une structure restreinte (sans le c1) = gain d'espace
+ //TODO on peut utiliser une structure restreinte (sans le c1)
extern __shared__ tcontribs scumuls[];
-
+
//calcul contribs du snake
+ scumuls[CFI(OPib)].cx = 0;
+ scumuls[CFI(OPib)].cx2 = 0;
for (int k=BPi ; k < OPi ; k++)
{
posG = (BPi+k)*l + BPj ;
posD = posG - BPj + OPj ;
- scumuls[CFI(OPib)].cx += cumul_x[ posD ] - cumul_x[ posG ] ;
- scumuls[CFI(OPib)].cx2 += cumul_x2[ posD ] - cumul_x2[ posG ];
+ scumuls[CFI(OPib)].cx += (cumul_x[ posD ] - cumul_x[ posG ]) ;
+ scumuls[CFI(OPib)].cx2 += (cumul_x2[ posD ] - cumul_x2[ posG ]);
}
//calcul de critère pour chaque snake
uint64 stat_sum_xe ; /* somme des xn region exterieure */
uint32 ne ; /* nombre de pixel region exterieure */
double sigi2, sige2; /* variance region interieure et exterieure */
- double criterion;
+ int index_crit;
/* variance des valeurs des niveaux de gris a l'interieur du snake */
sigi2 =
((double)scumuls[CFI(OPib)].cx/(uint64)C1)*((double)scumuls[CFI(OPib)].cx/(uint64)C1) ;
/* variance des valeurs des niveaux de gris a l'exterieur du snake */
- ne = SUM_1 - C1 ;
- stat_sum_xe = SUM_X - scumuls[CFI(OPib)].cx ;
+
+ ne = sums[3] - C1 ;
+ stat_sum_xe = sums[4] - scumuls[CFI(OPib)].cx ;
sige2 =
- ((double)SUM_X2-scumuls[CFI(OPib)].cx2)/(double)ne -
+ ((double)sums[5]-scumuls[CFI(OPib)].cx2)/(double)ne -
((double)stat_sum_xe/(uint64)ne)*((double)stat_sum_xe/(uint64)ne) ;
- if ((sigi2 > 0)|(sige2 > 0))
- criterion = 0.5*((double)C1*log(sigi2) + (double)ne*log(sige2)) ;
+ index_crit = blockSize*(BPi*gridDim.y + BPj) + OPib ;
+
- //tri meilleur snake du bloc ( necessite de passer SUM_1, SUM_X et SUM_X2 )
+ if ((sigi2 > 0)|(sige2 > 0))
+ {
+ critere[ index_crit ].bpi = BPi;
+ critere[ index_crit ].bpj = BPj;
+ critere[ index_crit ].opi = OPi;
+ critere[ index_crit ].opj = OPj;
+ critere[ index_crit ].crit = 0.5*((double)C1*log(sigi2) + (double)ne*log(sige2)) ;
+ }
+ else
+ {
+ critere[ index_crit ].bpi = BPi;
+ critere[ index_crit ].bpj = BPj;
+ critere[ index_crit ].opi = OPi;
+ critere[ index_crit ].opj = OPj;
+ critere[ index_crit ].crit = -1 ;
+ }
+ // identification meilleur snake du bloc
+ // laissé au CPU pour test mais le principe de ce kernel n'est pas efficace.
+
}