1 /*=========================================================================
2
3 Program: Visualization Toolkit
4 Module: vtkEncodedGradientEstimator.cxx
5
6 Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
7 All rights reserved.
8 See Copyright.txt or http://www.kitware.com/Copyright.htm for details.
9
10 This software is distributed WITHOUT ANY WARRANTY; without even
11 the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
12 PURPOSE. See the above copyright notice for more information.
13
14 =========================================================================*/
15 #include "vtkEncodedGradientEstimator.h"
16
17 #include "vtkGarbageCollector.h"
18 #include "vtkImageData.h"
19 #include "vtkMultiThreader.h"
20 #include "vtkRecursiveSphereDirectionEncoder.h"
21 #include "vtkTimerLog.h"
22
23 #include <cmath>
24
25
26 vtkCxxSetObjectMacro(vtkEncodedGradientEstimator, InputData, vtkImageData );
27
28 // Construct a vtkEncodedGradientEstimator with initial values of nullptr for
29 // the Input, EncodedNormal, and GradientMagnitude. Also,
30 // indicate that the IndexTable has not yet been initialized. The
31 // GradientMagnitudeRange and the GradientMangitudeTable are
32 // initialized to default values - these will change in the future
33 // when magnitude of gradient opacities are included
vtkEncodedGradientEstimator()34 vtkEncodedGradientEstimator::vtkEncodedGradientEstimator()
35 {
36 this->InputData = nullptr;
37 this->EncodedNormals = nullptr;
38 this->EncodedNormalsSize[0] = 0;
39 this->EncodedNormalsSize[1] = 0;
40 this->EncodedNormalsSize[2] = 0;
41 this->GradientMagnitudes = nullptr;
42 this->GradientMagnitudeScale = 1.0;
43 this->GradientMagnitudeBias = 0.0;
44 this->Threader = vtkMultiThreader::New();
45 this->NumberOfThreads = this->Threader->GetNumberOfThreads();
46 this->DirectionEncoder = vtkRecursiveSphereDirectionEncoder::New();
47 this->ComputeGradientMagnitudes = 1;
48 this->CylinderClip = 0;
49 this->CircleLimits = nullptr;
50 this->CircleLimitsSize = -1;
51 this->UseCylinderClip = 0;
52 this->LastUpdateTimeInSeconds = -1.0;
53 this->LastUpdateTimeInCPUSeconds = -1.0;
54 this->ZeroNormalThreshold = 0.0;
55 this->ZeroPad = 1;
56 this->BoundsClip = 0;
57 this->Bounds[0] =
58 this->Bounds[1] =
59 this->Bounds[2] =
60 this->Bounds[3] =
61 this->Bounds[4] =
62 this->Bounds[5] = 0;
63
64 }
65
66 // Destruct a vtkEncodedGradientEstimator - free up any memory used
~vtkEncodedGradientEstimator()67 vtkEncodedGradientEstimator::~vtkEncodedGradientEstimator()
68 {
69 this->SetInputData(nullptr);
70 this->Threader->Delete();
71 this->Threader = nullptr;
72
73 delete [] this->EncodedNormals;
74
75 delete [] this->GradientMagnitudes;
76
77 if ( this->DirectionEncoder )
78 {
79 this->DirectionEncoder->UnRegister( this );
80 }
81
82 delete [] this->CircleLimits;
83 }
84
SetZeroNormalThreshold(float v)85 void vtkEncodedGradientEstimator::SetZeroNormalThreshold( float v )
86 {
87 if ( this->ZeroNormalThreshold != v )
88 {
89 if ( v < 0.0 )
90 {
91 vtkErrorMacro( << "The ZeroNormalThreshold must be a value >= 0.0" );
92 return;
93 }
94
95 this->ZeroNormalThreshold = v;
96 this->Modified();
97 }
98 }
99
100 void
SetDirectionEncoder(vtkDirectionEncoder * direnc)101 vtkEncodedGradientEstimator::SetDirectionEncoder(vtkDirectionEncoder *direnc)
102 {
103 // If we are setting it to its current value, don't do anything
104 if ( this->DirectionEncoder == direnc )
105 {
106 return;
107 }
108
109 // If we already have a direction encoder, unregister it.
110 if ( this->DirectionEncoder )
111 {
112 this->DirectionEncoder->UnRegister(this);
113 this->DirectionEncoder = nullptr;
114 }
115
116 // If we are passing in a non-nullptr encoder, register it
117 if ( direnc )
118 {
119 direnc->Register( this );
120 }
121
122 // Actually set the encoder, and consider the object Modified
123 this->DirectionEncoder = direnc;
124 this->Modified();
125 }
126
GetEncodedNormalIndex(vtkIdType xyzIndex)127 int vtkEncodedGradientEstimator::GetEncodedNormalIndex( vtkIdType xyzIndex )
128 {
129 this->Update();
130 return *(this->EncodedNormals + xyzIndex);
131 }
132
GetEncodedNormalIndex(int xIndex,int yIndex,int zIndex)133 int vtkEncodedGradientEstimator::GetEncodedNormalIndex( int xIndex,
134 int yIndex,
135 int zIndex )
136 {
137 vtkIdType ystep, zstep;
138
139 this->Update();
140
141 // Compute steps through the volume in x, y, and z
142 ystep = this->InputSize[0];
143 zstep = ystep * this->InputSize[1];
144
145 return *(this->EncodedNormals + zIndex * zstep + yIndex * ystep + xIndex);
146 }
147
GetEncodedNormals()148 unsigned short *vtkEncodedGradientEstimator::GetEncodedNormals()
149 {
150 this->Update();
151
152 return this->EncodedNormals;
153 }
154
GetGradientMagnitudes()155 unsigned char *vtkEncodedGradientEstimator::GetGradientMagnitudes()
156 {
157 this->Update();
158
159 return this->GradientMagnitudes;
160 }
161
Update()162 void vtkEncodedGradientEstimator::Update( )
163 {
164 int scalarInputSize[3];
165 double scalarInputAspect[3];
166 double startSeconds, endSeconds;
167 double startCPUSeconds, endCPUSeconds;
168
169 if ( !this->InputData )
170 {
171 vtkErrorMacro(<< "No input in gradient estimator.");
172 return;
173 }
174
175 if ( this->GetMTime() > this->BuildTime ||
176 this->DirectionEncoder->GetMTime() > this->BuildTime ||
177 this->InputData->GetMTime() > this->BuildTime ||
178 !this->EncodedNormals )
179 {
180
181 startSeconds = vtkTimerLog::GetUniversalTime();
182 startCPUSeconds = vtkTimerLog::GetCPUTime();
183
184 // Get the dimensions of the data and its aspect ratio
185 this->InputData->GetDimensions( scalarInputSize );
186 this->InputData->GetSpacing( scalarInputAspect );
187
188 // If we previously have allocated space for the encoded normals,
189 // and this space is no longer the right size, delete it
190 if ( this->EncodedNormalsSize[0] != scalarInputSize[0] ||
191 this->EncodedNormalsSize[1] != scalarInputSize[1] ||
192 this->EncodedNormalsSize[2] != scalarInputSize[2] )
193 {
194 delete [] this->EncodedNormals;
195 this->EncodedNormals = nullptr;
196
197 delete [] this->GradientMagnitudes;
198 this->GradientMagnitudes = nullptr;
199 }
200
201 // Compute the number of encoded voxels
202 vtkIdType encodedSize = scalarInputSize[0];
203 encodedSize *= scalarInputSize[1];
204 encodedSize *= scalarInputSize[2];
205
206 // Allocate space for the encoded normals if necessary
207 if ( !this->EncodedNormals )
208 {
209 this->EncodedNormals = new unsigned short[ encodedSize ];
210 this->EncodedNormalsSize[0] = scalarInputSize[0];
211 this->EncodedNormalsSize[1] = scalarInputSize[1];
212 this->EncodedNormalsSize[2] = scalarInputSize[2];
213 }
214
215 if ( !this->GradientMagnitudes && this->ComputeGradientMagnitudes )
216 {
217 this->GradientMagnitudes = new unsigned char[ encodedSize ];
218 }
219
220 // Copy info that multi threaded function will need into temp variables
221 memcpy( this->InputSize, scalarInputSize, 3 * sizeof(int) );
222 // TODO cleanup when double changes are further along
223 this->InputAspect[0] = static_cast<float>(scalarInputAspect[0]);
224 this->InputAspect[1] = static_cast<float>(scalarInputAspect[1]);
225 this->InputAspect[2] = static_cast<float>(scalarInputAspect[2]);
226 // memcpy( this->InputAspect, scalarInputAspect, 3 * sizeof(float) );
227
228 if ( this->CylinderClip &&
229 (this->InputSize[0] == this->InputSize[1]) )
230 {
231 this->UseCylinderClip = 1;
232 this->ComputeCircleLimits( this->InputSize[0] );
233 }
234 else
235 {
236 this->UseCylinderClip = 0;
237 }
238 this->UpdateNormals();
239
240 this->BuildTime.Modified();
241
242 endSeconds = vtkTimerLog::GetUniversalTime();
243 endCPUSeconds = vtkTimerLog::GetCPUTime();
244
245 this->LastUpdateTimeInSeconds = static_cast<float>(endSeconds - startSeconds);
246 this->LastUpdateTimeInCPUSeconds = static_cast<float>(endCPUSeconds - startCPUSeconds);
247 }
248 }
249
ComputeCircleLimits(int size)250 void vtkEncodedGradientEstimator::ComputeCircleLimits( int size )
251 {
252 int *ptr, y;
253 double w, halfsize, length, start, end;
254
255 if ( this->CircleLimitsSize != size )
256 {
257 delete [] this->CircleLimits;
258 this->CircleLimits = new int[2*size];
259 this->CircleLimitsSize = size;
260 }
261
262 ptr = this->CircleLimits;
263
264 halfsize = (size-1)/2.0;
265
266 for ( y = 0; y < size; y++ )
267 {
268 w = halfsize - y;
269 length = static_cast<int>( sqrt( (halfsize*halfsize) - (w*w) ) + 0.5 );
270 start = halfsize - length - 1;
271 end = halfsize + length + 1;
272 start = (start<0)?(0):(start);
273 end = (end>(size-1))?(size-1):(end);
274
275 *(ptr++) = static_cast<int>(start);
276 *(ptr++) = static_cast<int>(end);
277 }
278 }
279
280 // Print the vtkEncodedGradientEstimator
PrintSelf(ostream & os,vtkIndent indent)281 void vtkEncodedGradientEstimator::PrintSelf(ostream& os, vtkIndent indent)
282 {
283 this->Superclass::PrintSelf(os, indent);
284
285 if ( this->InputData )
286 {
287 os << indent << "InputData: (" << this->InputData << ")\n";
288 }
289 else
290 {
291 os << indent << "Input: (none)\n";
292 }
293
294 if ( this->DirectionEncoder )
295 {
296 os << indent << "DirectionEncoder: (" << this->DirectionEncoder << ")\n";
297 }
298 else
299 {
300 os << indent << "DirectionEncoder: (none)\n";
301 }
302
303 os << indent << "Build Time: "
304 << this->BuildTime.GetMTime() << endl;
305
306 os << indent << "Gradient Magnitude Scale: "
307 << this->GradientMagnitudeScale << endl;
308
309 os << indent << "Gradient Magnitude Bias: "
310 << this->GradientMagnitudeBias << endl;
311
312 os << indent << "Zero Pad: "
313 << ((this->ZeroPad)?"On":"Off") << endl;
314
315 os << indent << "Bounds Clip: "
316 << ((this->BoundsClip)?"On":"Off") << endl;
317
318 os << indent << "Bounds: ("
319 << this->Bounds[0] << ", " << this->Bounds[1] << ", "
320 << this->Bounds[2] << ", " << this->Bounds[3] << ", "
321 << this->Bounds[4] << ", " << this->Bounds[5] << ")\n";
322
323 os << indent << "Zero Normal Threshold: "
324 << this->ZeroNormalThreshold << endl;
325
326 os << indent << "Compute Gradient Magnitudes: "
327 << ((this->ComputeGradientMagnitudes)?"On":"Off") << endl;
328
329 os << indent << "Cylinder Clip: "
330 << ((this->CylinderClip)?"On":"Off") << endl;
331
332 os << indent << "Number Of Threads: "
333 << this->NumberOfThreads << endl;
334
335 os << indent << "Last Update Time In Seconds: "
336 << this->LastUpdateTimeInSeconds << endl;
337
338 os << indent << "Last Update Time In CPU Seconds: "
339 << this->LastUpdateTimeInCPUSeconds << endl;
340
341 // I don't want to print out these variables - they are
342 // internal and the get methods are included only for access
343 // within the threaded function
344 // os << indent << "Use Cylinder Clip: "
345 // << this->UseCylinderClip << endl;
346 // os << indent << " Input Size: "
347 // << this->InputSize << endl;
348 // os << indent << " Input Aspect Clip: "
349 // << this->InputAspect << endl;
350
351 }
352
353 //----------------------------------------------------------------------------
354 void
ReportReferences(vtkGarbageCollector * collector)355 vtkEncodedGradientEstimator::ReportReferences(vtkGarbageCollector* collector)
356 {
357 this->Superclass::ReportReferences(collector);
358 vtkGarbageCollectorReport(collector, this->InputData, "Input");
359 }
360