[Scipy-svn] r4066 - trunk/scipy/ndimage/src/segment
scipy-svn@scip...
scipy-svn@scip...
Tue Apr 1 16:04:24 CDT 2008
Author: tom.waite
Date: 2008-04-01 16:04:22 -0500 (Tue, 01 Apr 2008)
New Revision: 4066
Modified:
trunk/scipy/ndimage/src/segment/Segmenter_IMPL.c
Log:
Added support for 3D blob extraction and improve the current 2D blob.
Modified: trunk/scipy/ndimage/src/segment/Segmenter_IMPL.c
===================================================================
--- trunk/scipy/ndimage/src/segment/Segmenter_IMPL.c 2008-04-01 21:04:08 UTC (rev 4065)
+++ trunk/scipy/ndimage/src/segment/Segmenter_IMPL.c 2008-04-01 21:04:22 UTC (rev 4066)
@@ -200,8 +200,289 @@
}
-int NI_GetBlobs(int samples, int rows, int cols, unsigned short *edges, unsigned short *connectedEdges, int *groups){
+int NI_GetBlobs3D(int samples, int layers, int rows, int cols, unsigned short *edges,
+ unsigned short *connectedEdges, int *groups, int mask){
+ int i, j, k, l, m;
+ int lOffset, rOffset, Label;
+ int lOffsetP, lOffsetN;
+ int rOffsetP, rOffsetN;
+ int Classes[4096];
+ int dwImageSize, ptr;
+ int HighZ, HighY, HighX, LowZ, LowX, LowY;
+ bool NewLabel;
+ bool Change;
+ bool connected;
+ int T[27];
+ int *ccompImage;
+ int layerSize;
+ int TargetArea;
+ int count;
+ int status;
+
+ layerSize = rows * cols;
+ dwImageSize = layers * rows * cols;
+ ccompImage = calloc(dwImageSize, sizeof(int ));
+
+ Label = 1;
+ for(i = 1; i < layers-1; ++i){
+ lOffset = i * layerSize;
+ lOffsetP = lOffset+layerSize;
+ lOffsetN = lOffset-layerSize;
+ for(j = 1; j < rows-1; ++j){
+ rOffset = j * cols;
+ rOffsetP = rOffset+cols;
+ rOffsetN = rOffset-cols;
+ for(k = 1; k < cols-1; ++k){
+ if(edges[lOffset+rOffset+k]){
+ /*
+ check 3x3x3 connectivity
+ */
+
+ T[0] = edges[lOffset+rOffset+k];
+ T[1] = edges[lOffset+rOffset+k+1];
+ T[2] = edges[lOffset+rOffsetN+k+1];
+ T[3] = edges[lOffset+rOffsetN+k];
+ T[4] = edges[lOffset+rOffsetN+k-1];
+ T[5] = edges[lOffset+rOffset+k-1];
+ T[6] = edges[lOffset+rOffsetP+k-1];
+ T[7] = edges[lOffset+rOffsetP+k];
+ T[8] = edges[lOffset+rOffsetP+k+1];
+
+ T[9] = edges[lOffsetN+rOffset+k];
+ T[10] = edges[lOffsetN+rOffset+k+1];
+ T[11] = edges[lOffsetN+rOffsetN+k+1];
+ T[12] = edges[lOffsetN+rOffsetN+k];
+ T[13] = edges[lOffsetN+rOffsetN+k-1];
+ T[14] = edges[lOffsetN+rOffset+k-1];
+ T[15] = edges[lOffsetN+rOffsetP+k-1];
+ T[16] = edges[lOffsetN+rOffsetP+k];
+ T[17] = edges[lOffsetN+rOffsetP+k+1];
+
+ T[18] = edges[lOffsetP+rOffset+k];
+ T[19] = edges[lOffsetP+rOffset+k+1];
+ T[20] = edges[lOffsetP+rOffsetN+k+1];
+ T[21] = edges[lOffsetP+rOffsetN+k];
+ T[22] = edges[lOffsetP+rOffsetN+k-1];
+ T[23] = edges[lOffsetP+rOffset+k-1];
+ T[24] = edges[lOffsetP+rOffsetP+k-1];
+ T[25] = edges[lOffsetP+rOffsetP+k];
+ T[26] = edges[lOffsetP+rOffsetP+k+1];
+
+ connected = FALSE;
+ if(mask == 1){
+ count = 0;
+ for(l = 1; l < 27; ++l){
+ count += T[l];
+ }
+ if(count){
+ connected = TRUE;
+ }
+ }
+ else if(mask == 6){
+ count = (T[2] + T[4] + T[6] + T[8] + T[9] + T[18]);
+ if(count == 6){
+ connected = TRUE;
+ }
+ }
+ else if(mask == 14){
+ count = (T[2] + T[4] + T[6] + T[8] + T[9] + T[18] + T[11] +
+ T[13] + T[15] + T[17] + T[20] + T[22] + T[24] + T[26]);
+ if(count == 14){
+ connected = TRUE;
+ }
+ }
+ else if(mask == 26){
+ count = 0;
+ for(l = 1; l < 27; ++l){
+ count += T[l];
+ }
+ if(count == 26){
+ connected = TRUE;
+ }
+ }
+ if(connected){
+ ccompImage[lOffset+rOffset+k] = Label++;
+ }
+ }
+ }
+ }
+ }
+
+
+ while(1){
+ Change = FALSE;
+ /*
+ // TOP-DOWN Pass for labeling
+ */
+ for(i = 1; i < layers-1; ++i){
+ lOffset = i * layerSize;
+ lOffsetP = lOffset+layerSize;
+ lOffsetN = lOffset-layerSize;
+ for(j = 1; j < rows-1; ++j){
+ rOffset = j * cols;
+ rOffsetP = rOffset+cols;
+ rOffsetN = rOffset-cols;
+ for(k = 1; k < cols-1; ++k){
+ if(ccompImage[lOffset+rOffset+k] != 0){
+
+ T[0] = ccompImage[lOffset+rOffset+k];
+ T[1] = ccompImage[lOffset+rOffset+k+1];
+ T[2] = ccompImage[lOffset+rOffsetN+k+1];
+ T[3] = ccompImage[lOffset+rOffsetN+k];
+ T[4] = ccompImage[lOffset+rOffsetN+k-1];
+ T[5] = ccompImage[lOffset+rOffset+k-1];
+ T[6] = ccompImage[lOffset+rOffsetP+k-1];
+ T[7] = ccompImage[lOffset+rOffsetP+k];
+ T[8] = ccompImage[lOffset+rOffsetP+k+1];
+
+ T[9] = ccompImage[lOffsetN+rOffset+k];
+ T[10] = ccompImage[lOffsetN+rOffset+k+1];
+ T[11] = ccompImage[lOffsetN+rOffsetN+k+1];
+ T[12] = ccompImage[lOffsetN+rOffsetN+k];
+ T[13] = ccompImage[lOffsetN+rOffsetN+k-1];
+ T[14] = ccompImage[lOffsetN+rOffset+k-1];
+ T[15] = ccompImage[lOffsetN+rOffsetP+k-1];
+ T[16] = ccompImage[lOffsetN+rOffsetP+k];
+ T[17] = ccompImage[lOffsetN+rOffsetP+k+1];
+
+ T[18] = ccompImage[lOffsetP+rOffset+k];
+ T[19] = ccompImage[lOffsetP+rOffset+k+1];
+ T[20] = ccompImage[lOffsetP+rOffsetN+k+1];
+ T[21] = ccompImage[lOffsetP+rOffsetN+k];
+ T[22] = ccompImage[lOffsetP+rOffsetN+k-1];
+ T[23] = ccompImage[lOffsetP+rOffset+k-1];
+ T[24] = ccompImage[lOffsetP+rOffsetP+k-1];
+ T[25] = ccompImage[lOffsetP+rOffsetP+k];
+ T[26] = ccompImage[lOffsetP+rOffsetP+k+1];
+
+ m = T[0];
+ for(l = 1; l < 27; ++l){
+ if(T[l] != 0){
+ if(T[l] < m) m = T[l];
+ }
+ }
+ if(m != ccompImage[lOffset+rOffset+k]){
+ Change = TRUE;
+ ccompImage[lOffset+rOffset+k] = m;
+ }
+ }
+ }
+ }
+ }
+ /*
+ // BOTTOM-UP Pass for labeling
+ */
+ for(i = layers-1; i > 0; --i){
+ lOffset = i * layerSize;
+ lOffsetP = lOffset+layerSize;
+ lOffsetN = lOffset-layerSize;
+ for(j = rows-1; j > 0; --j){
+ rOffset = j * cols;
+ rOffsetP = rOffset+cols;
+ rOffsetN = rOffset-cols;
+ for(k = cols-1; k > 0; --k){
+ if(ccompImage[lOffset+rOffset+k] != 0){
+
+ T[0] = ccompImage[lOffset+rOffset+k];
+ T[1] = ccompImage[lOffset+rOffset+k+1];
+ T[2] = ccompImage[lOffset+rOffsetN+k+1];
+ T[3] = ccompImage[lOffset+rOffsetN+k];
+ T[4] = ccompImage[lOffset+rOffsetN+k-1];
+ T[5] = ccompImage[lOffset+rOffset+k-1];
+ T[6] = ccompImage[lOffset+rOffsetP+k-1];
+ T[7] = ccompImage[lOffset+rOffsetP+k];
+ T[8] = ccompImage[lOffset+rOffsetP+k+1];
+
+ T[9] = ccompImage[lOffsetN+rOffset+k];
+ T[10] = ccompImage[lOffsetN+rOffset+k+1];
+ T[11] = ccompImage[lOffsetN+rOffsetN+k+1];
+ T[12] = ccompImage[lOffsetN+rOffsetN+k];
+ T[13] = ccompImage[lOffsetN+rOffsetN+k-1];
+ T[14] = ccompImage[lOffsetN+rOffset+k-1];
+ T[15] = ccompImage[lOffsetN+rOffsetP+k-1];
+ T[16] = ccompImage[lOffsetN+rOffsetP+k];
+ T[17] = ccompImage[lOffsetN+rOffsetP+k+1];
+
+ T[18] = ccompImage[lOffsetP+rOffset+k];
+ T[19] = ccompImage[lOffsetP+rOffset+k+1];
+ T[20] = ccompImage[lOffsetP+rOffsetN+k+1];
+ T[21] = ccompImage[lOffsetP+rOffsetN+k];
+ T[22] = ccompImage[lOffsetP+rOffsetN+k-1];
+ T[23] = ccompImage[lOffsetP+rOffset+k-1];
+ T[24] = ccompImage[lOffsetP+rOffsetP+k-1];
+ T[25] = ccompImage[lOffsetP+rOffsetP+k];
+ T[26] = ccompImage[lOffsetP+rOffsetP+k+1];
+
+ m = T[0];
+ for(l = 1; l < 27; ++l){
+ if(T[l] != 0){
+ if(T[l] < m) m = T[l];
+ }
+ }
+ if(m != ccompImage[lOffset+rOffset+k]){
+ Change = TRUE;
+ ccompImage[lOffset+rOffset+k] = m;
+ }
+ }
+ }
+ }
+ }
+
+ if(!Change) break;
+
+ } /* end while loop */
+
+ Label = 1;
+ Classes[0] = 0;
+ ptr = 0;
+ for(i = 0; i < layers; ++i){
+ for(j = 0; j < rows; ++j){
+ for(k = 0; k < cols; ++k){
+ m = ccompImage[ptr];
+ ++ptr;
+ if(m > 0){
+ NewLabel = TRUE;
+ for(l = 1; l < Label; ++l){
+ if(Classes[l] == m) NewLabel = FALSE;
+ }
+ if(NewLabel){
+ Classes[Label++] = m;
+ if(Label > 4000){
+ return 0;
+ }
+ }
+ }
+ }
+ }
+ }
+
+ ptr = 0;
+ for(i = 0; i < layers; ++i){
+ for(j = 0; j < rows; ++j){
+ for(k = 0; k < cols; ++k){
+ m = ccompImage[ptr];
+ for(l = 1; l < Label; ++l){
+ if(Classes[l] == m){
+ connectedEdges[ptr] = l;
+ break;
+ }
+ }
+ ++ptr;
+ }
+ }
+ }
+
+ free(ccompImage);
+
+ status = 1;
+ return(status);
+
+}
+
+int NI_GetBlobs2D(int samples, int rows, int cols, unsigned short *edges, unsigned short *connectedEdges,
+ int *groups, int mask){
+
int i, j, k, l, m;
int offset;
int Label;
@@ -209,10 +490,12 @@
int Classes[4096];
bool NewLabel;
bool Change;
+ bool connected;
+ int count;
unsigned short T[12];
/*
- // connected components labeling. pixels touch within 3x3 mask for edge connectedness.
+ // connected components labeling. pixels with 1, 4 or 8 connectedness.
*/
Label = 1;
offset = 0;
@@ -220,7 +503,34 @@
for(j = 0; j < cols; ++j){
connectedEdges[offset+j] = 0;
if(edges[offset+j] == 1){
- connectedEdges[offset+j] = Label++;
+ connected = FALSE;
+ if(mask == 1){
+ count = 0;
+ for(l = 1; l < 9; ++l){
+ count += T[l];
+ }
+ if(count){
+ connected = TRUE;
+ }
+ }
+ else if(mask == 4){
+ count = (T[2] + T[4] + T[6] + T[8]);
+ if(count == 4){
+ connected = TRUE;
+ }
+ }
+ else if(mask == 8){
+ count = 0;
+ for(l = 1; l < 9; ++l){
+ count += T[l];
+ }
+ if(count == 8){
+ connected = TRUE;
+ }
+ }
+ if(connected){
+ connectedEdges[offset+j] = Label++;
+ }
}
}
offset += cols;
More information about the Scipy-svn
mailing list