Blame view

scalarslice.cu 3.2 KB
3f56f1f9   dmayerich   initial commit
1
2
  #include "scalarslice.h"

  

d6f53e68   dmayerich   rts organization
3
  #include "rts/cuda/error.h"
3f56f1f9   dmayerich   initial commit
4
  #include "cublas_v2.h"
d6f53e68   dmayerich   rts organization
5
  #include "rts/envi/envi.h"

3f56f1f9   dmayerich   initial commit
6
7
8
9
10
11
12
13
  

  scalarslice::scalarslice(int x, int y)

  {

  	//set the resolution

  	R[0] = x;

  	R[1] = y;

  

  	//allocate memory on the GPU

3f36b18e   David Mayerich   Adding planewave ...
14
15
  	HANDLE_ERROR(cudaMalloc( (void**)&S, sizeof(ptype) * x * y ));
  	//std::cout<<"Scalerslice created."<<std::endl;

3f56f1f9   dmayerich   initial commit
16
17
18
19
20
21
  }
  
  scalarslice::scalarslice()
  {
      R[0] = R[1] = 0;
      S = NULL;
3f36b18e   David Mayerich   Adding planewave ...
22
23
  
      //std::cout<<"Scalerslice created (default)."<<std::endl;
3f56f1f9   dmayerich   initial commit
24
25
26
27
  }

  

  scalarslice::~scalarslice()

  {

51b6469a   dmayerich   added look-up tables
28
29
  	if(S != NULL)

  		HANDLE_ERROR(cudaFree(S));

3f36b18e   David Mayerich   Adding planewave ...
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
  	S = NULL;
  	R[0] = R[1] = 0;
  
  	//std::cout<<"Scalerslice destroyed."<<std::endl;

  }
  
  //assignment operator
  scalarslice & scalarslice::operator= (const scalarslice & rhs)
  {
      //de-allocate any existing GPU memory
      if(S != NULL)

  		HANDLE_ERROR(cudaFree(S));
  
      //copy the slice resolution
      R[0] = rhs.R[0];
      R[1] = rhs.R[1];
  
      //allocate the necessary memory
      HANDLE_ERROR(cudaMalloc(&S, sizeof(ptype) * R[0] * R[1]));
  
      //copy the slice
      HANDLE_ERROR(cudaMemcpy(S, rhs.S, sizeof(ptype) * R[0] * R[1], cudaMemcpyDeviceToDevice));
  
  
      std::cout<<"Assignment operator."<<std::endl;
  
  	return *this;
  
3f56f1f9   dmayerich   initial commit
58
59
  }

  

51b6469a   dmayerich   added look-up tables
60
  void scalarslice::toImage(std::string filename, ptype vmin, ptype vmax, rts::colormapType cmap)

3f56f1f9   dmayerich   initial commit
61
  {

51b6469a   dmayerich   added look-up tables
62
  	rts::gpu2image<ptype>(S, filename, R[0], R[1], vmin, vmax, cmap);

3f56f1f9   dmayerich   initial commit
63
64
  }
  
51b6469a   dmayerich   added look-up tables
65
  void scalarslice::toImage(std::string filename, bool positive, rts::colormapType cmap)
3f56f1f9   dmayerich   initial commit
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
  {
      cublasStatus_t stat;
      cublasHandle_t handle;
  
      //create a CUBLAS handle
      stat = cublasCreate(&handle);
      if(stat != CUBLAS_STATUS_SUCCESS)
      {
          std::cout<<"CUBLAS Error: initialization failed"<<std::endl;
          exit(1);
      }
  
      //find the index of the value with maximum magnitude
      int N = R[0] * R[1];
      int result;
e70a251f   dmayerich   fixed double prec...
81
  #ifdef PRECISION_SINGLE
3f56f1f9   dmayerich   initial commit
82
      stat = cublasIsamax(handle, N, S, 1, &result);
e70a251f   dmayerich   fixed double prec...
83
84
85
86
87
88
89
90
91
92
93
94
  #elif defined PRECISION_DOUBLE
  	stat = cublasIdamax(handle, N, S, 1, &result);
  #endif
  
  	//adjust for 1-based indexing
  	result -= 1;
  
  	if(stat != CUBLAS_STATUS_SUCCESS)
  	{
  		std::cout<<"CUBLAS Error: failure finding maximum value."<<std::endl;
  		exit(1);
  	}
3f56f1f9   dmayerich   initial commit
95
  
51b6469a   dmayerich   added look-up tables
96
  
3f56f1f9   dmayerich   initial commit
97
98
99
100
  
      //retrieve the maximum value
      ptype maxVal;
      HANDLE_ERROR(cudaMemcpy(&maxVal, S + result, sizeof(ptype), cudaMemcpyDeviceToHost));
3f56f1f9   dmayerich   initial commit
101
102
103
104
105
106
107
108
  
      //destroy the CUBLAS handle
      cublasDestroy(handle);
  
      //output the image
      if(positive)
          toImage(filename, 0, maxVal, cmap);
      else
51b6469a   dmayerich   added look-up tables
109
          toImage(filename, -abs(maxVal), abs(maxVal), cmap);
3f56f1f9   dmayerich   initial commit
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
  }
  
  void scalarslice::toEnvi(std::string filename, ptype wavelength, bool append)
  {
      std::string mode;
      if(append) mode = "a";
      else       mode = "w";
  
      //open the ENVI file
      EnviFile outfile(filename, mode);
  
      //get the scalar slice from the GPU to the CPU
      int memsize = sizeof(ptype) * R[0] * R[1];
      ptype* cpuData = (ptype*) malloc( memsize );
      HANDLE_ERROR(cudaMemcpy( cpuData, S, memsize, cudaMemcpyDeviceToHost));
  
3f56f1f9   dmayerich   initial commit
126
127
      //add a band to the ENVI file
      outfile.addBand(cpuData, R[0], R[1], wavelength);
3f56f1f9   dmayerich   initial commit
128
129
130
131
132
133
134
135
136
137
138
139
140
  
      outfile.close();
  
  
  }

  

  void scalarslice::clear()

  {

  	//this function sets the slice to zero

  	if(S != NULL)

  		//HANDLE_ERROR(cudaMalloc( (void**)&S, sizeof(ptype) * R[0] * R[1] ));

  		HANDLE_ERROR(cudaMemset(S, 0, sizeof(ptype) * R[0] * R[1]));

  }