Commit afb1883de8958319b58e681623a35d7e1191191a

Authored by Tianshu Cheng
1 parent 91d8912e

fixed problem about sigma

stim/cuda/bsds500/Pb.cpp 0 → 100644
  1 +#include <stim/image/image.h>
  2 +//#include <cmath>
  3 +#include <stim/visualization/colormap.h>
  4 +#include <stim/image/image_contour_detection.h>
  5 +
  6 +#define SIGMA_N 3
  7 +
  8 +void array_abs(float* img, unsigned int N);
  9 +void array_multiply(float* lhs, float rhs, unsigned int N);
  10 +void array_cos(float* ptr1, float* cpu_out, unsigned int N);
  11 +void array_sin(float* ptr1, float* cpu_out, unsigned int N);
  12 +void array_atan(float* ptr1, float* cpu_out, unsigned int N);
  13 +void array_divide(float* ptr1, float* ptr2,float* cpu_quotient, unsigned int N);
  14 +void array_multiply(float* ptr1, float* ptr2, float* product, unsigned int N);
  15 +void array_add(float* ptr1, float* ptr2, float* sum, unsigned int N);
  16 +
  17 +/// This function uses odd-symmetric gaussian derivative filter to evaluate
  18 +/// the max probability of a contour on one scale, given an one-channel image
  19 +
  20 +/// @param img is an one-channel image
  21 +/// @param r is an array of radii for different scaled discs(filters)
  22 +/// @param sigma_n is the number of standard deviations used to define the sigma
  23 +
  24 +stim::image<float> Pb(stim::image<float> image, int sigma){
  25 +
  26 + unsigned int w = image.width(); // get the width of picture
  27 + unsigned int h = image.height(); // get the height of picture
  28 + unsigned N = w * h; // get the number of pixels of picture
  29 +
  30 + int r = SIGMA_N * sigma; // set the radius of filter
  31 + int winsize = 2 * r + 1; // set the winsdow size of filter
  32 +
  33 + stim::image<float> I(w, h, 1, 2); // allocate space for return image of dG1_conv2
  34 + stim::image<float> theta(w, h); // allocate space for theta matrix
  35 + stim::image<float> cos(w, h); // allocate space for cos(theta)
  36 + stim::image<float> sin(w, h); // allocate space for sin(theta)
  37 + stim::image<float> temp(w, h); // allocate space for temp
  38 + stim::image<float> Ix(w, h); // allocate space for Ix
  39 + stim::image<float> Iy(w, h); // allocate space for Iy
  40 + stim::image<float> Pb(w, h); // allocate space for Pb
  41 +
  42 + I = dG1_conv2(image, sigma); // calculate the Ix, Iy
  43 + Ix = I.channel(0);
  44 + array_abs(Ix.data(), N); //get |Ix|;
  45 + //stim::cpu2image(Ix.data(), "data_output/Pb_Ix_0924.bmp", w, h, stim::cmBrewer);
  46 + Iy = I.channel(1);
  47 + array_abs(Iy.data(), N); //get |Iy|;
  48 + //stim::cpu2image(Iy.data(), "data_output/Pb_Iy_0924.bmp", w, h, stim::cmBrewer);
  49 +
  50 + array_divide(Iy.data(), Ix.data(), temp.data(), N); //temp = Iy./Ix
  51 + array_atan(temp.data(), theta.data(), N); //theta = atan(temp)
  52 + array_cos(theta.data(), cos.data(), N); //cos = cos(theta)
  53 + array_sin(theta.data(), sin.data(), N); //sin = sin(theta)
  54 + array_multiply(Ix.data(), cos.data(), Ix.data(), N); //Ix = Ix.*cos
  55 + array_multiply(Iy.data(), sin.data(), Iy.data(), N); //Iy = Iy.*sin
  56 + array_add(Ix.data(), Iy.data(), Pb.data(), N); //Pb = Ix + Iy;
  57 +
  58 + float max = Pb.maxv(); // get the maximum of Pb used for normalization
  59 + array_multiply(Pb.data(), 1/max, N); // normalize the Pb
  60 +
  61 + //stim::cpu2image(Pb.data(), "data_output/Pb_0924.bmp", w, h, stim::cmBrewer); show the Pb(optional)
  62 +
  63 + return Pb;
  64 +
  65 +}
0 \ No newline at end of file 66 \ No newline at end of file
stim/cuda/bsds500/cPb.cpp
@@ -15,7 +15,7 @@ void array_add(float* ptr1, float* ptr2, float* sum, unsigned int N); @@ -15,7 +15,7 @@ void array_add(float* ptr1, float* ptr2, float* sum, unsigned int N);
15 /// @param alpha is is an array of weights for different scaled discs(filters) 15 /// @param alpha is is an array of weights for different scaled discs(filters)
16 /// @param s is the number of scales 16 /// @param s is the number of scales
17 17
18 -stim::image<float> cPb(stim::image<float> img, int* r, float* alpha, int s){ 18 +stim::image<float> cPb(stim::image<float> img, int* sigma, float* alpha, int s){
19 19
20 unsigned int w = img.width(); // get the width of picture 20 unsigned int w = img.width(); // get the width of picture
21 unsigned int h = img.height(); // get the height of picture 21 unsigned int h = img.height(); // get the height of picture
@@ -28,20 +28,19 @@ stim::image&lt;float&gt; cPb(stim::image&lt;float&gt; img, int* r, float* alpha, int s){ @@ -28,20 +28,19 @@ stim::image&lt;float&gt; cPb(stim::image&lt;float&gt; img, int* r, float* alpha, int s){
28 28
29 29
30 unsigned int N = w * h; // get the number of pixels 30 unsigned int N = w * h; // get the number of pixels
31 - int sigma_n = 3; // set the number of standard deviations used to define the sigma  
32 31
33 - std::ostringstream ss; // (optional) set the stream to designate the test result file 32 + //std::ostringstream ss; // (optional) set the stream to designate the test result file
34 33
35 stim::image<float> temp; // set the temporary image to store the addtion result 34 stim::image<float> temp; // set the temporary image to store the addtion result
36 35
37 for (int i = 0; i < c; i++){ 36 for (int i = 0; i < c; i++){
38 for (int j = 0; j < s; j++){ 37 for (int j = 0; j < s; j++){
39 38
40 - ss << "data_output/cPb_slice"<< i*s + j << ".bmp"; // set the name for test result file (optional)  
41 - std::string sss = ss.str(); 39 + //ss << "data_output/cPb_slice"<< i*s + j << ".bmp"; // set the name for test result file (optional)
  40 + //std::string sss = ss.str();
42 41
43 // get the gaussian gradient by convolving each image slice with the mask 42 // get the gaussian gradient by convolving each image slice with the mask
44 - temp = Pb(img.channel(i), r[i*s + j], sigma_n); 43 + temp = Pb(img.channel(i), sigma[i*s + j]);
45 44
46 // output the test result of each slice (optional) 45 // output the test result of each slice (optional)
47 //stim::cpu2image(temp.data(), sss, w, h, stim::cmBrewer); 46 //stim::cpu2image(temp.data(), sss, w, h, stim::cmBrewer);
@@ -52,7 +51,7 @@ stim::image&lt;float&gt; cPb(stim::image&lt;float&gt; img, int* r, float* alpha, int s){ @@ -52,7 +51,7 @@ stim::image&lt;float&gt; cPb(stim::image&lt;float&gt; img, int* r, float* alpha, int s){
52 // add up all the weighted gaussian gradients 51 // add up all the weighted gaussian gradients
53 array_add(cPb.data(), temp.data(), cPb.data(), N); 52 array_add(cPb.data(), temp.data(), cPb.data(), N);
54 53
55 - ss.str(""); //(optional) clear the space for stream 54 + //ss.str(""); //(optional) clear the space for stream
56 55
57 } 56 }
58 } 57 }
stim/cuda/bsds500/dG1_conv2.cpp
@@ -2,25 +2,26 @@ @@ -2,25 +2,26 @@
2 //#include <cmath> 2 //#include <cmath>
3 #include <stim/visualization/colormap.h> 3 #include <stim/visualization/colormap.h>
4 #include <stim/image/image_contour_detection.h> 4 #include <stim/image/image_contour_detection.h>
  5 +#define SIGMA_N 3
5 6
6 /// This function generates the first-order gaussian derivative filter gx gy, 7 /// This function generates the first-order gaussian derivative filter gx gy,
7 /// convolves the image with gx gy, 8 /// convolves the image with gx gy,
8 /// and returns an image class which channel(0) is Ix and channel(1) is Iy 9 /// and returns an image class which channel(0) is Ix and channel(1) is Iy
9 10
10 /// @param img is the one-channel image 11 /// @param img is the one-channel image
11 -/// @param r is an array of radii for different scaled discs(filters)  
12 -/// @param sigma_n is the number of standard deviations used to define the sigma 12 +/// @param sigma is the parameter for gaussian function
13 13
14 void conv2_sep(float* img, unsigned int x, unsigned int y, float* kernel0, unsigned int k0, float* kernel1, unsigned int k1); 14 void conv2_sep(float* img, unsigned int x, unsigned int y, float* kernel0, unsigned int k0, float* kernel1, unsigned int k1);
15 //void array_abs(float* img, unsigned int N); 15 //void array_abs(float* img, unsigned int N);
16 16
17 -stim::image<float> Gd1(stim::image<float> image, int r, unsigned int sigma_n){ 17 +stim::image<float> dG1_conv2(stim::image<float> image, int sigma){
18 18
19 unsigned int w = image.width(); // get the width of picture 19 unsigned int w = image.width(); // get the width of picture
20 unsigned int h = image.height(); // get the height of picture 20 unsigned int h = image.height(); // get the height of picture
21 unsigned N = w * h; // get the number of pixels of picture 21 unsigned N = w * h; // get the number of pixels of picture
22 - int winsize = 2 * r + 1; // set the winsdow size of disc(filter)  
23 - float sigma = float(r)/float(sigma_n); // calculate the sigma used in gaussian function 22 +
  23 + int r = SIGMA_N * sigma;
  24 + int winsize = 2 * SIGMA_N * sigma + 1; // set the winsdow size of filter
24 25
25 stim::image<float> I(w, h, 1, 2); // allocate space for return image class 26 stim::image<float> I(w, h, 1, 2); // allocate space for return image class
26 stim::image<float> Ix(w, h); // allocate space for Ix 27 stim::image<float> Ix(w, h); // allocate space for Ix
stim/cuda/bsds500/dG1_theta_conv2.cpp
@@ -4,6 +4,7 @@ @@ -4,6 +4,7 @@
4 #include <stim/image/image_contour_detection.h> 4 #include <stim/image/image_contour_detection.h>
5 5
6 #define PI 3.1415926 6 #define PI 3.1415926
  7 +#define SIGMA_N 3
7 8
8 void array_multiply(float* lhs, float rhs, unsigned int N); 9 void array_multiply(float* lhs, float rhs, unsigned int N);
9 void array_add(float* ptr1, float* ptr2, float* sum, unsigned int N); 10 void array_add(float* ptr1, float* ptr2, float* sum, unsigned int N);
@@ -16,33 +17,36 @@ void array_abs(float* img, unsigned int N); @@ -16,33 +17,36 @@ void array_abs(float* img, unsigned int N);
16 /// @param sigma_n is the number of standard deviations used to define the sigma 17 /// @param sigma_n is the number of standard deviations used to define the sigma
17 /// @param theta is angle used for computing the gradient 18 /// @param theta is angle used for computing the gradient
18 19
19 -stim::image<float> Gd_odd(stim::image<float> image, int r, unsigned int sigma_n, float theta){ 20 +stim::image<float> dG1_theta_conv2(stim::image<float> image, int sigma, float theta){
20 21
21 float theta_r = (theta * PI)/180; //change angle unit from degree to rad 22 float theta_r = (theta * PI)/180; //change angle unit from degree to rad
22 23
23 unsigned int w = image.width(); // get the width of picture 24 unsigned int w = image.width(); // get the width of picture
24 unsigned int h = image.height(); // get the height of picture 25 unsigned int h = image.height(); // get the height of picture
  26 +
  27 + int r = SIGMA_N * sigma;
25 unsigned N = w * h; // get the number of pixels of picture 28 unsigned N = w * h; // get the number of pixels of picture
26 - int winsize = 2 * r + 1; // set the winsdow size of disc(filter)  
27 29
28 - stim::image<float> I(w, h, 1, 2); // allocate space for return image of Gd1  
29 - stim::image<float> Ix(w, h); // allocate space for Ix  
30 - stim::image<float> Iy(w, h); // allocate space for Iy  
31 - stim::image<float> Gd_odd_theta(w, h); // allocate space for Pb 30 + int winsize = 2 * SIGMA_N * sigma + 1; // set the winsdow size of filter
  31 +
  32 + stim::image<float> I(w, h, 1, 2); // allocate space for return image of dG1_conv2
  33 + stim::image<float> Ix(w, h); // allocate space for Ix
  34 + stim::image<float> Iy(w, h); // allocate space for Iy
  35 + stim::image<float> dG1_theta(w, h); // allocate space for Pb
32 36
33 - I = Gd1(image, r, sigma_n); // calculate the Ix, Iy 37 + I = dG1_conv2(image, sigma); // calculate the Ix, Iy
34 Ix = I.channel(0); 38 Ix = I.channel(0);
35 Iy = I.channel(1); 39 Iy = I.channel(1);
36 40
37 array_multiply(Ix.data(), cos(theta_r), N); //Ix = Ix*cos(theta_r) 41 array_multiply(Ix.data(), cos(theta_r), N); //Ix = Ix*cos(theta_r)
38 array_multiply(Iy.data(), sin(theta_r), N); //Iy = Iy*sin(theta_r) 42 array_multiply(Iy.data(), sin(theta_r), N); //Iy = Iy*sin(theta_r)
39 - array_add(Ix.data(), Iy.data(), Gd_odd_theta.data(), N); //Gd_odd_theta = Ix + Iy;  
40 - array_abs(Gd_odd_theta.data(), N); 43 + array_add(Ix.data(), Iy.data(), dG1_theta.data(), N); //dG1_theta = Ix + Iy;
  44 + array_abs(dG1_theta.data(), N);
41 45
42 - //stim::cpu2image(I.channel(0).data(), "data_output/Gd_odd_x_0919.bmp", w, h, stim::cmBrewer);  
43 - //stim::cpu2image(I.channel(1).data(), "data_output/Gd_odd_y_0919.bmp", w, h, stim::cmBrewer);  
44 - //stim::cpu2image(Gd_odd_theta.data(), "data_output/Gd_odd_theta_0919.bmp", w, h, stim::cmBrewer); 46 + //stim::cpu2image(I.channel(0).data(), "data_output/dG1_theta_x_0919.bmp", w, h, stim::cmBrewer);
  47 + //stim::cpu2image(I.channel(1).data(), "data_output/dG1_theta_y_0919.bmp", w, h, stim::cmBrewer);
  48 + //stim::cpu2image(dG1_theta.data(), "data_output/dG1_theta_0919.bmp", w, h, stim::cmBrewer);
45 49
46 - return Gd_odd_theta; 50 + return dG1_theta;
47 51
48 } 52 }
stim/cuda/bsds500/dG2_conv2.cpp
@@ -2,25 +2,26 @@ @@ -2,25 +2,26 @@
2 //#include <cmath> 2 //#include <cmath>
3 #include <stim/visualization/colormap.h> 3 #include <stim/visualization/colormap.h>
4 #include <stim/image/image_contour_detection.h> 4 #include <stim/image/image_contour_detection.h>
  5 +#define SIGMA_N 3
5 6
6 /// This function generates the second-order gaussian derivative filter gxx gyy, 7 /// This function generates the second-order gaussian derivative filter gxx gyy,
7 /// convolves the image with gxx gyy, 8 /// convolves the image with gxx gyy,
8 /// and returns an image class which channel(0) is Ixx and channel(1) is Iyy 9 /// and returns an image class which channel(0) is Ixx and channel(1) is Iyy
9 10
10 /// @param img is the one-channel image 11 /// @param img is the one-channel image
11 -/// @param r is an array of radii for different scaled discs(filters)  
12 -/// @param sigma_n is the number of standard deviations used to define the sigma 12 +/// @param sigma is the parameter for gaussian function
13 13
14 void conv2_sep(float* img, unsigned int x, unsigned int y, float* kernel0, unsigned int k0, float* kernel1, unsigned int k1); 14 void conv2_sep(float* img, unsigned int x, unsigned int y, float* kernel0, unsigned int k0, float* kernel1, unsigned int k1);
15 //void array_abs(float* img, unsigned int N); 15 //void array_abs(float* img, unsigned int N);
16 16
17 -stim::image<float> Gd2(stim::image<float> image, int r, unsigned int sigma_n){ 17 +stim::image<float> dG2_conv2(stim::image<float> image, int sigma){
18 18
19 unsigned int w = image.width(); // get the width of picture 19 unsigned int w = image.width(); // get the width of picture
20 unsigned int h = image.height(); // get the height of picture 20 unsigned int h = image.height(); // get the height of picture
21 unsigned N = w * h; // get the number of pixels of picture 21 unsigned N = w * h; // get the number of pixels of picture
22 - int winsize = 2 * r + 1; // set the winsdow size of disc(filter)  
23 - float sigma = float(r)/float(sigma_n); // calculate the sigma used in gaussian function 22 +
  23 + int winsize = 2 * SIGMA_N * sigma + 1; // set the winsdow size of filter
  24 + int r = SIGMA_N * sigma;
24 25
25 stim::image<float> I(w, h, 1, 2); // allocate space for return image class 26 stim::image<float> I(w, h, 1, 2); // allocate space for return image class
26 stim::image<float> Ixx(w, h); // allocate space for Ixx 27 stim::image<float> Ixx(w, h); // allocate space for Ixx
stim/cuda/bsds500/dG2_d2x_theta_conv2.cpp 0 → 100644
  1 +#include <stim/image/image.h>
  2 +#include <cmath>
  3 +#include <stim/visualization/colormap.h>
  4 +#include <stim/image/image_contour_detection.h>
  5 +
  6 +#define SIGMA_N 3
  7 +
  8 +/// This function evaluates the theta-dependent even-symmetric gaussian derivative gradient of an one-channel image
  9 +
  10 +/// @param img is the one-channel image
  11 +/// @param r is an array of radii for different scaled discs(filters)
  12 +/// @param sigma_n is the number of standard deviations used to define the sigma
  13 +/// @param theta is angle used for computing the gradient
  14 +
  15 +void conv2(float* img, float* mask, float* cpu_copy, unsigned int w, unsigned int h, unsigned int M);
  16 +void array_abs(float* img, unsigned int N);
  17 +
  18 +stim::image<float> dG2_d2x_theta_conv2(stim::image<float> image, int sigma, float theta){
  19 +
  20 + unsigned int w = image.width(); // get the width of picture
  21 + unsigned int h = image.height(); // get the height of picture
  22 + unsigned N = w * h; // get the number of pixels of picture
  23 +
  24 + int r = SIGMA_N * sigma; // set the radius of filter
  25 + int winsize = 2 * SIGMA_N * sigma + 1; // set the winsdow size of filter
  26 +
  27 + stim::image<float> I(w, h, 1, 2); // allocate space for return image class
  28 + stim::image<float> dG2_d2x_theta(w, h); // allocate space for dG2_d2x_theta
  29 + stim::image<float> mask_x(winsize, winsize); // allocate space for x-axis-oriented filter
  30 + stim::image<float> mask_r(winsize, winsize); // allocate space for theta-oriented filter
  31 +
  32 + for (int j = 0; j < winsize; j++){
  33 + for (int i = 0; i< winsize; i++){
  34 +
  35 + int x = i - r; //range of x
  36 + int y = j - r; //range of y
  37 +
  38 + // create the x-oriented gaussian derivative filter mask_x
  39 + mask_x.data()[j*winsize + i] = (-1) * (1 - pow(x, 2)) * exp((-1)*(pow(x, 2))/(2*pow(sigma, 2))) * exp((-1)*(pow(y, 2))/(2*pow(sigma, 2)));
  40 +
  41 + }
  42 + }
  43 +
  44 + mask_r = mask_x.rotate(theta, r, r);
  45 + //mask_r = mask_x.rotate(45, r, r);
  46 + //stim::cpu2image(mask_r.data(), "data_output/mask_r_0919.bmp", winsize, winsize, stim::cmBrewer);
  47 +
  48 + // do the 2D convolution with image and mask
  49 + conv2(image.data(), mask_r.data(), dG2_d2x_theta.data(), w, h, winsize);
  50 + array_abs(dG2_d2x_theta.data(), N);
  51 +
  52 + //stim::cpu2image(dG2_d2x_theta.data(), "data_output/dG2_d2x_theta_0919.bmp", w, h, stim::cmGrayscale);
  53 +
  54 + return dG2_d2x_theta;
  55 +}
0 \ No newline at end of file 56 \ No newline at end of file
stim/cuda/bsds500/kmeans.cpp
@@ -8,7 +8,7 @@ @@ -8,7 +8,7 @@
8 /// This function use cvkmeans to cluster given textons 8 /// This function use cvkmeans to cluster given textons
9 9
10 /// @param testons is a multi-channel image 10 /// @param testons is a multi-channel image
11 -/// @param k is the number of clusters 11 +/// @param K is the number of clusters
12 12
13 stim::image<float> kmeans(stim::image<float> textons, unsigned int K){ 13 stim::image<float> kmeans(stim::image<float> textons, unsigned int K){
14 14
stim/cuda/bsds500/laplacian_conv2.cpp
@@ -4,6 +4,7 @@ @@ -4,6 +4,7 @@
4 #include <stim/image/image_contour_detection.h> 4 #include <stim/image/image_contour_detection.h>
5 5
6 #define PI 3.1415926 6 #define PI 3.1415926
  7 +#define SIGMA_N 3
7 8
8 void array_multiply(float* lhs, float rhs, unsigned int N); 9 void array_multiply(float* lhs, float rhs, unsigned int N);
9 void array_add(float* ptr1, float* ptr2, float* sum, unsigned int N); 10 void array_add(float* ptr1, float* ptr2, float* sum, unsigned int N);
@@ -15,27 +16,28 @@ void array_abs(float* img, unsigned int N); @@ -15,27 +16,28 @@ void array_abs(float* img, unsigned int N);
15 /// @param r is an array of radii for different scaled discs(filters) 16 /// @param r is an array of radii for different scaled discs(filters)
16 /// @param sigma_n is the number of standard deviations used to define the sigma 17 /// @param sigma_n is the number of standard deviations used to define the sigma
17 18
18 -stim::image<float> Gd_center(stim::image<float> image, int r, unsigned int sigma_n){ 19 +stim::image<float> laplacian_conv2(stim::image<float> image, int sigma){
19 20
20 unsigned int w = image.width(); // get the width of picture 21 unsigned int w = image.width(); // get the width of picture
21 unsigned int h = image.height(); // get the height of picture 22 unsigned int h = image.height(); // get the height of picture
22 unsigned N = w * h; // get the number of pixels of picture 23 unsigned N = w * h; // get the number of pixels of picture
23 - int winsize = 2 * r + 1; // set the winsdow size of disc(filter) 24 +
  25 + int winsize = 2 * SIGMA_N * sigma + 1; // set the winsdow size of filter
24 26
25 - stim::image<float> I(w, h, 1, 2); // allocate space for return image of Gd2 27 + stim::image<float> I(w, h, 1, 2); // allocate space for return image of dG2_conv2
26 stim::image<float> Ixx(w, h); // allocate space for Ixx 28 stim::image<float> Ixx(w, h); // allocate space for Ixx
27 stim::image<float> Iyy(w, h); // allocate space for Iyy 29 stim::image<float> Iyy(w, h); // allocate space for Iyy
28 - stim::image<float> Gd_center(w, h); // allocate space for Pb 30 + stim::image<float> laplacian(w, h); // allocate space for Pb
29 31
30 - I = Gd2(image, r, sigma_n); // calculate the Ixx, Iyy 32 + I = dG2_conv2(image, sigma); // calculate the Ixx, Iyy
31 Ixx = I.channel(0); 33 Ixx = I.channel(0);
32 Iyy = I.channel(1); 34 Iyy = I.channel(1);
33 35
34 - array_add(Ixx.data(), Iyy.data(), Gd_center.data(), N); //Gd_center = Ixx + Iyy;  
35 - array_abs(Gd_center.data(), N); 36 + array_add(Ixx.data(), Iyy.data(), laplacian.data(), N); //laplacian = Ixx + Iyy;
  37 + array_abs(laplacian.data(), N);
36 38
37 - //stim::cpu2image(Gd_center.data(), "data_output/Gd_center_0919.bmp", w, h, stim::cmBrewer); 39 + //stim::cpu2image(laplacian.data(), "data_output/laplacian_0919.bmp", w, h, stim::cmBrewer);
38 40
39 - return Gd_center; 41 + return laplacian;
40 42
41 } 43 }
stim/cuda/bsds500/main.cpp 0 → 100644
  1 +#include <stim/image/image.h>
  2 +#include <cmath>
  3 +#include <stim/visualization/colormap.h>
  4 +#include <stim/image/image_contour_detection.h>
  5 +#include <iostream>
  6 +#include <stim/visualization/colormap.h>
  7 +#include <sstream>
  8 +/// calculate the cPb, tPb and mPb given a multi-channel image
  9 +
  10 +void array_multiply(float* lhs, float rhs, unsigned int N);
  11 +void array_add(float* ptr1, float* ptr2, float* sum, unsigned int N);
  12 +
  13 +int main()
  14 +{
  15 + stim::image<float> image; // generate an image object
  16 +
  17 + //std::ostringstream ss1,ss2; // (optional) set the stream to designate the test result file
  18 +
  19 + //for(unsigned int i = 27; i < 31; i++){
  20 +
  21 + //ss1 << "3d_sample/b0"<< i << ".bmp"; // (optional) set the name for test result file
  22 + //std::string sss1 = ss1.str(); // (optional)
  23 +
  24 + //image.load("bone1001_3.bmp");
  25 + //image.load(sss1);
  26 + //image.load("101085.bmp"); // load the input image
  27 + image.load("slice00_500_500.bmp"); // load the input image
  28 + image = image.channel(0); // get the first channel of black-and-white image
  29 +
  30 +
  31 + unsigned int w = image.width(); // get the width of picture
  32 + unsigned int h = image.height(); // get the height of picture
  33 + unsigned int N = w * h; // get the number of pixels
  34 + int c = image.channels(); // get the number if channels of picture
  35 + int s = 3; // set the number of scales
  36 + int sigma[3] = {1,2,3}; // set an array of sigmas for different scaled gaussian filters for cPb, the length of array is c*s
  37 + int r[3] = {5,10,20}; // set an array of radii for different scaled discs(filters) for tPb, the length of array is c*s
  38 + float alpha[3] = {1,1,1}; // set an array of weights for different scaled discs(filters)
  39 + unsigned int theta_n = 8; // set the number of angles used in filter orientation
  40 + unsigned int bin_n = 16; // set the number of bins used in chi-distance
  41 + unsigned int K = 16; // set the number of cludters, K should be multiple of bin_n
  42 +
  43 + stim::image<float> img_cPb(w, h); // allocate the space for cPb
  44 + stim::image<float> img_tPb(w, h); // allocate the space for tPb
  45 + stim::image<float> img_mPb(w, h); // allocate the space for tPb
  46 +
  47 + std::cout<<"imagesize: "<< w <<"*"<< h <<'\n';
  48 +
  49 + //*******************cPb********************
  50 + std::clock_t start1; // (optional) set the timer to calculate the total time
  51 + start1 = std::clock(); // (optional) set timer start point
  52 +
  53 + std::cout<<"begin cPb"<<'\n';
  54 +
  55 + img_cPb = cPb(image, sigma, alpha, s);
  56 +
  57 + // show the cPb (optional)
  58 + stim::cpu2image(img_cPb.data(), "data_output/img_cPb_1008.bmp", w, h, stim::cmBrewer);
  59 + //ss2 << "3d_sample/cPb0"<< i << ".bmp"; // (optional) set the name for test result file
  60 + //std::string sss2 = ss2.str();
  61 +
  62 + //stim::cpu2image(img_cPb.data(), "data_output/bone_cPb3.1.bmp", w, h, stim::cmBrewer);
  63 +
  64 + double duration1 = ( std::clock() - start1 ) / (double) CLOCKS_PER_SEC; // (optional) calculate the total time
  65 + std::cout<<"cPb time: "<< duration1 <<"s"<<'\n'; // (optional) show the total time
  66 +
  67 + //ss1.str(""); //(optional) clear the space for stream
  68 + //ss2.str(""); //(optional) clear the space for stream
  69 +
  70 +
  71 + //*******************tPb********************
  72 +
  73 + std::cout<<"begin tPb"<<'\n';
  74 +
  75 +
  76 + std::clock_t start2; // (optional) set the timer to calculate the total time
  77 + start2 = std::clock(); // (optional) set timer start point
  78 +
  79 + img_tPb = tPb(image, r, alpha, theta_n, bin_n, s, K);
  80 +
  81 + // show the tPb (optional)
  82 + // stim::cpu2image(img_tPb.data(), "data_output/bone_tPb_3.1.bmp", w, h, stim::cmBrewer);
  83 + stim::cpu2image(img_tPb.data(), "data_output/img_tPb_1008.bmp", w, h, stim::cmBrewer);
  84 +
  85 + double duration2 = ( std::clock() - start2 ) / (double) CLOCKS_PER_SEC; // (optional) calculate the total time
  86 + std::cout<<"tPb time: "<< duration2 <<"s"<<'\n'; // (optional) show the total time
  87 +
  88 +
  89 + //*******************************************
  90 +
  91 + double duration3 = ( std::clock() - start1 ) / (double) CLOCKS_PER_SEC; // (optional) calculate the total time
  92 + std::cout<<"total running time: "<< duration3 <<"s"<<'\n'; // (optional) show the total time
  93 + //}
  94 + //******************mPb**********************
  95 + // set parameters for linear combination
  96 +
  97 + float a = 1;
  98 + float b = 0.5;
  99 +
  100 + // create mPb by linearly combined cPb and tPb
  101 + array_multiply(img_cPb.data(), a, N);
  102 + array_multiply(img_tPb.data(), b, N);
  103 + array_add(img_cPb.data(), img_tPb.data(), img_mPb.data(), N);
  104 +
  105 + //ss2 << "3d_sample/mPb0"<< i << ".bmp"; // (optional) set the name for test result file
  106 + //std::string sss2 = ss2.str();
  107 + // show the mPb (optional)
  108 + stim::cpu2image(img_mPb.data(), "data_output/img_mPb_1008.bmp", w, h, stim::cmBrewer);
  109 +
  110 + //stim::cpu2image(img_mPb.data(), sss2, w, h, stim::cmBrewer);
  111 +
  112 + //ss1.str(""); //(optional) clear the space for stream
  113 + //ss2.str(""); //(optional) clear the space for stream
  114 +
  115 + //}
  116 + return 0;
  117 +
  118 +}
stim/cuda/bsds500/tPb.cpp
@@ -39,7 +39,7 @@ stim::image&lt;float&gt; tPb(stim::image&lt;float&gt; img, int* r, float* alpha, unsigned in @@ -39,7 +39,7 @@ stim::image&lt;float&gt; tPb(stim::image&lt;float&gt; img, int* r, float* alpha, unsigned in
39 39
40 img_texture = kmeans(img_textons, K); // changing kmeans result into float type is required 40 img_texture = kmeans(img_textons, K); // changing kmeans result into float type is required
41 41
42 - stim::cpu2image(img_texture.data(), "data_output/texture_0925.bmp", w, h, stim::cmBrewer); 42 + stim::cpu2image(img_texture.data(), "data_output/texture.bmp", w, h, stim::cmBrewer);
43 43
44 44
45 unsigned int max1 = img_texture.maxv(); // get the maximum of Pb used for normalization 45 unsigned int max1 = img_texture.maxv(); // get the maximum of Pb used for normalization
stim/cuda/bsds500/textons.cpp
@@ -19,12 +19,7 @@ stim::image&lt;float&gt; textons(stim::image&lt;float&gt; image, unsigned int theta_n){ @@ -19,12 +19,7 @@ stim::image&lt;float&gt; textons(stim::image&lt;float&gt; image, unsigned int theta_n){
19 stim::image<float> textons(w, h, 1, theta_n*2+1); // allocate space for textons 19 stim::image<float> textons(w, h, 1, theta_n*2+1); // allocate space for textons
20 stim::image<float> temp(w, h); // allocate space for temp 20 stim::image<float> temp(w, h); // allocate space for temp
21 21
22 - unsigned int r_odd = 3; // set disc radii for odd-symmetric filter  
23 - unsigned int sigma_n_odd = 3; // set sigma_n for odd-symmetric filter  
24 - unsigned int r_even = 3; // set disc radii for even-symmetric filter  
25 - unsigned int sigma_n_even = 3; // set sigma_n for even-symmetric filter  
26 - unsigned int r_center = 3; // set disc radii for center-surround filter  
27 - unsigned int sigma_n_center = 3; // set sigma_n for center-surround filter 22 + int sigma = 1; // set sigma for odd-symmetric, even-symmetric and center-surround filter filter
28 23
29 //std::ostringstream ss1, ss2; // (optional) set the stream to designate the test result file 24 //std::ostringstream ss1, ss2; // (optional) set the stream to designate the test result file
30 25
@@ -37,11 +32,11 @@ stim::image&lt;float&gt; textons(stim::image&lt;float&gt; image, unsigned int theta_n){ @@ -37,11 +32,11 @@ stim::image&lt;float&gt; textons(stim::image&lt;float&gt; image, unsigned int theta_n){
37 32
38 float theta = 180 * ((float)i/theta_n); // calculate the even-splited angle for each oriented filter 33 float theta = 180 * ((float)i/theta_n); // calculate the even-splited angle for each oriented filter
39 34
40 - temp = Gd_odd(image, r_odd, sigma_n_odd, theta); // return Gd_odd to temp 35 + temp = dG1_theta_conv2(image, sigma, theta); // return dG1_theta_conv2 to temp
41 //stim::cpu2image(temp.data(), sss1, w, h, stim::cmBrewer); 36 //stim::cpu2image(temp.data(), sss1, w, h, stim::cmBrewer);
42 textons.set_channel(i, temp.data()); // copy temp to ith channel of textons 37 textons.set_channel(i, temp.data()); // copy temp to ith channel of textons
43 38
44 - temp = Gd_even(image, r_even, sigma_n_even, theta); // return Gd_even to temp 39 + temp = dG2_d2x_theta_conv2(image, sigma, theta); // return dG2_d2x_theta_conv2 to temp
45 //stim::cpu2image(temp.data(), sss2, w, h, stim::cmBrewer); 40 //stim::cpu2image(temp.data(), sss2, w, h, stim::cmBrewer);
46 textons.set_channel(i + theta_n, temp.data()); // copy temp to (i+theta_n)th channel of textons 41 textons.set_channel(i + theta_n, temp.data()); // copy temp to (i+theta_n)th channel of textons
47 42
@@ -50,7 +45,7 @@ stim::image&lt;float&gt; textons(stim::image&lt;float&gt; image, unsigned int theta_n){ @@ -50,7 +45,7 @@ stim::image&lt;float&gt; textons(stim::image&lt;float&gt; image, unsigned int theta_n){
50 45
51 } 46 }
52 47
53 - temp = Gd_center(image, r_center, sigma_n_center); // return Gd_center to temp 48 + temp = laplacian_conv2(image, sigma); // return laplacian_conv2 to temp
54 //stim::cpu2image(temp.data(), "data_output/textons_channel_16.bmp", w, h, stim::cmBrewer); 49 //stim::cpu2image(temp.data(), "data_output/textons_channel_16.bmp", w, h, stim::cmBrewer);
55 textons.set_channel(theta_n*2, temp.data()); // copy temp to (theta_n*2)th channel of textons 50 textons.set_channel(theta_n*2, temp.data()); // copy temp to (theta_n*2)th channel of textons
56 51
stim/image/image_contour_detection.h
@@ -3,14 +3,14 @@ @@ -3,14 +3,14 @@
3 //stim::image<float> func_mPb_theta(stim::image<float> img, float theta, int* r, float* alpha, int s); 3 //stim::image<float> func_mPb_theta(stim::image<float> img, float theta, int* r, float* alpha, int s);
4 //stim::image<float> func_mPb(stim::image<float> img, unsigned int theta_n, int* r, float* alpha, int s); 4 //stim::image<float> func_mPb(stim::image<float> img, unsigned int theta_n, int* r, float* alpha, int s);
5 5
6 -stim::image<float> Gd1(stim::image<float> image, int r, unsigned int sigma_n);  
7 -stim::image<float> Gd2(stim::image<float> image, int r, unsigned int sigma_n);  
8 -stim::image<float> Gd_odd(stim::image<float> image, int r, unsigned int sigma_n, float theta);  
9 -stim::image<float> Gd_even(stim::image<float> image, int r, unsigned int sigma_n, float theta);  
10 -stim::image<float> Gd_center(stim::image<float> image, int r, unsigned int sigma_n); 6 +stim::image<float> dG1_conv2(stim::image<float> image, int sigma);
  7 +stim::image<float> dG2_conv2(stim::image<float> image, int sigma);
  8 +stim::image<float> dG1_theta_conv2(stim::image<float> image, int sigma, float theta);
  9 +stim::image<float> dG2_d2x_theta_conv2(stim::image<float> image, int sigma, float theta);
  10 +stim::image<float> laplacian_conv2(stim::image<float> image, int sigma);
11 11
12 stim::image<float> textons(stim::image<float> image, unsigned int theta_n); 12 stim::image<float> textons(stim::image<float> image, unsigned int theta_n);
13 stim::image<float> kmeans(stim::image<float> textons, unsigned int K); 13 stim::image<float> kmeans(stim::image<float> textons, unsigned int K);
14 -stim::image<float> Pb(stim::image<float> image, int r, unsigned int sigma_n);  
15 -stim::image<float> cPb(stim::image<float> img, int* r, float* alpha, int s); 14 +stim::image<float> Pb(stim::image<float> image, int sigma);
  15 +stim::image<float> cPb(stim::image<float> img, int* sigma, float* alpha, int s);
16 stim::image<float> tPb(stim::image<float> img, int* r, float* alpha, unsigned int theta_n, unsigned int bin_n, int s, unsigned int K); 16 stim::image<float> tPb(stim::image<float> img, int* r, float* alpha, unsigned int theta_n, unsigned int bin_n, int s, unsigned int K);
17 \ No newline at end of file 17 \ No newline at end of file