summaryrefslogtreecommitdiffstats
path: root/src
diff options
context:
space:
mode:
Diffstat (limited to 'src')
-rw-r--r--src/Core/regularisers_CPU/utils.c143
-rw-r--r--src/Core/regularisers_CPU/utils.h2
-rw-r--r--src/Matlab/mex_compile/supp_routines/Imscale2D.c50
3 files changed, 138 insertions, 57 deletions
diff --git a/src/Core/regularisers_CPU/utils.c b/src/Core/regularisers_CPU/utils.c
index 7a4e80b..5bb3a5c 100644
--- a/src/Core/regularisers_CPU/utils.c
+++ b/src/Core/regularisers_CPU/utils.c
@@ -1,21 +1,21 @@
/*
-This work is part of the Core Imaging Library developed by
-Visual Analytics and Imaging System Group of the Science Technology
-Facilities Council, STFC
-
-Copyright 2017 Daniil Kazanteev
-Copyright 2017 Srikanth Nagella, Edoardo Pasca
-
-Licensed under the Apache License, Version 2.0 (the "License");
-you may not use this file except in compliance with the License.
-You may obtain a copy of the License at
-http://www.apache.org/licenses/LICENSE-2.0
-Unless required by applicable law or agreed to in writing, software
-distributed under the License is distributed on an "AS IS" BASIS,
-WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
-See the License for the specific language governing permissions and
-limitations under the License.
-*/
+ * This work is part of the Core Imaging Library developed by
+ * Visual Analytics and Imaging System Group of the Science Technology
+ * Facilities Council, STFC
+ *
+ * Copyright 2017 Daniil Kazanteev
+ * Copyright 2017 Srikanth Nagella, Edoardo Pasca
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ * http://www.apache.org/licenses/LICENSE-2.0
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
#include "utils.h"
#include <math.h>
@@ -23,19 +23,19 @@ limitations under the License.
/* Copy Image (float) */
float copyIm(float *A, float *U, long dimX, long dimY, long dimZ)
{
- long j;
+ long j;
#pragma omp parallel for shared(A, U) private(j)
- for (j = 0; j<dimX*dimY*dimZ; j++) U[j] = A[j];
- return *U;
+ for (j = 0; j<dimX*dimY*dimZ; j++) U[j] = A[j];
+ return *U;
}
/* Copy Image */
unsigned char copyIm_unchar(unsigned char *A, unsigned char *U, int dimX, int dimY, int dimZ)
{
- int j;
+ int j;
#pragma omp parallel for shared(A, U) private(j)
- for (j = 0; j<dimX*dimY*dimZ; j++) U[j] = A[j];
- return *U;
+ for (j = 0; j<dimX*dimY*dimZ; j++) U[j] = A[j];
+ return *U;
}
/*Roll image symmetrically from top to bottom*/
@@ -59,59 +59,88 @@ float copyIm_roll(float *A, float *U, int dimX, int dimY, int roll_value, int sw
/* function that calculates TV energy
* type - 1: 2*lambda*min||\nabla u|| + ||u -u0||^2
- * type - 2: 2*lambda*min||\nabla u||
+ * type - 2: 2*lambda*min||\nabla u||
* */
float TV_energy2D(float *U, float *U0, float *E_val, float lambda, int type, int dimX, int dimY)
{
- int i, j, i1, j1, index;
- float NOMx_2, NOMy_2, E_Grad=0.0f, E_Data=0.0f;
-
- /* first calculate \grad U_xy*/
- for(j=0; j<dimY; j++) {
- for(i=0; i<dimX; i++) {
- index = j*dimX+i;
- /* boundary conditions */
- i1 = i + 1; if (i == dimX-1) i1 = i;
- j1 = j + 1; if (j == dimY-1) j1 = j;
-
- /* Forward differences */
- NOMx_2 = powf((float)(U[j1*dimX + i] - U[index]),2); /* x+ */
- NOMy_2 = powf((float)(U[j*dimX + i1] - U[index]),2); /* y+ */
- E_Grad += 2.0f*lambda*sqrtf((float)(NOMx_2) + (float)(NOMy_2)); /* gradient term energy */
- E_Data += powf((float)(U[index]-U0[index]),2); /* fidelity term energy */
- }
- }
- if (type == 1) E_val[0] = E_Grad + E_Data;
- if (type == 2) E_val[0] = E_Grad;
- return *E_val;
+ int i, j, i1, j1, index;
+ float NOMx_2, NOMy_2, E_Grad=0.0f, E_Data=0.0f;
+
+ /* first calculate \grad U_xy*/
+ for(j=0; j<dimY; j++) {
+ for(i=0; i<dimX; i++) {
+ index = j*dimX+i;
+ /* boundary conditions */
+ i1 = i + 1; if (i == dimX-1) i1 = i;
+ j1 = j + 1; if (j == dimY-1) j1 = j;
+
+ /* Forward differences */
+ NOMx_2 = powf((float)(U[j1*dimX + i] - U[index]),2); /* x+ */
+ NOMy_2 = powf((float)(U[j*dimX + i1] - U[index]),2); /* y+ */
+ E_Grad += 2.0f*lambda*sqrtf((float)(NOMx_2) + (float)(NOMy_2)); /* gradient term energy */
+ E_Data += powf((float)(U[index]-U0[index]),2); /* fidelity term energy */
+ }
+ }
+ if (type == 1) E_val[0] = E_Grad + E_Data;
+ if (type == 2) E_val[0] = E_Grad;
+ return *E_val;
}
float TV_energy3D(float *U, float *U0, float *E_val, float lambda, int type, int dimX, int dimY, int dimZ)
{
- long i, j, k, i1, j1, k1, index;
- float NOMx_2, NOMy_2, NOMz_2, E_Grad=0.0f, E_Data=0.0f;
-
- /* first calculate \grad U_xy*/
+ long i, j, k, i1, j1, k1, index;
+ float NOMx_2, NOMy_2, NOMz_2, E_Grad=0.0f, E_Data=0.0f;
+
+ /* first calculate \grad U_xy*/
for(j=0; j<(long)(dimY); j++) {
for(i=0; i<(long)(dimX); i++) {
for(k=0; k<(long)(dimZ); k++) {
- index = (dimX*dimY)*k + j*dimX+i;
+ index = (dimX*dimY)*k + j*dimX+i;
/* boundary conditions */
i1 = i + 1; if (i == (long)(dimX-1)) i1 = i;
j1 = j + 1; if (j == (long)(dimY-1)) j1 = j;
k1 = k + 1; if (k == (long)(dimZ-1)) k1 = k;
- /* Forward differences */
+ /* Forward differences */
NOMx_2 = powf((float)(U[(dimX*dimY)*k + j1*dimX+i] - U[index]),2); /* x+ */
NOMy_2 = powf((float)(U[(dimX*dimY)*k + j*dimX+i1] - U[index]),2); /* y+ */
NOMz_2 = powf((float)(U[(dimX*dimY)*k1 + j*dimX+i] - U[index]),2); /* z+ */
E_Grad += 2.0f*lambda*sqrtf((float)(NOMx_2) + (float)(NOMy_2) + (float)(NOMz_2)); /* gradient term energy */
E_Data += (powf((float)(U[index]-U0[index]),2)); /* fidelity term energy */
- }
- }
- }
- if (type == 1) E_val[0] = E_Grad + E_Data;
- if (type == 2) E_val[0] = E_Grad;
- return *E_val;
+ }
+ }
+ }
+ if (type == 1) E_val[0] = E_Grad + E_Data;
+ if (type == 2) E_val[0] = E_Grad;
+ return *E_val;
+}
+
+/* Down-Up scaling of 2D images using bilinear interpolation */
+float Im_scale2D(float *Input, float *Scaled, int w, int h, int w2, int h2)
+{
+ int x, y, index, i, j;
+ float x_ratio = ((float)(w-1))/w2;
+ float y_ratio = ((float)(h-1))/h2;
+ float A, B, C, D, x_diff, y_diff, gray;
+ #pragma omp parallel for shared (Input, Scaled) private(x, y, index, A, B, C, D, x_diff, y_diff, gray)
+ for (j=0;j<w2;j++) {
+ for (i=0;i<h2;i++) {
+ x = (int)(x_ratio * j);
+ y = (int)(y_ratio * i);
+ x_diff = (x_ratio * j) - x;
+ y_diff = (y_ratio * i) - y;
+ index = y*w+x ;
+
+ A = Input[index];
+ B = Input[index+1];
+ C = Input[index+w];
+ D = Input[index+w+1];
+
+ gray = (float)(A*(1.0 - x_diff)*(1.0 - y_diff) + B*(x_diff)*(1.0 - y_diff) +
+ C*(y_diff)*(1.0 - x_diff) + D*(x_diff*y_diff));
+
+ Scaled[i*w2+j] = gray;
+ }}
+ return *Scaled;
}
diff --git a/src/Core/regularisers_CPU/utils.h b/src/Core/regularisers_CPU/utils.h
index cfaf6d7..8f0ba59 100644
--- a/src/Core/regularisers_CPU/utils.h
+++ b/src/Core/regularisers_CPU/utils.h
@@ -29,6 +29,8 @@ CCPI_EXPORT unsigned char copyIm_unchar(unsigned char *A, unsigned char *U, int
CCPI_EXPORT float copyIm_roll(float *A, float *U, int dimX, int dimY, int roll_value, int switcher);
CCPI_EXPORT float TV_energy2D(float *U, float *U0, float *E_val, float lambda, int type, int dimX, int dimY);
CCPI_EXPORT float TV_energy3D(float *U, float *U0, float *E_val, float lambda, int type, int dimX, int dimY, int dimZ);
+CCPI_EXPORT float TV_energy3D(float *U, float *U0, float *E_val, float lambda, int type, int dimX, int dimY, int dimZ);
+CCPI_EXPORT float Im_scale2D(float *Input, float *Scaled, int w, int h, int w2, int h2);
#ifdef __cplusplus
}
#endif
diff --git a/src/Matlab/mex_compile/supp_routines/Imscale2D.c b/src/Matlab/mex_compile/supp_routines/Imscale2D.c
new file mode 100644
index 0000000..4576cce
--- /dev/null
+++ b/src/Matlab/mex_compile/supp_routines/Imscale2D.c
@@ -0,0 +1,50 @@
+/*
+ * This work is part of the Core Imaging Library developed by
+ * Visual Analytics and Imaging System Group of the Science Technology
+ * Facilities Council, STFC and Diamond Light Source Ltd.
+ *
+ * Copyright 2019 Daniil Kazantsev
+ * Copyright 2019 Srikanth Nagella, Edoardo Pasca
+ * Copyright 2019 Diamond Light Source Ltd.
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ * http://www.apache.org/licenses/LICENSE-2.0
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+
+#include "matrix.h"
+#include "mex.h"
+#include "utils.h"
+
+/**************************************************/
+void mexFunction(
+ int nlhs, mxArray *plhs[],
+ int nrhs, const mxArray *prhs[])
+{
+ int number_of_dims, X_new, Y_new;
+ mwSize dimX, dimY, dimZ;
+ const mwSize *dim_array;
+ float *A, *B;
+ mwSize dim_array2[2]; /* for scaled 2D data */
+
+ dim_array = mxGetDimensions(prhs[0]);
+ number_of_dims = mxGetNumberOfDimensions(prhs[0]);
+
+ /*Handling Matlab input data*/
+ A = (float *) mxGetData(prhs[0]); /* a 2D or 3D image/volume */
+ X_new = (int) mxGetScalar(prhs[1]); /* new size for image */
+ Y_new = (int) mxGetScalar(prhs[2]); /* new size for image */
+
+ dimX = dim_array[0]; dimY = dim_array[1]; dimZ = dim_array[2];
+ dim_array2[0] = X_new; dim_array2[1] = Y_new;
+
+ B = (float*)mxGetPr(plhs[0] = mxCreateNumericArray(2, dim_array2, mxSINGLE_CLASS, mxREAL));
+
+ Im_scale2D(A, B, dimX, dimY, X_new, Y_new);
+ }
generated by cgit v1.2.3 (git 2.25.1) at 2024-11-05 16:33:13 +0000