1 /*=========================================================================
2
3 Program: Visualization Toolkit
4 Module: vtkImageContinuousDilate3D.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 "vtkImageContinuousDilate3D.h"
16
17 #include "vtkDataArray.h"
18 #include "vtkImageData.h"
19 #include "vtkImageEllipsoidSource.h"
20 #include "vtkInformation.h"
21 #include "vtkInformationVector.h"
22 #include "vtkObjectFactory.h"
23 #include "vtkPointData.h"
24 #include "vtkStreamingDemandDrivenPipeline.h"
25
26 vtkStandardNewMacro(vtkImageContinuousDilate3D);
27
28 //------------------------------------------------------------------------------
29 // Construct an instance of vtkImageContinuousDilate3D filter.
30 // By default zero values are dilated.
vtkImageContinuousDilate3D()31 vtkImageContinuousDilate3D::vtkImageContinuousDilate3D()
32 {
33 this->HandleBoundaries = 1;
34 this->KernelSize[0] = 0;
35 this->KernelSize[1] = 0;
36 this->KernelSize[2] = 0;
37
38 this->Ellipse = vtkImageEllipsoidSource::New();
39 // Setup the Ellipse to default size
40 this->SetKernelSize(1, 1, 1);
41 }
42
43 //------------------------------------------------------------------------------
~vtkImageContinuousDilate3D()44 vtkImageContinuousDilate3D::~vtkImageContinuousDilate3D()
45 {
46 if (this->Ellipse)
47 {
48 this->Ellipse->Delete();
49 this->Ellipse = nullptr;
50 }
51 }
52
53 //------------------------------------------------------------------------------
PrintSelf(ostream & os,vtkIndent indent)54 void vtkImageContinuousDilate3D::PrintSelf(ostream& os, vtkIndent indent)
55 {
56 this->Superclass::PrintSelf(os, indent);
57 }
58
59 //------------------------------------------------------------------------------
60 // This method sets the size of the neighborhood. It also sets the
61 // default middle of the neighborhood and computes the elliptical foot print.
SetKernelSize(int size0,int size1,int size2)62 void vtkImageContinuousDilate3D::SetKernelSize(int size0, int size1, int size2)
63 {
64 int modified = 0;
65
66 if (this->KernelSize[0] != size0)
67 {
68 modified = 1;
69 this->KernelSize[0] = size0;
70 this->KernelMiddle[0] = size0 / 2;
71 }
72 if (this->KernelSize[1] != size1)
73 {
74 modified = 1;
75 this->KernelSize[1] = size1;
76 this->KernelMiddle[1] = size1 / 2;
77 }
78 if (this->KernelSize[2] != size2)
79 {
80 modified = 1;
81 this->KernelSize[2] = size2;
82 this->KernelMiddle[2] = size2 / 2;
83 }
84
85 if (modified)
86 {
87 this->Modified();
88 this->Ellipse->SetWholeExtent(
89 0, this->KernelSize[0] - 1, 0, this->KernelSize[1] - 1, 0, this->KernelSize[2] - 1);
90 this->Ellipse->SetCenter(static_cast<float>(this->KernelSize[0] - 1) * 0.5,
91 static_cast<float>(this->KernelSize[1] - 1) * 0.5,
92 static_cast<float>(this->KernelSize[2] - 1) * 0.5);
93 this->Ellipse->SetRadius(static_cast<float>(this->KernelSize[0]) * 0.5,
94 static_cast<float>(this->KernelSize[1]) * 0.5, static_cast<float>(this->KernelSize[2]) * 0.5);
95
96 // make sure scalars have been allocated (needed if multithreaded is used)
97 vtkInformation* ellipseOutInfo = this->Ellipse->GetExecutive()->GetOutputInformation(0);
98 ellipseOutInfo->Set(vtkStreamingDemandDrivenPipeline::UPDATE_EXTENT(), 0,
99 this->KernelSize[0] - 1, 0, this->KernelSize[1] - 1, 0, this->KernelSize[2] - 1);
100 this->Ellipse->Update();
101 }
102 }
103
104 //------------------------------------------------------------------------------
105 // This templated function executes the filter on any region,
106 // whether it needs boundary checking or not.
107 // If the filter needs to be faster, the function could be duplicated
108 // for strictly center (no boundary ) processing.
109 template <class T>
vtkImageContinuousDilate3DExecute(vtkImageContinuousDilate3D * self,vtkImageData * mask,vtkImageData * inData,T * inPtr,vtkImageData * outData,int * outExt,T * outPtr,int id,vtkDataArray * inArray,vtkInformation * inInfo)110 void vtkImageContinuousDilate3DExecute(vtkImageContinuousDilate3D* self, vtkImageData* mask,
111 vtkImageData* inData, T* inPtr, vtkImageData* outData, int* outExt, T* outPtr, int id,
112 vtkDataArray* inArray, vtkInformation* inInfo)
113 {
114 int *kernelMiddle, *kernelSize;
115 // For looping though output (and input) pixels.
116 int outMin0, outMax0, outMin1, outMax1, outMin2, outMax2;
117 int outIdx0, outIdx1, outIdx2;
118 vtkIdType inInc0, inInc1, inInc2;
119 vtkIdType outInc0, outInc1, outInc2;
120 T *inPtr0, *inPtr1, *inPtr2;
121 T *outPtr0, *outPtr1, *outPtr2;
122 int numComps, outIdxC;
123 // For looping through hood pixels
124 int hoodMin0, hoodMax0, hoodMin1, hoodMax1, hoodMin2, hoodMax2;
125 int hoodIdx0, hoodIdx1, hoodIdx2;
126 T *hoodPtr0, *hoodPtr1, *hoodPtr2;
127 // For looping through the mask.
128 unsigned char *maskPtr, *maskPtr0, *maskPtr1, *maskPtr2;
129 vtkIdType maskInc0, maskInc1, maskInc2;
130 // The extent of the whole input image
131 int inImageMin0, inImageMin1, inImageMin2;
132 int inImageMax0, inImageMax1, inImageMax2;
133 int inImageExt[6];
134 // to compute the range
135 T pixelMax;
136 unsigned long count = 0;
137 unsigned long target;
138 int* inExt;
139
140 inExt = inData->GetExtent();
141
142 // Get information to march through data
143 inData->GetIncrements(inInc0, inInc1, inInc2);
144 inInfo->Get(vtkStreamingDemandDrivenPipeline::UPDATE_EXTENT(), inImageExt);
145 inImageMin0 = inImageExt[0];
146 inImageMax0 = inImageExt[1];
147 inImageMin1 = inImageExt[2];
148 inImageMax1 = inImageExt[3];
149 inImageMin2 = inImageExt[4];
150 inImageMax2 = inImageExt[5];
151 outData->GetIncrements(outInc0, outInc1, outInc2);
152 outMin0 = outExt[0];
153 outMax0 = outExt[1];
154 outMin1 = outExt[2];
155 outMax1 = outExt[3];
156 outMin2 = outExt[4];
157 outMax2 = outExt[5];
158 numComps = outData->GetNumberOfScalarComponents();
159
160 // Get ivars of this object (easier than making friends)
161 kernelSize = self->GetKernelSize();
162 kernelMiddle = self->GetKernelMiddle();
163 hoodMin0 = -kernelMiddle[0];
164 hoodMin1 = -kernelMiddle[1];
165 hoodMin2 = -kernelMiddle[2];
166 hoodMax0 = hoodMin0 + kernelSize[0] - 1;
167 hoodMax1 = hoodMin1 + kernelSize[1] - 1;
168 hoodMax2 = hoodMin2 + kernelSize[2] - 1;
169
170 // Setup mask info
171 maskPtr = static_cast<unsigned char*>(mask->GetScalarPointer());
172 mask->GetIncrements(maskInc0, maskInc1, maskInc2);
173
174 // in and out should be marching through corresponding pixels.
175 inPtr = static_cast<T*>(inArray->GetVoidPointer(
176 (outMin0 - inExt[0]) * inInc0 + (outMin1 - inExt[2]) * inInc1 + (outMin2 - inExt[4]) * inInc2));
177
178 target =
179 static_cast<unsigned long>(numComps * (outMax2 - outMin2 + 1) * (outMax1 - outMin1 + 1) / 50.0);
180 target++;
181
182 // loop through components
183 for (outIdxC = 0; outIdxC < numComps; ++outIdxC)
184 {
185 // loop through pixels of output
186 outPtr2 = outPtr;
187 inPtr2 = inPtr;
188 for (outIdx2 = outMin2; outIdx2 <= outMax2; ++outIdx2)
189 {
190 outPtr1 = outPtr2;
191 inPtr1 = inPtr2;
192 for (outIdx1 = outMin1; !self->AbortExecute && outIdx1 <= outMax1; ++outIdx1)
193 {
194 if (!id)
195 {
196 if (!(count % target))
197 {
198 self->UpdateProgress(count / (50.0 * target));
199 }
200 count++;
201 }
202 outPtr0 = outPtr1;
203 inPtr0 = inPtr1;
204 for (outIdx0 = outMin0; outIdx0 <= outMax0; ++outIdx0)
205 {
206 // Find min
207 pixelMax = *inPtr0;
208 // loop through neighborhood pixels
209 // as sort of a hack to handle boundaries,
210 // input pointer will be marching through data that does not exist.
211 hoodPtr2 =
212 inPtr0 - kernelMiddle[0] * inInc0 - kernelMiddle[1] * inInc1 - kernelMiddle[2] * inInc2;
213 maskPtr2 = maskPtr;
214 for (hoodIdx2 = hoodMin2; hoodIdx2 <= hoodMax2; ++hoodIdx2)
215 {
216 hoodPtr1 = hoodPtr2;
217 maskPtr1 = maskPtr2;
218 for (hoodIdx1 = hoodMin1; hoodIdx1 <= hoodMax1; ++hoodIdx1)
219 {
220 hoodPtr0 = hoodPtr1;
221 maskPtr0 = maskPtr1;
222 for (hoodIdx0 = hoodMin0; hoodIdx0 <= hoodMax0; ++hoodIdx0)
223 {
224 // A quick but rather expensive way to handle boundaries
225 if (outIdx0 + hoodIdx0 >= inImageMin0 && outIdx0 + hoodIdx0 <= inImageMax0 &&
226 outIdx1 + hoodIdx1 >= inImageMin1 && outIdx1 + hoodIdx1 <= inImageMax1 &&
227 outIdx2 + hoodIdx2 >= inImageMin2 && outIdx2 + hoodIdx2 <= inImageMax2)
228 {
229 if (*maskPtr0)
230 {
231 if (*hoodPtr0 > pixelMax)
232 {
233 pixelMax = *hoodPtr0;
234 }
235 }
236 }
237
238 hoodPtr0 += inInc0;
239 maskPtr0 += maskInc0;
240 }
241 hoodPtr1 += inInc1;
242 maskPtr1 += maskInc1;
243 }
244 hoodPtr2 += inInc2;
245 maskPtr2 += maskInc2;
246 }
247 *outPtr0 = pixelMax;
248
249 inPtr0 += inInc0;
250 outPtr0 += outInc0;
251 }
252 inPtr1 += inInc1;
253 outPtr1 += outInc1;
254 }
255 inPtr2 += inInc2;
256 outPtr2 += outInc2;
257 }
258 ++inPtr;
259 ++outPtr;
260 }
261 }
262
263 //------------------------------------------------------------------------------
264 // This method contains the first switch statement that calls the correct
265 // templated function for the input and output Data types.
266 // It handles image boundaries, so the image does not shrink.
ThreadedRequestData(vtkInformation * vtkNotUsed (request),vtkInformationVector ** inputVector,vtkInformationVector * vtkNotUsed (outputVector),vtkImageData *** inData,vtkImageData ** outData,int outExt[6],int id)267 void vtkImageContinuousDilate3D::ThreadedRequestData(vtkInformation* vtkNotUsed(request),
268 vtkInformationVector** inputVector, vtkInformationVector* vtkNotUsed(outputVector),
269 vtkImageData*** inData, vtkImageData** outData, int outExt[6], int id)
270 {
271 // return if nothing to do
272 if (outExt[1] < outExt[0] || outExt[3] < outExt[2] || outExt[5] < outExt[4])
273 {
274 return;
275 }
276
277 int inExt[6], wholeExt[6];
278
279 vtkInformation* inInfo = inputVector[0]->GetInformationObject(0);
280 inInfo->Get(vtkStreamingDemandDrivenPipeline::WHOLE_EXTENT(), wholeExt);
281 this->InternalRequestUpdateExtent(inExt, outExt, wholeExt);
282 void* inPtr;
283 void* outPtr = outData[0]->GetScalarPointerForExtent(outExt);
284 vtkImageData* mask;
285
286 vtkDataArray* inArray = this->GetInputArrayToProcess(0, inputVector);
287 // Reset later.
288 inPtr = inArray->GetVoidPointer(0);
289
290 // Error checking on mask
291 mask = this->Ellipse->GetOutput();
292 if (mask->GetScalarType() != VTK_UNSIGNED_CHAR)
293 {
294 vtkErrorMacro(<< "Execute: mask has wrong scalar type");
295 return;
296 }
297
298 // this filter expects the output type to be same as input
299 if (outData[0]->GetScalarType() != inArray->GetDataType())
300 {
301 vtkErrorMacro(<< "Execute: output ScalarType, "
302 << vtkImageScalarTypeNameMacro(outData[0]->GetScalarType())
303 << " must match input array data type");
304 return;
305 }
306
307 switch (inArray->GetDataType())
308 {
309 vtkTemplateMacro(
310 vtkImageContinuousDilate3DExecute(this, mask, inData[0][0], static_cast<VTK_TT*>(inPtr),
311 outData[0], outExt, static_cast<VTK_TT*>(outPtr), id, inArray, inInfo));
312 default:
313 vtkErrorMacro(<< "Execute: Unknown ScalarType");
314 return;
315 }
316 }
317
318 //------------------------------------------------------------------------------
RequestData(vtkInformation * request,vtkInformationVector ** inputVector,vtkInformationVector * outputVector)319 int vtkImageContinuousDilate3D::RequestData(
320 vtkInformation* request, vtkInformationVector** inputVector, vtkInformationVector* outputVector)
321 {
322 this->Ellipse->Update();
323 return this->Superclass::RequestData(request, inputVector, outputVector);
324 }
325