OilfieldToolsNet/LBPM
4
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Using median within regions lbpm_uCT_pp

Browse Source
This commit is contained in:
James E McClure 2016-04-25 19:19:42 -04:00
parent ca97cd9133
commit 09b8c452db
1 changed files with 45 additions and 13 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

58
tests/lbpm_uCT_pp.cpp
View File

@ -697,32 +697,64 @@ int main(int argc, char **argv)
//..........................................
// Compute the means for each region
float mean_plus,mean_minus;
float count_plus,count_minus;
int count_plus,count_minus;
float mean_plus_global,mean_minus_global;
float count_plus_global,count_minus_global;
int count_plus_global,count_minus_global;
float *TmpMed;
TmpMed = new float[nsx*nsy*nsz];
// Compute median for regions of distance function
count_plus=count_minus=0;
mean_plus=mean_minus=0;
for (k=1;k<nsz-1;k++){
for (j=1;j<nsy-1;j++){
for (i=1;i<nsx-1;i++){
if (spDist(i,j,k) > 0.0){
mean_plus += spM(i,j,k);
count_plus += 1.0;
}
else{
mean_minus += spM(i,j,k);
count_minus += 1.0; }
if (spDist(i,j,k) > 0.0)
TmpMed[count_plus++]= spM(i,j,k);
}
}
}
for (int ii=0; ii<count_plus/2+1; ii++){
for (int jj=ii+1; jj<count_plus; jj++){
if (TmpMed[jj] < TmpMed[ii]){
float tmp = TmpMed[ii];
TmpMed[ii] = TmpMed[jj];
TmpMed[jj] = tmp;
}
}
}
mean_plus = TmpMed[count_plus/2];
for (k=1;k<nsz-1;k++){
for (j=1;j<nsy-1;j++){
for (i=1;i<nsx-1;i++){
if (spDist(i,j,k) < 0.0)
TmpMed[count_minus++]= spM(i,j,k);
}
}
}
for (int ii=0; ii<count_minus/2+1; ii++){
for (int jj=ii+1; jj<count_minus; jj++){
if (TmpMed[jj] < TmpMed[ii]){
float tmp = TmpMed[ii];
TmpMed[ii] = TmpMed[jj];
TmpMed[jj] = tmp;
}
}
}
mean_minus = TmpMed[count_minus/2];
MPI_Allreduce(&mean_plus,&mean_plus_global,1,MPI_FLOAT,MPI_SUM,Dm.Comm);
MPI_Allreduce(&count_plus,&count_plus_global,1,MPI_FLOAT,MPI_SUM,Dm.Comm);
MPI_Allreduce(&mean_minus,&mean_minus_global,1,MPI_FLOAT,MPI_SUM,Dm.Comm);
MPI_Allreduce(&count_minus,&count_minus_global,1,MPI_FLOAT,MPI_SUM,Dm.Comm);
mean_plus_global /= count_plus_global;
mean_minus_global /= count_minus_global;
mean_plus_global /= nprocs;
mean_minus_global /= nprocs;
if (rank==0) printf(" Region 1 mean (+): %f, Region 2 mean (-): %f \n",mean_plus_global, mean_minus_global);
delete[] TmpMed;
//..........................................
// intialize distance based on segmentation