Blame view

stim/cuda/bsds500/dG2_conv2.cuh 3.74 KB
018b00d6   David Mayerich   adapted bsds500 f...
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
  #ifndef STIM_CUDA_DG2_CONV2_CUH
  #define STIM_CUDA_DG2_CONV2_CUH
  
  #include <stim/image/image.h>
  //#include <cmath>
  #include <stim/visualization/colormap.h>
  //#include <stim/image/image_contour_detection.h>
  #define SIGMA_N 3
  
  /// This function generates the second-order gaussian derivative filter gxx gyy, 
  /// convolves the image with gxx gyy, 
  /// and returns an image class which channel(0) is Ixx and channel(1) is Iyy 
  
  /// @param img is the one-channel image
  /// @param sigma is the parameter for gaussian function
  
  //void conv2_sep(float* img, unsigned int x, unsigned int y, float* kernel0, unsigned int k0, float* kernel1, unsigned int k1);
  //void array_abs(float* img, unsigned int N);
  
  stim::image<float> dG2_conv2(stim::image<float> image, int sigma){
  
  	unsigned int w = image.width();    // get the width of picture
  	unsigned int h = image.height();   // get the height of picture
  	unsigned N = w * h;				   // get the number of pixels of picture
  	
  	int winsize = 2 * SIGMA_N * sigma + 1;           // set the winsdow size of filter
  	int r = SIGMA_N * sigma;  
  
  	stim::image<float> I(w, h, 1, 2);      // allocate space for return image class
  	stim::image<float> Ixx(w, h);      // allocate space for Ixx
  	stim::image<float> Iyy(w, h);      // allocate space for Iyy
  	Ixx = image;  // initialize Ixx
  	Iyy = image;  // initialize Iyy
  
  	float* array_x1;   
  	array_x1 = new float[winsize];  //allocate space for the 1D x-oriented gaussian derivative filter array_x1 for gxx
  	float* array_y1;   
  	array_y1 = new float[winsize];  //allocate space for the 1D y-oriented gaussian derivative filter array_y1 for gxx
  	float* array_x2;   
  	array_x2 = new float[winsize];  //allocate space for the 1D x-oriented gaussian derivative filter array_x2 for gyy
  	float* array_y2;   
  	array_y2 = new float[winsize];  //allocate space for the 1D y-oriented gaussian derivative filter array_y2 for gyy
  
  
  	for (int i = 0; i < winsize; i++){	
  
  		int x = i - r;          //range of x
  		int y = i - r;          //range of y
  
  		// create the 1D x-oriented gaussian derivative filter array_x1 for gxx
  		array_x1[i] = (-1) * (1 - pow(x, 2)) * exp((-1)*(pow(x, 2))/(2*pow(sigma, 2)));
  		// create the 1D y-oriented gaussian derivative filter array_y1  for gxx
  		array_y1[i] = exp((-1)*(pow(y, 2))/(2*pow(sigma, 2)));
  		// create the 1D x-oriented gaussian derivative filter array_x2 for gyy
  		array_x2[i] = exp((-1)*(pow(x, 2))/(2*pow(sigma, 2)));
  		// create the 1D y-oriented gaussian derivative filter array_y2  for gyy
  		array_y2[i] = (-1) * (1 - pow(y, 2)) * exp((-1)*(pow(y, 2))/(2*pow(sigma, 2)));
  	}
  
  	//stim::cpu2image(array_x1, "data_output/array_x1_0915.bmp", winsize, 1, stim::cmBrewer); // (optional) show the mask result
  	//stim::cpu2image(array_y1, "data_output/array_y1_0915.bmp", winsize, 1, stim::cmBrewer); // (optional) show the mask result
  	//stim::cpu2image(array_x2, "data_output/array_x2_0915.bmp", winsize, 1, stim::cmBrewer); // (optional) show the mask result
  	//stim::cpu2image(array_y2, "data_output/array_y2_0915.bmp", winsize, 1, stim::cmBrewer); // (optional) show the mask result
  
  	// get Ixx by convolving the image with gxx
  	conv2_sep(Ixx.data(), w, h, array_x1, winsize, array_y1, winsize);
  	
  	//stim::cpu2image(Ixx.data(), "data_output/Ixx_0915.bmp", w, h, stim::cmBrewer); 
  	// get Iyy by convolving the image with gyy
  	conv2_sep(Iyy.data(), w, h, array_x2, winsize, array_y2, winsize);
  	
  	//stim::cpu2image(Iyy.data(), "data_output/Iyy_0915.bmp", w, h, stim::cmBrewer); 
  
  	delete [] array_x1;            //free the memory of array_x1
  	delete [] array_y1;            //free the memory of array_y1
  	delete [] array_x2;            //free the memory of array_x2
  	delete [] array_y2;            //free the memory of array_y2
  
  	I.set_channel(0, Ixx.data());
  	I.set_channel(1, Iyy.data());
  
  	return I;
  
  }
  
  #endif