Done... Bug #2234<br>
<br>
<a href="http://www.vtk.org/Bug/bug.php?op=show&amp;bugid=2234&amp;pos=2">http://www.vtk.org/Bug/bug.php?op=show&amp;bugid=2234&amp;pos=2</a><br><br><div><span class="gmail_quote">On 9/14/05, <b class="gmail_sendername">
Mathieu Malaterre</b> &lt;<a href="mailto:mathieu.malaterre@kitware.com">mathieu.malaterre@kitware.com</a>&gt; wrote:</span><blockquote class="gmail_quote" style="border-left: 1px solid rgb(204, 204, 204); margin: 0pt 0pt 0pt 0.8ex; padding-left: 1ex;">
Ok this is a bug:<br><br>i ========== End of Data Statistics ===============<br>==30679==<br>==30679== ERROR SUMMARY: 1 errors from 1 contexts (suppressed: 9 from 4)<br>==30679==<br>==30679== 1 errors in context 1 of 1:<br>
==30679== Invalid read of size 8<br>==30679==&nbsp;&nbsp;&nbsp;&nbsp;at 0x12A4B8F5: void<br>vtkImageFFTExecute&lt;double&gt;(vtkImageFFT*, vtkImageData*, int*, double*,<br>vtkImageData*, int*, double*, int) (vtkImageFFT.cxx:137)<br>==30679==&nbsp;&nbsp;&nbsp;&nbsp;by 0x12A4A645: vtkImageFFT::ThreadedExecute(vtkImageData*,
<br>vtkImageData*, int*, int) (vtkImageFFT.cxx:203)<br>==30679==&nbsp;&nbsp;&nbsp;&nbsp;by 0x1353AEA2:<br>vtkThreadedImageAlgorithm::ThreadedRequestData(vtkInformation*,<br>vtkInformationVector**, vtkInformationVector*, vtkImageData***,<br>vtkImageData**, int*, int) (
vtkThreadedImageAlgorithm.cxx:327)<br>==30679==&nbsp;&nbsp;&nbsp;&nbsp;by 0x1353AA8E:<br>vtkThreadedImageAlgorithmThreadedExecute(void*)<br>(vtkThreadedImageAlgorithm.cxx:200)<br>==30679==&nbsp;&nbsp;&nbsp;&nbsp;by 0x13BA30F9: start_thread (in /lib64/tls/libpthread-
<a href="http://2.3.5.so">2.3.5.so</a>)<br>==30679==&nbsp;&nbsp;&nbsp;&nbsp;by 0x1561AA62: clone (in /lib64/tls/libc-<a href="http://2.3.5.so">2.3.5.so</a>)<br>==30679==&nbsp;&nbsp;Address 0x19992030 is 0 bytes after a block of size 16777216<br>alloc'd
<br>==30679==&nbsp;&nbsp;&nbsp;&nbsp;at 0x11B1A8E2: operator new[](unsigned long) (in<br>/usr/lib64/valgrind/vgpreload_memcheck.so)<br>==30679==&nbsp;&nbsp;&nbsp;&nbsp;by 0x1380B136: vtkDataArrayTemplate&lt;double&gt;::Allocate(long<br>long, long long) (vtkDataArrayTemplate.txx
:99)<br>==30679==&nbsp;&nbsp;&nbsp;&nbsp;by 0x1380BF2E:<br>vtkDataArrayTemplate&lt;double&gt;::SetNumberOfValues(long long)<br>(vtkDataArrayTemplate.txx:536)<br>==30679==&nbsp;&nbsp;&nbsp;&nbsp;by 0x1380B409:<br>vtkDataArrayTemplate&lt;double&gt;::SetNumberOfTuples(long long)
<br>(vtkDataArrayTemplate.txx:330)<br>==30679==&nbsp;&nbsp;&nbsp;&nbsp;by 0x12154566:<br>vtkArrayCalculator::RequestData(vtkInformation*, vtkInformationVector**,<br>vtkInformationVector*) (vtkArrayCalculator.cxx:279)<br>==30679==&nbsp;&nbsp;&nbsp;&nbsp;by 0x1341A5AA:
<br>vtkDataSetAlgorithm::ProcessRequest(vtkInformation*,<br>vtkInformationVector**, vtkInformationVector*)<br>(vtkDataSetAlgorithm.cxx:175)==30679==&nbsp;&nbsp;&nbsp;&nbsp;by 0x1342CBC9:<br>vtkExecutive::CallAlgorithm(vtkInformation*, int,<br>
vtkInformationVector**, vtkInformationVector*) (vtkExecutive.cxx:679)<br>==30679==&nbsp;&nbsp;&nbsp;&nbsp;by 0x13427AF7:<br>vtkDemandDrivenPipeline::ExecuteData(vtkInformation*,<br>vtkInformationVector**, vtkInformationVector*)<br>(vtkDemandDrivenPipeline.cxx
:464)<br>==30679==&nbsp;&nbsp;&nbsp;&nbsp;by 0x1342732B:<br>vtkDemandDrivenPipeline::ProcessRequest(vtkInformation*, int,<br>vtkInformationVector**, vtkInformationVector*)<br>(vtkDemandDrivenPipeline.cxx:244)<br>==30679==&nbsp;&nbsp;&nbsp;&nbsp;by 0x1352487E:<br>
vtkStreamingDemandDrivenPipeline::ProcessRequest(vtkInformation*, int,<br>vtkInformationVector**, vtkInformationVector*)<br>(vtkStreamingDemandDrivenPipeline.cxx:144)<br>==30679==&nbsp;&nbsp;&nbsp;&nbsp;by 0x13427927: vtkDemandDrivenPipeline::UpdateData(int)
<br>(vtkDemandDrivenPipeline.cxx:406)<br>==30679==&nbsp;&nbsp;&nbsp;&nbsp;by 0x13524A04:<br>vtkStreamingDemandDrivenPipeline::Update(int)<br>(vtkStreamingDemandDrivenPipeline.cxx:190)<br><br><br>Could you open a bug in the bugtracker(*) for this. Add your program as
<br>an attachment.<br><br>Thanks<br>Mathieu<br>(*) <a href="http://vtk.org/Bug">http://vtk.org/Bug</a><br><br><br><br>Randall Hand wrote:<br>&gt; Ok, I posted this once before, but here it is again with some more<br>&gt; information.
<br>&gt;<br>&gt; I've written a piece of code to create a synthetic dataset.&nbsp;&nbsp;It's a<br>&gt; 128^3 dataset, each point with a 3-component double value consisting of<br>&gt; (sin(x), junk, morejunk).&nbsp;&nbsp;I extract just the first component, and do an
<br>&gt; FFT.&nbsp;&nbsp;You'll find the code at the end of this message, along with my<br>&gt; &quot;PrintStatistics&quot; function attached.<br>&gt;<br>&gt; You'll see from the #if/#endif blocks that there's two possible ways to<br>
&gt; do it.&nbsp;&nbsp;The first is with a vtkArrayCalculator, where I set a very<br>&gt; trivial function to simply copy out the one value.&nbsp;&nbsp;The other way is<br>&gt; with a vtkImageExtractComponents.&nbsp;&nbsp;When doing it with the<br>&gt; vtkArrayCalculator, the filter seems to work just fine but then the FFT
<br>&gt; yields values like :<br>&gt; p&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; Array [0]: vtkDoubleArray &quot;result&quot;, 2,097,152 points, 8 bytes<br>&gt; per point<br>&gt;
p&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;
Magnitude
Range:&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;0.0000&nbsp;&nbsp;1482912.5793<br>&gt;
p&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;
Component 0
Range:&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;-1048575.9848&nbsp;&nbsp;
1048575.8760<br>&gt;
p&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;
Component 1
Range:&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;-1048572.9243&nbsp;&nbsp;
1048579.0968<br>&gt; With a simple sin wave, the result should be almost purely imaginary, so<br>&gt; the Component 0 range should be almost 0 to 0.<br>&gt;<br>&gt; When I switch to using the vtkImageExtractComponents, the results look
<br>&gt; better:<br>&gt; p&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; Array [0]: vtkDoubleArray &quot;image&quot;, 2,097,152 points, 8 bytes per<br>&gt; point<br>&gt;
p&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;
Magnitude
Range:&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;0.0000&nbsp;&nbsp;1048576.0091<br>&gt;
p&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;
Component 0 Range:&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;-0.0557 0.0011<br>&gt;
p&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;
Component 1
Range:&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;-1048576.0089&nbsp;&nbsp;
1048576.0091<br>&gt;<br>&gt;<br>&gt; Unfortunately, when working with the vtkImageExtractComponents filter, I<br>&gt; lose the other two fields if I had planned on using them for coloring or<br>&gt; something.&nbsp;&nbsp;Also, my results aren't quite what I expected.&nbsp;&nbsp;Much better,
<br>&gt; but on more complex datasets I seem to get results that are a little off<br>&gt; from the expected result.&nbsp;&nbsp;I don't know if the cause of all this is<br>&gt; related, but I'm chasing down what I can.<br>&gt;<br>&gt; So, any ideas?
<br>&gt;<br>&gt; ----- Code starts here -------<br>&gt; #include &lt;stdio.h&gt;<br>&gt; #include &lt;string.h&gt;<br>&gt; #include &lt;stdlib.h&gt;<br>&gt; #include &lt;unistd.h&gt;<br>&gt; #include &quot;stat.h&quot;<br>
&gt; #include &lt;sys/types.h&gt;<br>&gt; #include &lt;sys/stat.h&gt;<br>&gt; #include &lt;sys/mman.h&gt;<br>&gt; #include &lt;fcntl.h&gt;<br>&gt; #include &lt;vtkObject.h&gt;<br>&gt; #include &quot;vtkImageData.h&quot;<br>
&gt; #include &quot;vtkPointData.h&quot;<br>&gt; #include &lt;vtkDataArray.h&gt;<br>&gt; #include &lt;vtkUnsignedCharArray.h&gt;<br>&gt; #include &lt;vtkDoubleArray.h&gt;<br>&gt; #include &lt;vtkImageData.h&gt;<br>&gt; #include &lt;
vtkImageFFT.h&gt;<br>&gt; #include &lt;vtkArrayCalculator.h&gt;<br>&gt; #include &lt;vtkImageExtractComponents.h&gt;<br>&gt;<br>&gt; int main (int argc, char *argv[])<br>&gt; {<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; double *data;<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; int xdim, ydim, zdim;
<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; double kx;<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; int x,y,z;<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; long index;<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; xdim = 128;<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; ydim = 128;<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; zdim = 128;<br>&gt;<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; if ((data = (double*)malloc(sizeof(double)*xdim*ydim*zdim*3)) == NULL) {
<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; perror(&quot;Unable to allocate memory!\n\t&quot;);<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; return -1;<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; }<br>&gt;<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; printf(&quot;Constructing sin(kx * x) dataset...\n&quot;);<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; kx = 2.0 * 3.1415926
 / (double)xdim;<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; for(x=0; x&lt;xdim; x++) {<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; for(y=0; y&lt;ydim; y++) {<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; for(z=0; z&lt;zdim; z++) {<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;
index = z + (y*zdim) + (x * ydim * zdim);<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;
data[(index*3)+0] = sin(kx * (double)x);<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;
data[(index*3)+1] = cos(kx * (double)x);<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;
data[(index*3)+2] = tan(kx * (double)x);<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; }<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; }<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; }<br>&gt;<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; vtkDoubleArray *ucPointer = vtkDoubleArray::New();<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; ucPointer-&gt;SetNumberOfComponents(3);<br>
&gt;&nbsp;&nbsp;&nbsp;&nbsp; ucPointer-&gt;SetArray(data, xdim*ydim*zdim*3, 1);<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; ucPointer-&gt;SetName(&quot;image&quot;);<br>&gt;<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; vtkImageData *image = vtkImageData::New();<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; image-&gt;Initialize();<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; image-&gt;SetDimensions(xdim, ydim, zdim);
<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; image-&gt;SetExtent(1, xdim, 1, ydim, 1, zdim);<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; image-&gt;SetScalarTypeToDouble();<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; image-&gt;SetNumberOfScalarComponents(3);<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; image-&gt;GetPointData()-&gt;SetScalars(ucPointer);
<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; PrintStatistics(image);<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; printf(&quot;Extracting Velocity X-Component...\n&quot;);<br>&gt;<br>&gt; #if 0<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; vtkArrayCalculator *extract = vtkArrayCalculator::New();<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; extract-&gt;SetInput(image);
<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; extract-&gt;SetResultArrayName(&quot;result&quot;);<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; extract-&gt;AddScalarVariable(&quot;x&quot;, &quot;image&quot;, 0);<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; extract-&gt;SetFunction(&quot;(x*1.0)+0.0&quot;);<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; extract-&gt;ReplaceInvalidValuesOff();
<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; extract-&gt;Update();<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; extract-&gt;GetOutput()-&gt;GetPointData()-&gt;RemoveArray(&quot;image&quot;);<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; extract-&gt;GetOutput()-&gt;GetPointData()-&gt;SetActiveScalars(&quot;result&quot;);
<br>&gt; #else<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; vtkImageExtractComponents *extract = vtkImageExtractComponents::New();<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; extract-&gt;SetInput(image);<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; extract-&gt;SetComponents(0);<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; extract-&gt;Update();<br>&gt; #endif
<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; PrintStatistics(extract-&gt;GetOutput());<br>&gt;<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; printf(&quot;Computing FFT...\n&quot;);<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; vtkImageFFT *fft = vtkImageFFT::New();<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; fft-&gt;SetInput(extract-&gt;GetOutput());<br>
&gt;&nbsp;&nbsp;&nbsp;&nbsp; fft-&gt;Update();<br>&gt;&nbsp;&nbsp;&nbsp;&nbsp; PrintStatistics(fft-&gt;GetOutput());<br>&gt;<br>&gt; }<br>&gt;<br>&gt;<br>&gt; --<br>&gt; Randall Hand<br>&gt; Visualization Scientist,<br>&gt; ERDC-MSRC Vicksburg, MS<br>&gt; Homepage: 
<a href="http://www.yeraze.com">http://www.yeraze.com</a><br>&gt;<br>&gt;<br>&gt; ------------------------------------------------------------------------<br>&gt;<br>&gt; _______________________________________________<br>
&gt; This is the private VTK discussion list.<br>&gt; Please keep messages on-topic. Check the FAQ at: <a href="http://www.vtk.org/Wiki/VTK_FAQ">http://www.vtk.org/Wiki/VTK_FAQ</a><br>&gt; Follow this link to subscribe/unsubscribe:
<br>&gt; <a href="http://www.vtk.org/mailman/listinfo/vtkusers">http://www.vtk.org/mailman/listinfo/vtkusers</a><br><br></blockquote></div><br><br clear="all"><br>-- <br>Randall Hand<br>Visualization Scientist, <br>ERDC-MSRC Vicksburg, MS
<br>Homepage: <a href="http://www.yeraze.com">http://www.yeraze.com</a>