wassname/scikit-image
1
0
Fork 0
mirror of https://github.com/wassname/scikit-image.git synced 2026-08-13 12:40:24 +08:00
Code Issues Packages Projects Releases Wiki Activity

unwrap: Fix whitespace in c code.

Browse Source
This commit is contained in:
Jostein Bø Fløystad
2013-11-22 10:42:28 +01:00
parent e72e5de06e
commit f67b03643a
2 changed files with 167 additions and 167 deletions
Download Patch File Download Diff File Expand all files Collapse all files
skimage/exposure/Hussein_3D_unwrapper_with_mask_and_wrap_around_option.c
+104 -104
View File
@@ -1,17 +1,17 @@
// 3D phase unwrapping, modified for inclusion in scipy by Gregor Thalhammer
//This program was written by Hussein Abdul-Rahman and Munther Gdeisat to program the three-dimensional phase unwrapper
//entitled "Fast three-dimensional phase-unwrapping algorithm based on sorting by
//entitled "Fast three-dimensional phase-unwrapping algorithm based on sorting by
//reliability following a noncontinuous path"
//by Hussein Abdul-Rahman, Munther A. Gdeisat, David R. Burton, and Michael J. Lalor,
//published in the Proceedings of SPIE -
//by Hussein Abdul-Rahman, Munther A. Gdeisat, David R. Burton, and Michael J. Lalor,
//published in the Proceedings of SPIE -
//The International Society for Optical Engineering, Vol. 5856, No. 1, 2005, pp. 32-40
//This program was written by Munther Gdeisat, Liverpool John Moores University, United Kingdom.
//Date 31st August 2007
//The wrapped phase volume is assumed to be of floating point data type. The resultant unwrapped phase volume is also of floating point type.
//Read the data from the file frame by frame
//The mask is of byte data type.
//When the mask is 255 this means that the voxel is valid
//The mask is of byte data type.
//When the mask is 255 this means that the voxel is valid
//When the mask is 0 this means that the voxel is invalid (noisy or corrupted voxel)
//This program takes into consideration the image wrap around problem encountered in MRI imaging.
@@ -32,7 +32,7 @@ typedef struct
int y_connectivity;
int z_connectivity;
int no_of_edges;
} params_t;
} params_t;
//VOXELM information
struct VOXELM
@@ -55,7 +55,7 @@ typedef struct VOXELM VOXELM;
//the EDGE is the line that connects two voxels.
//if we have S voxels, then we have S horizontal edges and S vertical edges
struct EDGE
{
{
float reliab; //reliabilty of the edge and it depends on the two voxels
VOXELM *pointer_1; //pointer to the first voxel
VOXELM *pointer_2; //pointer to the second voxel
@@ -157,14 +157,14 @@ void initialiseVOXELs(float *WrappedVolume, unsigned char *input_mask, unsigned
for (j=0; j < volume_width; j++)
{
voxel_pointer->increment = 0;
voxel_pointer->number_of_voxels_in_group = 1;
voxel_pointer->number_of_voxels_in_group = 1;
voxel_pointer->value = *wrapped_volume_pointer;
voxel_pointer->reliability = 9999999.f + rand();
voxel_pointer->input_mask = *input_mask_pointer;
voxel_pointer->extended_mask = *extended_mask_pointer;
voxel_pointer->head = voxel_pointer;
voxel_pointer->last = voxel_pointer;
voxel_pointer->next = NULL;
voxel_pointer->next = NULL;
voxel_pointer->new_group = 0;
voxel_pointer->group = -1;
voxel_pointer++;
@@ -190,7 +190,7 @@ float wrap(float voxel_value)
// voxelL_value is the left voxel, voxelR_value is the right voxel
int find_wrap(float voxelL_value, float voxelR_value)
{
float difference;
float difference;
int wrap_value;
difference = voxelL_value - voxelR_value;
@@ -199,14 +199,14 @@ int find_wrap(float voxelL_value, float voxelR_value)
else wrap_value = 0;
return wrap_value;
}
}
void extend_mask(unsigned char *input_mask, unsigned char *extended_mask, int volume_width, int volume_height, int volume_depth, params_t *params)
{
int n, i, j;
int vw = volume_width, vh = volume_height, vd = volume_depth;
int fs = volume_width * volume_height; //frame size
int frame_size = volume_width * volume_height;
int frame_size = volume_width * volume_height;
int volume_size = volume_width * volume_height * volume_depth; //volume size
int vs = volume_size;
unsigned char *IMP = input_mask + frame_size + volume_width + 1; //input mask pointer
@@ -219,7 +219,7 @@ void extend_mask(unsigned char *input_mask, unsigned char *extended_mask, int vo
{
for (j=1; j < volume_width - 1; j++)
{
if( (*IMP) == NOMASK && (*(IMP - 1) == NOMASK) && (*(IMP + 1) == NOMASK) &&
if( (*IMP) == NOMASK && (*(IMP - 1) == NOMASK) && (*(IMP + 1) == NOMASK) &&
(*(IMP + vw) == NOMASK) && (*(IMP + vw - 1) == NOMASK) && (*(IMP + vw + 1) == NOMASK) &&
(*(IMP - vw) == NOMASK) && (*(IMP - vw - 1) == NOMASK) && (*(IMP - vw + 1) == NOMASK) &&
(*(IMP + fs) == NOMASK) && (*(IMP + fs - 1) == NOMASK) && (*(IMP + fs + 1) == NOMASK) &&
@@ -228,7 +228,7 @@ void extend_mask(unsigned char *input_mask, unsigned char *extended_mask, int vo
(*(IMP - fs) == NOMASK) && (*(IMP - fs - 1) == NOMASK) && (*(IMP - fs + 1) == NOMASK) &&
(*(IMP - fs - vw) == NOMASK) && (*(IMP - fs - vw - 1) == NOMASK) && (*(IMP - fs - vw + 1) == NOMASK) &&
(*(IMP - fs + vw) == NOMASK) && (*(IMP - fs + vw - 1) == NOMASK) && (*(IMP - fs + vw + 1) == NOMASK))
{
{
*EMP = NOMASK;
}
++EMP;
@@ -239,7 +239,7 @@ void extend_mask(unsigned char *input_mask, unsigned char *extended_mask, int vo
}
EMP += 2 * volume_width;
IMP += 2 * volume_width;
}
}
if (params->x_connectivity == 1)
{
@@ -251,19 +251,19 @@ void extend_mask(unsigned char *input_mask, unsigned char *extended_mask, int vo
for (i=1; i < volume_height - 1; i++)
{
if( (*IMP) == NOMASK && (*(IMP + vw - 1) == NOMASK) && (*(IMP + 1) == NOMASK) &&
(*(IMP - vw) == NOMASK) && (*(IMP + vw) == NOMASK) &&
(*(IMP - vw) == NOMASK) && (*(IMP + vw) == NOMASK) &&
(*(IMP - fs) == NOMASK) && (*(IMP + fs) == NOMASK) &&
(*(IMP - 1) == NOMASK) && (*(IMP + vw + 1) == NOMASK) &&
(*(IMP - vw + 1) == NOMASK) && (*(IMP + 2 * vw - 1) == NOMASK) &&
(*(IMP - fs - 1) == NOMASK) && (*(IMP + fs + vw + 1) == NOMASK) &&
(*(IMP - vw + 1) == NOMASK) && (*(IMP + 2 * vw - 1) == NOMASK) &&
(*(IMP - fs - 1) == NOMASK) && (*(IMP + fs + vw + 1) == NOMASK) &&
(*(IMP - fs - vw) == NOMASK) && (*(IMP + fs + vw) == NOMASK) &&
(*(IMP - fs - vw + 1) == NOMASK) && (*(IMP + fs + 2 * vw - 1) == NOMASK) &&
(*(IMP - fs + vw - 1) == NOMASK) && (*(IMP + fs + 1) == NOMASK) &&
(*(IMP - fs - vw + 1) == NOMASK) && (*(IMP + fs + 2 * vw - 1) == NOMASK) &&
(*(IMP - fs + vw - 1) == NOMASK) && (*(IMP + fs + 1) == NOMASK) &&
(*(IMP - fs + 1) == NOMASK) && (*(IMP + fs + vw - 1) == NOMASK) &&
(*(IMP - fs + 2 * vw - 1) == NOMASK) && (*(IMP + fs - vw + 1) == NOMASK) &&
(*(IMP - fs + vw) == NOMASK) && (*(IMP + fs - vw) == NOMASK) &&
(*(IMP - fs + vw) == NOMASK) && (*(IMP + fs - vw) == NOMASK) &&
(*(IMP - fs + vw + 1) == NOMASK) && (*(IMP + fs - 1) == NOMASK) )
{
{
*EMP = NOMASK;
}
EMP += vw;
@@ -281,19 +281,19 @@ void extend_mask(unsigned char *input_mask, unsigned char *extended_mask, int vo
for (i=1; i < volume_height - 1; i++)
{
if( (*IMP) == NOMASK && (*(IMP - vw + 1) == NOMASK) && (*(IMP - 1) == NOMASK) &&
(*(IMP - vw) == NOMASK) && (*(IMP + vw) == NOMASK) &&
(*(IMP - vw) == NOMASK) && (*(IMP + vw) == NOMASK) &&
(*(IMP - fs) == NOMASK) && (*(IMP + fs) == NOMASK) &&
(*(IMP - vw - 1) == NOMASK) && (*(IMP + 1) == NOMASK) &&
(*(IMP + vw - 1) == NOMASK) && (*(IMP - 2 * vw + 1) == NOMASK) &&
(*(IMP - fs - vw - 1) == NOMASK) && (*(IMP + fs + 1) == NOMASK) &&
(*(IMP + vw - 1) == NOMASK) && (*(IMP - 2 * vw + 1) == NOMASK) &&
(*(IMP - fs - vw - 1) == NOMASK) && (*(IMP + fs + 1) == NOMASK) &&
(*(IMP - fs - 2 * vw + 1) == NOMASK) && (*(IMP + fs + vw - 1) == NOMASK) &&
(*(IMP - fs - 1) == NOMASK) && (*(IMP + fs - vw + 1) == NOMASK) &&
(*(IMP - fs - vw + 1) == NOMASK) && (*(IMP + fs - 1) == NOMASK) &&
(*(IMP - fs - 1) == NOMASK) && (*(IMP + fs - vw + 1) == NOMASK) &&
(*(IMP - fs - vw + 1) == NOMASK) && (*(IMP + fs - 1) == NOMASK) &&
(*(IMP - fs - vw) == NOMASK) && (*(IMP + fs + vw) == NOMASK) &&
(*(IMP - fs + vw - 1) == NOMASK) && (*(IMP + fs - 2 * vw + 1) == NOMASK) &&
(*(IMP - fs + vw) == NOMASK) && (*(IMP + fs - vw) == NOMASK) &&
(*(IMP - fs + 1) == NOMASK) && (*(IMP + fs - vw - 1) == NOMASK) )
{
(*(IMP - fs + vw) == NOMASK) && (*(IMP + fs - vw) == NOMASK) &&
(*(IMP - fs + 1) == NOMASK) && (*(IMP + fs - vw - 1) == NOMASK) )
{
*EMP = NOMASK;
}
EMP += vw;
@@ -302,31 +302,31 @@ void extend_mask(unsigned char *input_mask, unsigned char *extended_mask, int vo
EMP += 2 * vw;
IMP += 2 *vw;
}
}
}
if (params->y_connectivity == 1)
{
//extend the mask to the left side of the phase volume
IMP = input_mask + frame_size + 1;
IMP = input_mask + frame_size + 1;
EMP = extended_mask + frame_size + 1;
for (n=1; n < volume_depth - 1; n++)
{
for (j=1; j < volume_width - 1; j++)
{
if( (*IMP) == NOMASK && (*(IMP - 1) == NOMASK) && (*(IMP + 1) == NOMASK) &&
(*(IMP + fs - vw) == NOMASK) && (*(IMP + vw) == NOMASK) &&
(*(IMP + fs - vw) == NOMASK) && (*(IMP + vw) == NOMASK) &&
(*(IMP - fs) == NOMASK) && (*(IMP + fs) == NOMASK) &&
(*(IMP + fs - vw - 1) == NOMASK) && (*(IMP + vw + 1) == NOMASK) &&
(*(IMP + fs - vw + 1) == NOMASK) && (*(IMP + vw - 1) == NOMASK) &&
(*(IMP - vw - 1) == NOMASK) && (*(IMP + fs + vw + 1) == NOMASK) &&
(*(IMP + fs - vw + 1) == NOMASK) && (*(IMP + vw - 1) == NOMASK) &&
(*(IMP - vw - 1) == NOMASK) && (*(IMP + fs + vw + 1) == NOMASK) &&
(*(IMP - vw) == NOMASK) && (*(IMP + fs + vw) == NOMASK) &&
(*(IMP - vw + 1) == NOMASK) && (*(IMP + fs + vw - 1) == NOMASK) &&
(*(IMP - fs - 1) == NOMASK) && (*(IMP + fs + 1) == NOMASK) &&
(*(IMP - vw + 1) == NOMASK) && (*(IMP + fs + vw - 1) == NOMASK) &&
(*(IMP - fs - 1) == NOMASK) && (*(IMP + fs + 1) == NOMASK) &&
(*(IMP - fs + 1) == NOMASK) && (*(IMP + fs - 1) == NOMASK) &&
(*(IMP - fs + vw - 1) == NOMASK) && (*(IMP + 2 * fs - vw + 1) == NOMASK) &&
(*(IMP - fs + vw) == NOMASK) && (*(IMP + 2 * fs - vw) == NOMASK) &&
(*(IMP - fs + vw) == NOMASK) && (*(IMP + 2 * fs - vw) == NOMASK) &&
(*(IMP - fs + vw + 1) == NOMASK) && (*(IMP + 2 * fs - vw - 1) == NOMASK) )
{
{
*EMP = NOMASK;
}
EMP++;
@@ -344,19 +344,19 @@ void extend_mask(unsigned char *input_mask, unsigned char *extended_mask, int vo
for (j=1; j < volume_width - 1; j++)
{
if( (*IMP) == NOMASK && (*(IMP + 1) == NOMASK) && (*(IMP - 1) == NOMASK) &&
(*(IMP - vw) == NOMASK) && (*(IMP - fs + vw) == NOMASK) &&
(*(IMP - vw) == NOMASK) && (*(IMP - fs + vw) == NOMASK) &&
(*(IMP - fs) == NOMASK) && (*(IMP + fs) == NOMASK) &&
(*(IMP - vw - 1) == NOMASK) && (*(IMP - fs + vw + 1) == NOMASK) &&
(*(IMP - vw + 1) == NOMASK) && (*(IMP - fs + vw - 1) == NOMASK) &&
(*(IMP - fs - vw - 1) == NOMASK) && (*(IMP + vw + 1) == NOMASK) &&
(*(IMP - vw + 1) == NOMASK) && (*(IMP - fs + vw - 1) == NOMASK) &&
(*(IMP - fs - vw - 1) == NOMASK) && (*(IMP + vw + 1) == NOMASK) &&
(*(IMP - fs - vw + 1) == NOMASK) && (*(IMP + vw - 1) == NOMASK) &&
(*(IMP - fs - vw) == NOMASK) && (*(IMP + vw) == NOMASK) &&
(*(IMP - fs - 1) == NOMASK) && (*(IMP + fs + 1) == NOMASK) &&
(*(IMP - fs - vw) == NOMASK) && (*(IMP + vw) == NOMASK) &&
(*(IMP - fs - 1) == NOMASK) && (*(IMP + fs + 1) == NOMASK) &&
(*(IMP - fs + 1) == NOMASK) && (*(IMP + fs - 1) == NOMASK) &&
(*(IMP - 2 * fs + vw - 1) == NOMASK) && (*(IMP + fs - vw + 1) == NOMASK) &&
(*(IMP - 2 * fs + vw) == NOMASK) && (*(IMP + fs - vw) == NOMASK) &&
(*(IMP - 2 * fs + vw + 1) == NOMASK) && (*(IMP + fs - vw - 1) == NOMASK) )
{
(*(IMP - 2 * fs + vw) == NOMASK) && (*(IMP + fs - vw) == NOMASK) &&
(*(IMP - 2 * fs + vw + 1) == NOMASK) && (*(IMP + fs - vw - 1) == NOMASK) )
{
*EMP = NOMASK;
}
EMP++;
@@ -370,26 +370,26 @@ void extend_mask(unsigned char *input_mask, unsigned char *extended_mask, int vo
if (params->z_connectivity == 1)
{
//extend the mask to the bottom side of the phase volume
IMP = input_mask + volume_width + 1;
IMP = input_mask + volume_width + 1;
EMP = extended_mask + volume_width + 1;
for (i=1; i < volume_height - 1; ++i)
{
for (j=1; j < volume_width - 1; ++j)
{
if( (*IMP) == NOMASK && (*(IMP - 1) == NOMASK) && (*(IMP + 1) == NOMASK) &&
(*(IMP - vw) == NOMASK) && (*(IMP + vw) == NOMASK) &&
(*(IMP - vw) == NOMASK) && (*(IMP + vw) == NOMASK) &&
(*(IMP + fs) == NOMASK) && (*(IMP + vs - fs) == NOMASK) &&
(*(IMP - vw - 1) == NOMASK) && (*(IMP + vw + 1) == NOMASK) &&
(*(IMP - vw + 1) == NOMASK) && (*(IMP + vw - 1) == NOMASK) &&
(*(IMP + vs - fs - vw - 1) == NOMASK) && (*(IMP + fs + vw + 1) == NOMASK) &&
(*(IMP - vw + 1) == NOMASK) && (*(IMP + vw - 1) == NOMASK) &&
(*(IMP + vs - fs - vw - 1) == NOMASK) && (*(IMP + fs + vw + 1) == NOMASK) &&
(*(IMP + vs - fs - vw) == NOMASK) && (*(IMP + fs + vw) == NOMASK) &&
(*(IMP + vs - fs - vw + 1) == NOMASK) && (*(IMP + fs + vw - 1) == NOMASK) &&
(*(IMP + vs - fs - 1) == NOMASK) && (*(IMP + fs + 1) == NOMASK) &&
(*(IMP + vs - fs - vw + 1) == NOMASK) && (*(IMP + fs + vw - 1) == NOMASK) &&
(*(IMP + vs - fs - 1) == NOMASK) && (*(IMP + fs + 1) == NOMASK) &&
(*(IMP + vs - fs + 1) == NOMASK) && (*(IMP + fs - 1) == NOMASK) &&
(*(IMP + vs - fs + vw - 1) == NOMASK) && (*(IMP + fs - vw + 1) == NOMASK) &&
(*(IMP + vs - fs + vw) == NOMASK) && (*(IMP + fs - vw) == NOMASK) &&
(*(IMP + vs - fs + vw) == NOMASK) && (*(IMP + fs - vw) == NOMASK) &&
(*(IMP + vs - fs + vw + 1) == NOMASK) && (*(IMP + fs - vw - 1) == NOMASK) )
{
{
*EMP = NOMASK;
}
EMP++;
@@ -400,26 +400,26 @@ void extend_mask(unsigned char *input_mask, unsigned char *extended_mask, int vo
}
//extend the mask to the top side of the phase volume
IMP = input_mask + volume_size - frame_size + volume_width + 1;
IMP = input_mask + volume_size - frame_size + volume_width + 1;
EMP = extended_mask + volume_size - frame_size + volume_width + 1;
for (i=1; i < volume_height - 1; ++i)
{
for (j=1; j < volume_width - 1; ++j)
{
if( (*IMP) == NOMASK && (*(IMP + 1) == NOMASK) && (*(IMP - 1) == NOMASK) &&
(*(IMP - vw) == NOMASK) && (*(IMP - fs + vw) == NOMASK) &&
(*(IMP - vw) == NOMASK) && (*(IMP - fs + vw) == NOMASK) &&
(*(IMP - fs) == NOMASK) && (*(IMP - vs + fs) == NOMASK) &&
(*(IMP - vw - 1) == NOMASK) && (*(IMP + vw + 1) == NOMASK) &&
(*(IMP - vw + 1) == NOMASK) && (*(IMP + vw - 1) == NOMASK) &&
(*(IMP - fs - vw - 1) == NOMASK) && (*(IMP - vs + fs + vw + 1) == NOMASK) &&
(*(IMP - vw + 1) == NOMASK) && (*(IMP + vw - 1) == NOMASK) &&
(*(IMP - fs - vw - 1) == NOMASK) && (*(IMP - vs + fs + vw + 1) == NOMASK) &&
(*(IMP - fs - vw + 1) == NOMASK) && (*(IMP - vs + fs + vw - 1) == NOMASK) &&
(*(IMP - fs - vw) == NOMASK) && (*(IMP - vs + fs + vw) == NOMASK) &&
(*(IMP - fs - 1) == NOMASK) && (*(IMP - vs + fs + 1) == NOMASK) &&
(*(IMP - fs - vw) == NOMASK) && (*(IMP - vs + fs + vw) == NOMASK) &&
(*(IMP - fs - 1) == NOMASK) && (*(IMP - vs + fs + 1) == NOMASK) &&
(*(IMP - fs + 1) == NOMASK) && (*(IMP - vs + fs - 1) == NOMASK) &&
(*(IMP - fs + vw - 1) == NOMASK) && (*(IMP - vs + fs - vw + 1) == NOMASK) &&
(*(IMP - fs + vw) == NOMASK) && (*(IMP - vs + fs - vw) == NOMASK) &&
(*(IMP - fs + vw + 1) == NOMASK) && (*(IMP - vs + fs - vw - 1) == NOMASK) )
{
(*(IMP - fs + vw) == NOMASK) && (*(IMP - vs + fs - vw) == NOMASK) &&
(*(IMP - fs + vw + 1) == NOMASK) && (*(IMP - vs + fs - vw - 1) == NOMASK) )
{
*EMP = NOMASK;
}
EMP++;
@@ -432,14 +432,14 @@ void extend_mask(unsigned char *input_mask, unsigned char *extended_mask, int vo
}
void calculate_reliability(float *wrappedVolume, VOXELM *voxel, int volume_width, int volume_height, int volume_depth, params_t *params)
{
{
int frame_size = volume_width * volume_height;
int volume_size = volume_width * volume_height * volume_depth;
VOXELM *voxel_pointer;
float H, V, N, D1, D2, D3, D4, D5, D6, D7, D8, D9, D10;
float *WVP;
int n, i, j;
WVP = wrappedVolume + frame_size + volume_width + 1;
voxel_pointer = voxel + frame_size + volume_width + 1;
for (n=1; n < volume_depth - 1; n++)
@@ -449,7 +449,7 @@ void calculate_reliability(float *wrappedVolume, VOXELM *voxel, int volume_width
for (j=1; j < volume_width - 1; j++)
{
if (voxel_pointer->extended_mask == NOMASK)
{
{
H = wrap(*(WVP - 1) - *WVP) - wrap(*WVP - *(WVP + 1));
V = wrap(*(WVP - volume_width) - *WVP) - wrap(*WVP - *(WVP + volume_width));
N = wrap(*(WVP - frame_size) - *WVP) - wrap(*WVP - *(WVP + frame_size));
@@ -463,7 +463,7 @@ void calculate_reliability(float *wrappedVolume, VOXELM *voxel, int volume_width
D8 = wrap(*(WVP - frame_size + volume_width - 1) - *WVP) - wrap(*WVP - *(WVP + frame_size - volume_width + 1));
D9 = wrap(*(WVP - frame_size + volume_width) - *WVP) - wrap(*WVP - *(WVP + frame_size - volume_width));
D10 = wrap(*(WVP - frame_size + volume_width + 1) - *WVP) - wrap(*WVP - *(WVP + frame_size - volume_width - 1));
voxel_pointer->reliability = H*H + V*V + N*N + D1*D1 + D2*D2 + D3*D3 + D4*D4 + D5*D5 + D6*D6
voxel_pointer->reliability = H*H + V*V + N*N + D1*D1 + D2*D2 + D3*D3 + D4*D4 + D5*D5 + D6*D6
+ D7*D7 + D8*D8 + D9*D9 + D10*D10;
}
voxel_pointer++;
@@ -486,7 +486,7 @@ void calculate_reliability(float *wrappedVolume, VOXELM *voxel, int volume_width
for (i=1; i < volume_height - 1; ++i)
{
if (voxel_pointer->extended_mask == NOMASK)
{
{
H = wrap(*(WVP + volume_width - 1) - *WVP) - wrap(*WVP - *(WVP + 1));
V = wrap(*(WVP - volume_width) - *WVP) - wrap(*WVP - *(WVP + volume_width));
N = wrap(*(WVP - frame_size) - *WVP) - wrap(*WVP - *(WVP + frame_size));
@@ -500,7 +500,7 @@ void calculate_reliability(float *wrappedVolume, VOXELM *voxel, int volume_width
D8 = wrap(*(WVP - frame_size + 2 * volume_width - 1) - *WVP) - wrap(*WVP - *(WVP + frame_size - volume_width + 1));
D9 = wrap(*(WVP - frame_size + volume_width) - *WVP) - wrap(*WVP - *(WVP + frame_size - volume_width));
D10 = wrap(*(WVP - frame_size + volume_width + 1) - *WVP) - wrap(*WVP - *(WVP + frame_size - 1));
voxel_pointer->reliability = H*H + V*V + N*N + D1*D1 + D2*D2 + D3*D3 + D4*D4 + D5*D5 + D6*D6
voxel_pointer->reliability = H*H + V*V + N*N + D1*D1 + D2*D2 + D3*D3 + D4*D4 + D5*D5 + D6*D6
+ D7*D7 + D8*D8 + D9*D9 + D10*D10;
}
voxel_pointer += volume_width;
@@ -509,7 +509,7 @@ void calculate_reliability(float *wrappedVolume, VOXELM *voxel, int volume_width
voxel_pointer += 2 * volume_width;
WVP += 2 * volume_width;
}
//calculating reliability for the rear side of the phase volume..... subtract volume_width
WVP = wrappedVolume + frame_size + 2 * volume_width - 1;
voxel_pointer = voxel + frame_size + 2 * volume_width - 1;
@@ -518,7 +518,7 @@ void calculate_reliability(float *wrappedVolume, VOXELM *voxel, int volume_width
for (i=1; i < volume_height - 1; ++i)
{
if (voxel_pointer->extended_mask == NOMASK)
{
{
H = wrap(*(WVP - volume_width + 1) - *WVP) - wrap(*WVP - *(WVP - 1));
V = wrap(*(WVP - volume_width) - *WVP) - wrap(*WVP - *(WVP + volume_width));
N = wrap(*(WVP - frame_size) - *WVP) - wrap(*WVP - *(WVP + frame_size));
@@ -532,7 +532,7 @@ void calculate_reliability(float *wrappedVolume, VOXELM *voxel, int volume_width
D8 = wrap(*(WVP - frame_size + volume_width - 1) - *WVP) - wrap(*WVP - *(WVP + frame_size - 2 * volume_width + 1));
D9 = wrap(*(WVP - frame_size + volume_width) - *WVP) - wrap(*WVP - *(WVP + frame_size - volume_width));
D10 = wrap(*(WVP - frame_size + 1) - *WVP) - wrap(*WVP - *(WVP + frame_size - volume_width - 1));
voxel_pointer->reliability = H*H + V*V + N*N + D1*D1 + D2*D2 + D3*D3 + D4*D4 + D5*D5 + D6*D6
voxel_pointer->reliability = H*H + V*V + N*N + D1*D1 + D2*D2 + D3*D3 + D4*D4 + D5*D5 + D6*D6
+ D7*D7 + D8*D8 + D9*D9 + D10*D10;
}
voxel_pointer += volume_width;
@@ -540,7 +540,7 @@ void calculate_reliability(float *wrappedVolume, VOXELM *voxel, int volume_width
}
voxel_pointer += 2 * volume_width;
WVP += 2 * volume_width;
}
}
}
if (params->y_connectivity == 1)
@@ -553,7 +553,7 @@ void calculate_reliability(float *wrappedVolume, VOXELM *voxel, int volume_width
for (j=1; j < volume_width - 1; ++j)
{
if (voxel_pointer->extended_mask == NOMASK)
{
{
H = wrap(*(WVP - 1) - *WVP) - wrap(*WVP - *(WVP + 1));
V = wrap(*(WVP + frame_size - volume_width) - *WVP) - wrap(*WVP - *(WVP + volume_width));
N = wrap(*(WVP - frame_size) - *WVP) - wrap(*WVP - *(WVP + frame_size));
@@ -567,7 +567,7 @@ void calculate_reliability(float *wrappedVolume, VOXELM *voxel, int volume_width
D8 = wrap(*(WVP - frame_size + volume_width - 1) - *WVP) - wrap(*WVP - *(WVP + 2 * frame_size - volume_width + 1));
D9 = wrap(*(WVP - frame_size + volume_width) - *WVP) - wrap(*WVP - *(WVP + 2 * frame_size - volume_width));
D10 = wrap(*(WVP - frame_size + volume_width + 1) - *WVP) - wrap(*WVP - *(WVP + 2 * frame_size - volume_width - 1));
voxel_pointer->reliability = H*H + V*V + N*N + D1*D1 + D2*D2 + D3*D3 + D4*D4 + D5*D5 + D6*D6
voxel_pointer->reliability = H*H + V*V + N*N + D1*D1 + D2*D2 + D3*D3 + D4*D4 + D5*D5 + D6*D6
+ D7*D7 + D8*D8 + D9*D9 + D10*D10;
}
voxel_pointer++;
@@ -575,7 +575,7 @@ void calculate_reliability(float *wrappedVolume, VOXELM *voxel, int volume_width
}
voxel_pointer += frame_size - volume_width + 2;
WVP += frame_size - volume_width + 2;
}
}
//calculating reliability for the right side of the phase volume...subtract frame_size
WVP = wrappedVolume + 2 * frame_size - volume_width + 1;
@@ -585,7 +585,7 @@ void calculate_reliability(float *wrappedVolume, VOXELM *voxel, int volume_width
for (j=1; j < volume_width - 1; ++j)
{
if (voxel_pointer->extended_mask == NOMASK)
{
{
H = wrap(*(WVP + 1) - *WVP) - wrap(*WVP - *(WVP - 1));
V = wrap(*(WVP - volume_width) - *WVP) - wrap(*WVP - *(WVP - frame_size + volume_width));
N = wrap(*(WVP - frame_size) - *WVP) - wrap(*WVP - *(WVP + frame_size));
@@ -599,7 +599,7 @@ void calculate_reliability(float *wrappedVolume, VOXELM *voxel, int volume_width
D8 = wrap(*(WVP - 2 * frame_size + volume_width - 1) - *WVP) - wrap(*WVP - *(WVP + frame_size - volume_width + 1));
D9 = wrap(*(WVP - 2 * frame_size + volume_width) - *WVP) - wrap(*WVP - *(WVP + frame_size - volume_width));
D10 = wrap(*(WVP - 2 * frame_size + volume_width + 1) - *WVP) - wrap(*WVP - *(WVP + frame_size - volume_width - 1));
voxel_pointer->reliability = H*H + V*V + N*N + D1*D1 + D2*D2 + D3*D3 + D4*D4 + D5*D5 + D6*D6
voxel_pointer->reliability = H*H + V*V + N*N + D1*D1 + D2*D2 + D3*D3 + D4*D4 + D5*D5 + D6*D6
+ D7*D7 + D8*D8 + D9*D9 + D10*D10;
}
voxel_pointer++;
@@ -607,7 +607,7 @@ void calculate_reliability(float *wrappedVolume, VOXELM *voxel, int volume_width
}
voxel_pointer += frame_size - volume_width + 2;
WVP += frame_size - volume_width + 2;
}
}
}
if (params->z_connectivity == 1)
@@ -620,7 +620,7 @@ void calculate_reliability(float *wrappedVolume, VOXELM *voxel, int volume_width
for (j=1; j < volume_width - 1; ++j)
{
if (voxel_pointer->extended_mask == NOMASK)
{
{
H = wrap(*(WVP - 1) - *WVP) - wrap(*WVP - *(WVP + 1));
V = wrap(*(WVP - volume_width) - *WVP) - wrap(*WVP - *(WVP + volume_width));
N = wrap(*(WVP + frame_size) - *WVP) - wrap(*WVP - *(WVP + volume_size - frame_size));
@@ -634,7 +634,7 @@ void calculate_reliability(float *wrappedVolume, VOXELM *voxel, int volume_width
D8 = wrap(*(WVP + volume_size - frame_size + volume_width - 1) - *WVP) - wrap(*WVP - *(WVP + frame_size - volume_width + 1));
D9 = wrap(*(WVP + volume_size - frame_size + volume_width) - *WVP) - wrap(*WVP - *(WVP + frame_size - volume_width));
D10 = wrap(*(WVP + volume_size - frame_size + volume_width + 1) - *WVP) - wrap(*WVP - *(WVP + frame_size - volume_width - 1));
voxel_pointer->reliability = H*H + V*V + N*N + D1*D1 + D2*D2 + D3*D3 + D4*D4 + D5*D5 + D6*D6
voxel_pointer->reliability = H*H + V*V + N*N + D1*D1 + D2*D2 + D3*D3 + D4*D4 + D5*D5 + D6*D6
+ D7*D7 + D8*D8 + D9*D9 + D10*D10;
}
voxel_pointer++;
@@ -642,7 +642,7 @@ void calculate_reliability(float *wrappedVolume, VOXELM *voxel, int volume_width
}
voxel_pointer += 2;
WVP += 2;
}
}
//calculating reliability for the top side of the phase volume...subtract volume_size
WVP = wrappedVolume + volume_size - frame_size + volume_width + 1;
@@ -652,7 +652,7 @@ void calculate_reliability(float *wrappedVolume, VOXELM *voxel, int volume_width
for (j=1; j < volume_width - 1; ++j)
{
if (voxel_pointer->extended_mask == NOMASK)
{
{
H = wrap(*(WVP + 1) - *WVP) - wrap(*WVP - *(WVP - 1));
V = wrap(*(WVP - volume_width) - *WVP) - wrap(*WVP - *(WVP + volume_width));
N = wrap(*(WVP - frame_size) - *WVP) - wrap(*WVP - *(WVP - volume_size + frame_size));
@@ -660,13 +660,13 @@ void calculate_reliability(float *wrappedVolume, VOXELM *voxel, int volume_width
D2 = wrap(*(WVP - volume_width + 1) - *WVP) - wrap(*WVP - *(WVP + volume_width - 1));
D3 = wrap(*(WVP - frame_size - volume_width - 1) - *WVP) - wrap(*WVP - *(WVP - volume_size + frame_size + volume_width + 1));
D4 = wrap(*(WVP - frame_size - volume_width + 1) - *WVP) - wrap(*WVP - *(WVP - volume_size + frame_size + volume_width - 1));
D5 = wrap(*(WVP - frame_size - volume_width) - *WVP) - wrap(*WVP - *(WVP - volume_size + frame_size + volume_width));
D5 = wrap(*(WVP - frame_size - volume_width) - *WVP) - wrap(*WVP - *(WVP - volume_size + frame_size + volume_width));
D6 = wrap(*(WVP - frame_size - 1) - *WVP) - wrap(*WVP - *(WVP - volume_size + frame_size + 1));
D7 = wrap(*(WVP - frame_size + 1) - *WVP) - wrap(*WVP - *(WVP - volume_size + frame_size - 1));
D8 = wrap(*(WVP - frame_size + volume_width - 1) - *WVP) - wrap(*WVP - *(WVP - volume_size + frame_size - volume_width + 1));
D9 = wrap(*(WVP - frame_size + volume_width) - *WVP) - wrap(*WVP - *(WVP - volume_size + frame_size - volume_width));
D10 = wrap(*(WVP - frame_size + volume_width + 1) - *WVP) - wrap(*WVP - *(WVP - volume_size + frame_size - volume_width - 1));
voxel_pointer->reliability = H*H + V*V + N*N + D1*D1 + D2*D2 + D3*D3 + D4*D4 + D5*D5 + D6*D6
voxel_pointer->reliability = H*H + V*V + N*N + D1*D1 + D2*D2 + D3*D3 + D4*D4 + D5*D5 + D6*D6
+ D7*D7 + D8*D8 + D9*D9 + D10*D10;
}
voxel_pointer++;
@@ -674,7 +674,7 @@ void calculate_reliability(float *wrappedVolume, VOXELM *voxel, int volume_width
}
voxel_pointer += 2;
WVP += 2;
}
}
}
}
@@ -688,12 +688,12 @@ void horizontalEDGEs(VOXELM *voxel, EDGE *edge, int volume_width, int volume_he
EDGE *edge_pointer = edge;
VOXELM *voxel_pointer = voxel;
int no_of_edges = params->no_of_edges;
for (n=0; n < volume_depth; n++)
{
for (i = 0; i < volume_height; i++)
{
for (j = 0; j < volume_width - 1; j++)
for (j = 0; j < volume_width - 1; j++)
{
if (voxel_pointer->input_mask == NOMASK && (voxel_pointer + 1)->input_mask == NOMASK )
{
@@ -737,7 +737,7 @@ void verticalEDGEs(VOXELM *voxel, EDGE *edge, int volume_width, int volume_heig
int n, i, j;
int no_of_edges = params->no_of_edges;
VOXELM *voxel_pointer = voxel;
EDGE *edge_pointer = edge + no_of_edges;
EDGE *edge_pointer = edge + no_of_edges;
int frame_size = volume_width * volume_height;
int next_voxel = frame_size - volume_width;
@@ -745,7 +745,7 @@ void verticalEDGEs(VOXELM *voxel, EDGE *edge, int volume_width, int volume_heig
{
for (i=0; i<volume_height - 1; i++)
{
for (j=0; j < volume_width; j++)
for (j=0; j < volume_width; j++)
{
if (voxel_pointer->input_mask == NOMASK && (voxel_pointer + volume_width)->input_mask == NOMASK )
{
@@ -760,7 +760,7 @@ void verticalEDGEs(VOXELM *voxel, EDGE *edge, int volume_width, int volume_heig
}
}
voxel_pointer += volume_width;
}
}
if (params->y_connectivity == 1)
{
@@ -800,7 +800,7 @@ void normalEDGEs(VOXELM *voxel, EDGE *edge, int volume_width, int volume_height
{
for (i=0; i<volume_height; i++)
{
for (j=0; j < volume_width; j++)
for (j=0; j < volume_width; j++)
{
if (voxel_pointer->input_mask == NOMASK && (voxel_pointer + frame_size)->input_mask == NOMASK )
{
@@ -816,7 +816,7 @@ void normalEDGEs(VOXELM *voxel, EDGE *edge, int volume_width, int volume_height
}
}
if (params->z_connectivity == 1)
{
voxel_pointer = voxel + next_voxel;
@@ -840,11 +840,11 @@ void normalEDGEs(VOXELM *voxel, EDGE *edge, int volume_width, int volume_height
params->no_of_edges = no_of_edges;
}
//gather the voxels of the volume into groups
//gather the voxels of the volume into groups
void gatherVOXELs(EDGE *edge, params_t *params)
{
int k;
VOXELM *VOXEL1;
VOXELM *VOXEL1;
VOXELM *VOXEL2;
VOXELM *group1;
VOXELM *group2;
@@ -862,7 +862,7 @@ void gatherVOXELs(EDGE *edge, params_t *params)
if (VOXEL2->head != VOXEL1->head)
{
//VOXELM 2 is alone in its group
//merge this voxel with VOXELM 1 group and find the number of 2 pi to add
//merge this voxel with VOXELM 1 group and find the number of 2 pi to add
//to or subtract to unwrap it
if ((VOXEL2->next == NULL) && (VOXEL2->head == VOXEL2))
{
@@ -874,7 +874,7 @@ void gatherVOXELs(EDGE *edge, params_t *params)
}
//VOXELM 1 is alone in its group
//merge this voxel with VOXELM 2 group and find the number of 2 pi to add
//merge this voxel with VOXELM 2 group and find the number of 2 pi to add
//to or subtract to unwrap it
else if ((VOXEL1->next == NULL) && (VOXEL1->head == VOXEL1))
{
@@ -883,7 +883,7 @@ void gatherVOXELs(EDGE *edge, params_t *params)
(VOXEL2->head->number_of_voxels_in_group)++;
VOXEL1->head = VOXEL2->head;
VOXEL1->increment = VOXEL2->increment+pointer_edge->increment;
}
}
//VOXELM 1 and VOXELM 2 both have groups
else
@@ -909,7 +909,7 @@ void gatherVOXELs(EDGE *edge, params_t *params)
group2->increment += incremento;
group2 = group2->next;
}
}
}
//if the no. of voxels in VOXELM 2 group is larger than the no. of voxels
//in VOXELM 1 group. Merge VOXELM 1 group to VOXELM 2 group
@@ -938,7 +938,7 @@ void gatherVOXELs(EDGE *edge, params_t *params)
}
}
//unwrap the volume
//unwrap the volume
void unwrapVolume(VOXELM *voxel, int volume_width, int volume_height, int volume_depth)
{
int i;
@@ -969,7 +969,7 @@ void maskVolume(VOXELM *voxel, unsigned char *input_mask, int volume_width, int
//find the minimum of the unwrapped phase
for (i = 0; i < volume_size; i++)
{
if ((pointer_voxel->value < min) && (*IMP == NOMASK))
if ((pointer_voxel->value < min) && (*IMP == NOMASK))
min = pointer_voxel->value;
pointer_voxel++;
@@ -977,7 +977,7 @@ void maskVolume(VOXELM *voxel, unsigned char *input_mask, int volume_width, int
}
pointer_voxel = voxel;
IMP = input_mask;
IMP = input_mask;
//set the masked voxels to minimum
for (i = 0; i < volume_size; i++)
@@ -1002,7 +1002,7 @@ void returnVolume(VOXELM *voxel, float *unwrappedVolume, int volume_width, int
float *unwrappedVolume_pointer = unwrappedVolume;
VOXELM *voxel_pointer = voxel;
for (i=0; i < volume_size; i++)
for (i=0; i < volume_size; i++)
{
*unwrappedVolume_pointer = voxel_pointer->value;
voxel_pointer++;
skimage/exposure/Miguel_2D_unwrapper_with_mask_and_wrap_around_option.c
+63 -63
View File
@@ -1,15 +1,15 @@
// 2D phase unwrapping, modified for inclusion in scipy by Gregor Thalhammer
//This program was written by Munther Gdeisat and Miguel Arevallilo Herraez to program the two-dimensional unwrapper
//entitled "Fast two-dimensional phase-unwrapping algorithm based on sorting by
//entitled "Fast two-dimensional phase-unwrapping algorithm based on sorting by
//reliability following a noncontinuous path"
//by Miguel Arevallilo Herraez, David R. Burton, Michael J. Lalor, and Munther A. Gdeisat
//published in the Journal Applied Optics, Vol. 41, No. 35, pp. 7437, 2002.
//This program was written by Munther Gdeisat, Liverpool John Moores University, United Kingdom.
//Date 26th August 2007
//The wrapped phase map is assumed to be of floating point data type. The resultant unwrapped phase map is also of floating point type.
//The mask is of byte data type.
//When the mask is 255 this means that the pixel is valid
//The mask is of byte data type.
//When the mask is 255 this means that the pixel is valid
//When the mask is 0 this means that the pixel is invalid (noisy or corrupted pixel)
//This program takes into consideration the image wrap around problem encountered in MRI imaging.
@@ -55,13 +55,13 @@ typedef struct PIXELM PIXELM;
//the EDGE is the line that connects two pixels.
//if we have S pixels, then we have S horizontal edges and S vertical edges
struct EDGE
{
{
float reliab; //reliabilty of the edge and it depends on the two pixels
PIXELM *pointer_1; //pointer to the first pixel
PIXELM *pointer_2; //pointer to the second pixel
int increment; //No. of 2*pi to add to one of the pixels to
//unwrap it with respect to the second
};
};
typedef struct EDGE EDGE;
@@ -147,20 +147,20 @@ void initialisePIXELs(float *wrapped_image, unsigned char *input_mask, unsigned
unsigned char *input_mask_pointer = input_mask;
unsigned char *extended_mask_pointer = extended_mask;
int i, j;
for (i=0; i < image_height; i++)
{
for (j=0; j < image_width; j++)
{
pixel_pointer->increment = 0;
pixel_pointer->number_of_pixels_in_group = 1;
pixel_pointer->number_of_pixels_in_group = 1;
pixel_pointer->value = *wrapped_image_pointer;
pixel_pointer->reliability = 9999999.f + rand();
pixel_pointer->input_mask = *input_mask_pointer;
pixel_pointer->extended_mask = *extended_mask_pointer;
pixel_pointer->head = pixel_pointer;
pixel_pointer->last = pixel_pointer;
pixel_pointer->next = NULL;
pixel_pointer->next = NULL;
pixel_pointer->new_group = 0;
pixel_pointer->group = -1;
pixel_pointer++;
@@ -185,7 +185,7 @@ float wrap(float pixel_value)
// pixelL_value is the left pixel, pixelR_value is the right pixel
int find_wrap(float pixelL_value, float pixelR_value)
{
float difference;
float difference;
int wrap_value;
difference = pixelL_value - pixelR_value;
@@ -194,10 +194,10 @@ int find_wrap(float pixelL_value, float pixelR_value)
else wrap_value = 0;
return wrap_value;
}
}
void extend_mask(unsigned char *input_mask, unsigned char *extended_mask,
int image_width, int image_height,
void extend_mask(unsigned char *input_mask, unsigned char *extended_mask,
int image_width, int image_height,
params_t *params)
{
int i,j;
@@ -211,11 +211,11 @@ void extend_mask(unsigned char *input_mask, unsigned char *extended_mask,
{
for (j=1; j < image_width - 1; ++j)
{
if ( (*IMP) == NOMASK && (*(IMP + 1) == NOMASK) && (*(IMP - 1) == NOMASK) &&
if ( (*IMP) == NOMASK && (*(IMP + 1) == NOMASK) && (*(IMP - 1) == NOMASK) &&
(*(IMP + image_width) == NOMASK) && (*(IMP - image_width) == NOMASK) &&
(*(IMP - image_width_minus_one) == NOMASK) && (*(IMP - image_width_plus_one) == NOMASK) &&
(*(IMP + image_width_minus_one) == NOMASK) && (*(IMP + image_width_plus_one) == NOMASK) )
{
{
*EMP = NOMASK;
}
++EMP;
@@ -224,7 +224,7 @@ void extend_mask(unsigned char *input_mask, unsigned char *extended_mask,
EMP += 2;
IMP += 2;
}
if (params->x_connectivity == 1)
{
//extend the mask for the right border of the image
@@ -236,7 +236,7 @@ void extend_mask(unsigned char *input_mask, unsigned char *extended_mask,
(*(IMP + image_width) == NOMASK) && (*(IMP - image_width) == NOMASK) &&
(*(IMP - image_width - 1) == NOMASK) && (*(IMP - image_width + 1) == NOMASK) &&
(*(IMP + image_width - 1) == NOMASK) && (*(IMP - 2 * image_width + 1) == NOMASK) )
{
{
*EMP = NOMASK;
}
EMP += image_width;
@@ -248,11 +248,11 @@ void extend_mask(unsigned char *input_mask, unsigned char *extended_mask,
EMP = extended_mask + image_width;
for (i=1; i < image_height - 1; ++i)
{
if ( (*IMP) == NOMASK && (*(IMP - 1) == NOMASK) && (*(IMP + 1) == NOMASK) &&
if ( (*IMP) == NOMASK && (*(IMP - 1) == NOMASK) && (*(IMP + 1) == NOMASK) &&
(*(IMP + image_width) == NOMASK) && (*(IMP - image_width) == NOMASK) &&
(*(IMP - image_width + 1) == NOMASK) && (*(IMP + image_width + 1) == NOMASK) &&
(*(IMP + image_width - 1) == NOMASK) && (*(IMP + 2 * image_width - 1) == NOMASK) )
{
{
*EMP = NOMASK;
}
EMP += image_width;
@@ -267,11 +267,11 @@ void extend_mask(unsigned char *input_mask, unsigned char *extended_mask,
EMP = extended_mask + 1;
for (i=1; i < image_width - 1; ++i)
{
if ( (*IMP) == NOMASK && (*(IMP - 1) == NOMASK) && (*(IMP + 1) == NOMASK) &&
if ( (*IMP) == NOMASK && (*(IMP - 1) == NOMASK) && (*(IMP + 1) == NOMASK) &&
(*(IMP + image_width) == NOMASK) && (*(IMP + image_width * (image_height - 1)) == NOMASK) &&
(*(IMP + image_width + 1) == NOMASK) && (*(IMP + image_width - 1) == NOMASK) &&
(*(IMP + image_width * (image_height - 1) - 1) == NOMASK) && (*(IMP + image_width * (image_height - 1) + 1) == NOMASK) )
{
{
*EMP = NOMASK;
}
EMP++;
@@ -283,21 +283,21 @@ void extend_mask(unsigned char *input_mask, unsigned char *extended_mask,
EMP = extended_mask + image_width * (image_height - 1) + 1;
for (i=1; i < image_width - 1; ++i)
{
if ( (*IMP) == NOMASK && (*(IMP - 1) == NOMASK) && (*(IMP + 1) == NOMASK) &&
if ( (*IMP) == NOMASK && (*(IMP - 1) == NOMASK) && (*(IMP + 1) == NOMASK) &&
(*(IMP - image_width) == NOMASK) && (*(IMP - image_width - 1) == NOMASK) && (*(IMP - image_width + 1) == NOMASK) &&
(*(IMP - image_width * (image_height - 1) ) == NOMASK) &&
(*(IMP - image_width * (image_height - 1) - 1) == NOMASK) &&
(*(IMP - image_width * (image_height - 1) ) == NOMASK) &&
(*(IMP - image_width * (image_height - 1) - 1) == NOMASK) &&
(*(IMP - image_width * (image_height - 1) + 1) == NOMASK) )
{
{
*EMP = NOMASK;
}
EMP++;
IMP++;
}
}
}
}
void calculate_reliability(float *wrappedImage, PIXELM *pixel,
void calculate_reliability(float *wrappedImage, PIXELM *pixel,
int image_width, int image_height,
params_t *params)
{
@@ -307,7 +307,7 @@ void calculate_reliability(float *wrappedImage, PIXELM *pixel,
float *WIP = wrappedImage + image_width_plus_one; //WIP is the wrapped image pointer
float H, V, D1, D2;
int i, j;
for (i = 1; i < image_height -1; ++i)
{
for (j = 1; j < image_width - 1; ++j)
@@ -331,8 +331,8 @@ void calculate_reliability(float *wrappedImage, PIXELM *pixel,
{
//calculating the reliability for the left border of the image
PIXELM *pixel_pointer = pixel + image_width;
float *WIP = wrappedImage + image_width;
float *WIP = wrappedImage + image_width;
for (i = 1; i < image_height - 1; ++i)
{
if (pixel_pointer->extended_mask == NOMASK)
@@ -349,8 +349,8 @@ void calculate_reliability(float *wrappedImage, PIXELM *pixel,
//calculating the reliability for the right border of the image
pixel_pointer = pixel + 2 * image_width - 1;
WIP = wrappedImage + 2 * image_width - 1;
WIP = wrappedImage + 2 * image_width - 1;
for (i = 1; i < image_height - 1; ++i)
{
if (pixel_pointer->extended_mask == NOMASK)
@@ -370,8 +370,8 @@ void calculate_reliability(float *wrappedImage, PIXELM *pixel,
{
//calculating the reliability for the top border of the image
PIXELM *pixel_pointer = pixel + 1;
float *WIP = wrappedImage + 1;
float *WIP = wrappedImage + 1;
for (i = 1; i < image_width - 1; ++i)
{
if (pixel_pointer->extended_mask == NOMASK)
@@ -388,8 +388,8 @@ void calculate_reliability(float *wrappedImage, PIXELM *pixel,
//calculating the reliability for the bottom border of the image
pixel_pointer = pixel + (image_height - 1) * image_width + 1;
WIP = wrappedImage + (image_height - 1) * image_width + 1;
WIP = wrappedImage + (image_height - 1) * image_width + 1;
for (i = 1; i < image_width - 1; ++i)
{
if (pixel_pointer->extended_mask == NOMASK)
@@ -407,10 +407,10 @@ void calculate_reliability(float *wrappedImage, PIXELM *pixel,
}
//calculate the reliability of the horizontal edges of the image
//it is calculated by adding the reliability of pixel and the relibility of
//it is calculated by adding the reliability of pixel and the relibility of
//its right-hand neighbour
//edge is calculated between a pixel and its next neighbour
void horizontalEDGEs(PIXELM *pixel, EDGE *edge,
void horizontalEDGEs(PIXELM *pixel, EDGE *edge,
int image_width, int image_height,
params_t *params)
{
@@ -418,10 +418,10 @@ void horizontalEDGEs(PIXELM *pixel, EDGE *edge,
EDGE *edge_pointer = edge;
PIXELM *pixel_pointer = pixel;
int no_of_edges = params->no_of_edges;
for (i = 0; i < image_height; i++)
{
for (j = 0; j < image_width - 1; j++)
for (j = 0; j < image_width - 1; j++)
{
if (pixel_pointer->input_mask == NOMASK && (pixel_pointer + 1)->input_mask == NOMASK)
{
@@ -458,20 +458,20 @@ void horizontalEDGEs(PIXELM *pixel, EDGE *edge,
}
//calculate the reliability of the vertical edges of the image
//it is calculated by adding the reliability of pixel and the relibility of
//it is calculated by adding the reliability of pixel and the relibility of
//its lower neighbour in the image.
void verticalEDGEs(PIXELM *pixel, EDGE *edge,
void verticalEDGEs(PIXELM *pixel, EDGE *edge,
int image_width, int image_height,
params_t *params)
{
int i, j;
int no_of_edges = params->no_of_edges;
PIXELM *pixel_pointer = pixel;
EDGE *edge_pointer = edge + no_of_edges;
EDGE *edge_pointer = edge + no_of_edges;
for (i=0; i < image_height - 1; i++)
{
for (j=0; j < image_width; j++)
for (j=0; j < image_width; j++)
{
if (pixel_pointer->input_mask == NOMASK && (pixel_pointer + image_width)->input_mask == NOMASK)
{
@@ -507,11 +507,11 @@ void verticalEDGEs(PIXELM *pixel, EDGE *edge,
params->no_of_edges = no_of_edges;
}
//gather the pixels of the image into groups
//gather the pixels of the image into groups
void gatherPIXELs(EDGE *edge, params_t *params)
{
int k;
PIXELM *PIXEL1;
PIXELM *PIXEL1;
PIXELM *PIXEL2;
PIXELM *group1;
PIXELM *group2;
@@ -529,7 +529,7 @@ void gatherPIXELs(EDGE *edge, params_t *params)
if (PIXEL2->head != PIXEL1->head)
{
//PIXELM 2 is alone in its group
//merge this pixel with PIXELM 1 group and find the number of 2 pi to add
//merge this pixel with PIXELM 1 group and find the number of 2 pi to add
//to or subtract to unwrap it
if ((PIXEL2->next == NULL) && (PIXEL2->head == PIXEL2))
{
@@ -541,7 +541,7 @@ void gatherPIXELs(EDGE *edge, params_t *params)
}
//PIXELM 1 is alone in its group
//merge this pixel with PIXELM 2 group and find the number of 2 pi to add
//merge this pixel with PIXELM 2 group and find the number of 2 pi to add
//to or subtract to unwrap it
else if ((PIXEL1->next == NULL) && (PIXEL1->head == PIXEL1))
{
@@ -550,7 +550,7 @@ void gatherPIXELs(EDGE *edge, params_t *params)
(PIXEL2->head->number_of_pixels_in_group)++;
PIXEL1->head = PIXEL2->head;
PIXEL1->increment = PIXEL2->increment+pointer_edge->increment;
}
}
//PIXELM 1 and PIXELM 2 both have groups
else
@@ -577,7 +577,7 @@ void gatherPIXELs(EDGE *edge, params_t *params)
group2->increment += incremento;
group2 = group2->next;
}
}
}
//if the no. of pixels in PIXELM 2 group is larger than the
//no. of pixels in PIXELM 1 group. Merge PIXELM 1 group to
@@ -607,7 +607,7 @@ void gatherPIXELs(EDGE *edge, params_t *params)
}
}
//unwrap the image
//unwrap the image
void unwrapImage(PIXELM *pixel, int image_width, int image_height)
{
int i;
@@ -638,7 +638,7 @@ void maskImage(PIXELM *pixel, unsigned char *input_mask, int image_width, int i
//find the minimum of the unwrapped phase
for (i = 0; i < image_size; i++)
{
if ((pointer_pixel->value < min) && (*IMP == NOMASK))
if ((pointer_pixel->value < min) && (*IMP == NOMASK))
min = pointer_pixel->value;
pointer_pixel++;
@@ -646,7 +646,7 @@ void maskImage(PIXELM *pixel, unsigned char *input_mask, int image_width, int i
}
pointer_pixel = pixel;
IMP = input_mask;
IMP = input_mask;
//set the masked pixels to minimum
for (i = 0; i < image_size; i++)
@@ -670,8 +670,8 @@ void returnImage(PIXELM *pixel, float *unwrapped_image, int image_width, int im
int image_size = image_width * image_height;
float *unwrapped_image_pointer = unwrapped_image;
PIXELM *pixel_pointer = pixel;
for (i=0; i < image_size; i++)
for (i=0; i < image_size; i++)
{
*unwrapped_image_pointer = pixel_pointer->value;
pixel_pointer++;
@@ -681,8 +681,8 @@ void returnImage(PIXELM *pixel, float *unwrapped_image, int image_width, int im
//the main function of the unwrapper
void
unwrap2D(float* wrapped_image, float* UnwrappedImage, unsigned char* input_mask,
int image_width, int image_height,
unwrap2D(float* wrapped_image, float* UnwrappedImage, unsigned char* input_mask,
int image_width, int image_height,
int wrap_around_x, int wrap_around_y)
{
params_t params = {TWOPI, wrap_around_x, wrap_around_y, 0};
@@ -691,32 +691,32 @@ unwrap2D(float* wrapped_image, float* UnwrappedImage, unsigned char* input_mask,
EDGE *edge;
int image_size = image_height * image_width;
int No_of_Edges_initially = 2 * image_width * image_height;
extended_mask = (unsigned char *) calloc(image_size, sizeof(unsigned char));
pixel = (PIXELM *) calloc(image_size, sizeof(PIXELM));
edge = (EDGE *) calloc(No_of_Edges_initially, sizeof(EDGE));
extend_mask(input_mask, extended_mask, image_width, image_height, &params);
initialisePIXELs(wrapped_image, input_mask, extended_mask, pixel, image_width, image_height);
calculate_reliability(wrapped_image, pixel, image_width, image_height, &params);
horizontalEDGEs(pixel, edge, image_width, image_height, &params);
verticalEDGEs(pixel, edge, image_width, image_height, &params);
//sort the EDGEs depending on their reiability. The PIXELs with higher
//relibility (small value) first
quicker_sort(edge, edge + params.no_of_edges - 1);
//gather PIXELs into groups
gatherPIXELs(edge, &params);
unwrapImage(pixel, image_width, image_height);
maskImage(pixel, input_mask, image_width, image_height);
//copy the image from PIXELM structure to the unwrapped phase array
//passed to this function
//TODO: replace by (cython?) function to directly write into numpy array ?
returnImage(pixel, UnwrappedImage, image_width, image_height);
free(edge);
free(pixel);
free(extended_mask);