Hi,<br><br><br>than may be I set my origin wrong. Here the sample code with comments.<br><br>    typedef float RealType;<br>    <br>// It&#39;s a volume<br>    typedef itk::Vector&lt;RealType, 3&gt; VectorType;<br>// The deformation is over time<br>
    typedef itk ::Image &lt;VectorType , 4 &gt; TimeVaryingDeformationFieldType;<br>    typedef itk ::PointSet &lt;VectorType , 4 &gt; TimeVaryingDeformationFieldPointSetType;<br><br>    TimeVaryingDeformationFieldPointSetType::Pointer fieldPoints = TimeVaryingDeformationFieldPointSetType:: New ();<br>
    fieldPoints-&gt;Initialize();<br><br>    TimeVaryingDeformationFieldType::PointType        origin;<br>    TimeVaryingDeformationFieldType::SpacingType    spacing;<br>    TimeVaryingDeformationFieldType::SizeType        size;<br>
// I want to have a sapcing of 10 mm in each dimension except in the temporal there it is 2 <br>    // Define the parametric domain of BSpline Field<br>    spacing[0] = 10.0f;<br>    spacing[1] = 10.0f;<br>    spacing[2] = 10.0f;<br>
    spacing[3] = 2.0f;<br>// Here the size of the BSpline field, I have 133 numberOfTimePoints<br>    size[0] = 128;<br>    size[1] = 128;<br>    size[2] = 128;<br>    size[3] = 67;<br>    <br>    <br>// The origin ( where I&#39;m not sure, I also tried 0,0,0)<br>
    float O[3] ={0.0f,0.0f,0.0f};<br>    for (int i=0; i&lt;3; i++)<br>            O[i] = m_configuration.origin[i] - spacing[i] * (size[i] - 1) *0.5f;<br><br>    origin[0] = O[0];<br>    origin[1] = O[1];<br>    origin[2] = O[2];<br>
    origin[3] = 0.0f;<br><br>    // Set the sample points<br>    for ( int t = 0; t &lt; numberOfTimePoints ; t++ ) <br>    {<br>        int num = 0;<br>        for( int idx = 0; idx &lt; <a href="http://endo_points.at">endo_points.at</a>(t).size()-3; idx+=3 )<br>
        {<br>            <br>            VectorType v;<br>// The sampling points range from -70 to +70 for example<br>            v[0] = <a href="http://endo_points.at">endo_points.at</a>(t)[idx]   - <a href="http://endo_points.at">endo_points.at</a>(0)[idx]  ;<br>
            v[1] = <a href="http://endo_points.at">endo_points.at</a>(t)[idx+1] - <a href="http://endo_points.at">endo_points.at</a>(0)[idx+1];<br>            v[2] = <a href="http://endo_points.at">endo_points.at</a>(t)[idx+2] - <a href="http://endo_points.at">endo_points.at</a>(0)[idx+2];<br>
<br>            <br>            TimeVaryingDeformationFieldPointSetType::PointType point;<br>            point[0] = <a href="http://endo_points.at">endo_points.at</a>(0)[idx];<br>            point[1] = <a href="http://endo_points.at">endo_points.at</a>(0)[idx+1];<br>
            point[2] = <a href="http://endo_points.at">endo_points.at</a>(0)[idx+2];<br>            point[3] = t;<br><br>            fieldPoints-&gt;SetPoint( num, point );<br>            fieldPoints-&gt;SetPointData(num, v);<br>
            num++;<br><br>        }<br>    }<br><br><br>    <br><br>    // Define the filter and set the parameters<br>    typedef itk::myBSplineScatteredDataPointSetToImageFilter&lt;TimeVaryingDeformationFieldPointSetType, TimeVaryingDeformationFieldType&gt; BSplineFilterType;<br>
    BSplineFilterType::Pointer filter = BSplineFilterType::New();<br><br>    filter -&gt;SetSize( size );<br>    filter -&gt;SetOrigin( origin );<br>    filter -&gt;SetSpacing( spacing );<br>    filter -&gt;SetGenerateOutputImage ( false );<br>
    filter -&gt;SetNumberOfLevels ( 2 );<br>    filter -&gt;SetSplineOrder ( 3 );<br><br>    BSplineFilterType :: ArrayType ncps;<br>    ncps.SetElement(0,8);<br>    ncps.SetElement(1,8);<br>    ncps.SetElement(2,8);<br>    ncps.SetElement(3,5);<br>
    filter -&gt; SetNumberOfControlPoints ( ncps );<br>    <br>    filter -&gt;SetInput( fieldPoints );<br><br>    filter-&gt;SetDebug(true);<br><br>    try<br>    {<br>        filter -&gt;Update ();<br>    }<br>    catch ( itk::ExceptionObject &amp; excp )<br>
    {<br>        std :: cerr &lt;&lt; &quot;Test 2: itkBSplineScatteredDataImageFilter exception thrown&quot; &lt;&lt; std:: endl;<br>        std::cerr&lt;&lt; excp &lt;&lt;std::endl;<br>    }<br><br>That&#39;s the Spline definition and afterwards I wanna resample the BSpline at arbitrary positions:<br>
<br>float fPoint[3] = {0.0f,0.0f,0.0f};<br><br>    float sample_O[3] ={0.0f,0.0f,0.0f};<br>    for (int i=0; i&lt;3; i++)<br>        sample_O[i] = m_configuration.origin[i] - m_configuration.vs * (m_configuration.size[i] - 1) *0.5f;<br>
<br>    int sample_size[3];<br>    sample_size[0] = m_configuration.size[0];<br>    sample_size[1] = m_configuration.size[1];<br>    sample_size[2] = m_configuration.size[2];<br>    <br>    for (int iz=0, idx=0; iz&lt;sample_size[2]; iz++)<br>
    {<br>        fPoint[2] = sample_O[2] + iz*m_configuration.vs;<br>        for (int iy=0; iy&lt;sample_size[1]; iy++)<br>        {<br>            fPoint[1] = sample_O[1] + iy*m_configuration.vs;<br>            for (int ix=0; ix&lt;sample_size[0]; ix++)<br>
            {<br>                fPoint[0] = sample_O[0] + ix*m_configuration.vs;<br>                <br>                TimeVaryingDeformationFieldPointSetType :: PointType point;<br>                point [0] = fPoint[0];<br>
                point [1] = fPoint[1];<br>                point [2] = fPoint[2];<br><br>                VectorType V;<br>                filter-&gt;EvaluateAtPoint(point, V);<br>            }<br>            //cout &lt;&lt; iy &lt;&lt; endl;<br>
        }<br>        cout &lt;&lt; iz &lt;&lt; endl;<br>    }<br><br>Thans a lot for your help!<br><br>Best,<br><br>Kerstin<br><div class="gmail_quote">2012/1/10 Nicholas Tustison <span dir="ltr">&lt;<a href="mailto:ntustison@gmail.com">ntustison@gmail.com</a>&gt;</span><br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div style="word-wrap:break-word">That&#39;s fine.  I&#39;m assuming, though, that they reside<div>in a bounding box of some sort, correct?  Just </div>
<div>assign those points to their continuous coordinates</div><div>inside that bounding box.</div><div><div class="h5"><div><br></div><div><br><div><div><br></div><div><br><div><div>On Jan 10, 2012, at 5:31 AM, Kerstin Müller wrote:</div>
<br><blockquote type="cite">But my sample points are on continous world coordinates between for example -70.0 and 70.0.<br><br>Best,<br><br>Kerstin<br><br><div class="gmail_quote">2012/1/10 Nicholas Tustison <span dir="ltr">&lt;<a href="mailto:ntustison@gmail.com" target="_blank">ntustison@gmail.com</a>&gt;</span><br>

<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div style="word-wrap:break-word">Since your points reside in an image, an easy<div>parameterization is to assign their location to </div>

<div>their x,y,z location in the image.  Does that make</div><div>sense?<div><div><br><div><div><br><div><div>On Jan 10, 2012, at 5:19 AM, Kerstin Müller wrote:</div><br><blockquote type="cite"><br>Hi,<br><br>yes I also used your sample code, but still I got a heap corruption error. That&#39;s why I just wanted to try if my sample points are the problem.<br>

Do I need to scale them to a range between [0 1]?<br><br>
Best Regards,<br><br>Kerstin<br><br><div class="gmail_quote">2012/1/10 Nicholas Tustison <span dir="ltr">&lt;<a href="mailto:ntustison@gmail.com" target="_blank">ntustison@gmail.com</a>&gt;</span><br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">


<div style="word-wrap:break-word">Hi Kerstin,<div><br></div><div>Please keep the correspondence on the user&#39;s list.  </div><div><br></div><div>Why are you setting your vector to all zeros? </div><div><div><br>
</div><div><blockquote type="cite">            VectorType v;<br>            v[0] = 0.0;<br>            v[1] = 0.0;<br>            v[2] = 0.0;<br></blockquote><br></div></div><div>and why are you setting the location of your points to </div>


<div>(idx,idx,idx,t)?</div><div><div><br></div><div><blockquote type="cite">            point[0] = idx;<br>            point[1] = idx;<br>            point[2] = idx;<br>            point[3] = t;<br></blockquote>
<br></div></div><div>That can&#39;t possibly be right.  Did you look at the sample</div><div>code I gave you?</div><div><br></div><div>Nick</div><div><div><div><br></div><div><br></div><div><br><div><div>On Jan 9, 2012, at 2:04 AM, Kerstin Müller wrote:</div>


<br><blockquote type="cite">Hi,<br><br>thank you very much!<br>However I&#39;m still struggeling with a heap corruption problem. I&#39;m not sure if I set my parameters correct.<br><br>    typedef float RealType;<br>    <br>


    typedef itk::Vector&lt;RealType, 3&gt; VectorType;<br>
    typedef itk ::Image &lt;VectorType , 4 &gt; TimeVaryingDeformationFieldType;<br>    typedef itk :: PointSet &lt;VectorType , 4 &gt; TimeVaryingDeformationFieldPointSetType;<br><br>    TimeVaryingDeformationFieldPointSetType::Pointer fieldPoints = TimeVaryingDeformationFieldPointSetType:: New ();<br>



    fieldPoints-&gt;Initialize();<br><br>    // Set the sample points just for testing 10 linear points no movement<br>    for ( int t = 0; t &lt; m_configuration.num_projections ; t++ ) <br>    {<br>        int num = 0;<br>



        for( int idx = 0; idx &lt;= 10; idx++ )<br>        {<br>            <br>            VectorType v;<br>            v[0] = 0.0;<br>            v[1] = 0.0;<br>            v[2] = 0.0;<br><br>            <br>            TimeVaryingDeformationFieldPointSetType::PointType point;<br>



            point[0] = idx;<br>            point[1] = idx;<br>            point[2] = idx;<br>            point[3] = t;<br><br>            fieldPoints-&gt;SetPoint( num, point );<br>            fieldPoints-&gt;SetPointData(num, v);<br>



            num++;<br><br>        }<br>    }<br><br><br>    TimeVaryingDeformationFieldType::PointType        origin;<br>    TimeVaryingDeformationFieldType::SpacingType    spacing;<br>    TimeVaryingDeformationFieldType::SizeType        size;<br>



<br>    // Define the parametric domain<br>    spacing[0] = 1.0f;<br>    spacing[1] = 1.0f;<br>    spacing[2] = 1.0f;<br>    spacing[3] = 2.0f;<br><br>    size[0] = 128;<br>    size[1] = 128;<br>    size[2] = 128;<br>    size[3] = 66;<br>



    <br>    origin[0] = 0.0f;<br>    origin[1] = 0.0f;<br>    origin[2] = 0.0f;<br>    origin[3] = 0.0f;<br><br>    // Define the filter and set the parameters<br>    typedef itk::myBSplineScatteredDataPointSetToImageFilter&lt;TimeVaryingDeformationFieldPointSetType, TimeVaryingDeformationFieldType&gt; BSplineFilterType;<br>



    BSplineFilterType::Pointer filter = BSplineFilterType::New();<br><br>    filter -&gt;SetSize( size );<br>    filter -&gt;SetOrigin( origin );<br>    filter -&gt;SetSpacing( spacing );<br>    filter -&gt; SetGenerateOutputImage ( false );<br>



    filter -&gt; SetNumberOfLevels ( 2 );<br>    filter -&gt; SetSplineOrder ( 3 );<br><br>    <br>    BSplineFilterType :: ArrayType ncps;<br>    ncps.SetElement(0,8);<br>    ncps.SetElement(1,8);<br>    ncps.SetElement(2,8);<br>



    ncps.SetElement(3,34);<br>    filter -&gt; SetNumberOfControlPoints ( ncps );<br>    <br>    filter -&gt;SetInput( fieldPoints );<br><br>    filter-&gt;SetDebug(true);<br><br>    try<br>    {<br>        filter -&gt;Update ();<br>



    }<br>    catch ( itk::ExceptionObject &amp; excp )<br>    {<br>        std :: cerr &lt;&lt; &quot;Test 2: itkBSplineScatteredDataImageFilter exception thrown&quot; &lt;&lt; std:: endl;<br>        std::cerr&lt;&lt; excp &lt;&lt;std::endl;<br>



    }<br><br>Thank you very much!<br><br>Best,<br><br>Kerstin<br><div class="gmail_quote">2012/1/4 Nicholas Tustison <span dir="ltr">&lt;<a href="mailto:ntustison@gmail.com" target="_blank">ntustison@gmail.com</a>&gt;</span><br>




<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div style="word-wrap:break-word"><div>Okay, here&#39;s some pseudocode which hasn&#39;t been tested</div><div>since you have to provide it the interface to your values </div>




<div>but you should get the general idea.</div><div><br></div><div>  typedef float RealType;</div><div><div>  typedef itk::Image&lt;RealType, ImageDimension&gt; RealImageType;</div><div><br></div><div>  // Read in your reference image which will define the domain of </div>




<div>  // your B-spline displacement field.  Assume, for simplicity, that it</div><div>  // has identity direction</div><div><br></div><div>  typedef itk::ImageFileReader&lt;RealImageType&gt; ImageReaderType;</div><div>  typename ImageReaderType::Pointer reader = ImageReaderType::New();</div>




<div>  reader-&gt;SetFileName( XXXX );</div><div>  reader-&gt;Update();</div><div><br></div><div>  typedef itk::Vector&lt;RealType, ImageDimension&gt; VectorType;</div><div>  typedef itk::Image&lt;VectorType, ImageDimension+1&gt; TimeVaryingDeformationFieldType;</div>




<div><br></div><div>  typedef itk::PointSet&lt;VectorType, ImageDimension+1&gt; TimeVaryingDeformationFieldPointSetType;</div></div><div>  typename TimeVaryingDeformationFieldPointSetType::Pointer fieldPoints =</div><div>




    TimeVaryingDeformationFieldPointSetType::New();</div><div>  fieldPoints-&gt;Initialize();</div><div>  unsigned long count = 0;</div><div><br></div><div>  // Assume points are stored in a 2-D matrix called &#39;mySamplePoints&#39; where the column</div>




<div>  // is the time point and the row is the point</div><div><br></div><div>  for( unsigned int t = 0; t &lt; numberOfTimePoints; t++ )</div><div>    {</div><div>    for( unsigned int n = 0;  n &lt; numberOfPointsPerTimePoint; n++ )</div>




<div>      {</div><div>      VectorType displacement = mySamplePoints[n][t] - mySamplePoints[n][0];</div><div>      </div><div>      TimeVaryingDeformationFieldPointSetType::PointType parametric point;</div><div>      for( unsigned int d = 0; d &lt; ImageDimension; d++ )</div>




<div>        {</div><div>        point[d] = ( mySamplePoints[n][0] )[d];</div><div>        }</div><div>      point[ImageDimension] = t;</div><div><br></div><div>      fieldPoints-&gt;SetPoint( count, point );</div><div>      fieldPoints-&gt;SetPointData( count, displacement );</div>




<div>      count++;</div><div>      }</div><div>    }   </div><div><br></div><div><div>  TimeVaryingDeformationFieldType::PointType origin;</div></div><div><div>  TimeVaryingDeformationFieldType::SpacingType spacing;</div>




</div><div><div>  TimeVaryingDeformationFieldType::SizeType size;</div></div><div><br></div><div>   for( unsigned int d = 0; d &lt; ImageDimension; d++ )</div><div>      {</div><div>      origin[d] = reader-&gt;GetOutput()-&gt;GetOrigin()[d];</div>




<div>      size[d] = reader-&gt;GetOutput()-&gt;GetLargestPossibleRegion().GetSize()[d];</div><div>      spacing[d] = reader-&gt;GetOutput()-&gt;GetSpacing()[d];</div><div>      }</div><div>  // Now include the temporal information.  You can change this to whatever</div>




<div>  // resolution you like.  You just need to make sure that </div><div>  //  ( size[ImageDimension] - 1 ) * spacing[ImageDimension] = ( numberOfTimePoints - 1)</div><div>  // </div><div>  origin[ImageDimension] = 0;</div>




<div>  size[ImageDimension]  = numberOfTimePoints - 1;</div><div>  spacing[ImageDimension] = 1;</div><div><br></div><div><div>  typedef itk::BSplineScatteredDataPointSetToImageFilter</div><div>    &lt;TimeVaryingDeformationFieldPointSetType, TimeVaryingDeformationFieldType&gt; BSplineFilterType;</div>




<div>  typename BSplineFilterType::Pointer bspliner = BSplineFilterType::New();</div></div><div><div>  bspliner-&gt;SetOrigin( origin );</div><div>  bspliner-&gt;SetSpacing( spacing );</div><div>  bspliner-&gt;SetSize( size );</div>




<div>  bspliner-&gt;SetGenerateOutputImage( XXXX );</div><div>  bspliner-&gt;SetNumberOfLevels( XXXX );</div><div>  bspliner-&gt;SetSplineOrder( XXXX );</div><div>  bspliner-&gt;SetNumberOfControlPoints( XXXX );</div><div>




  bspliner-&gt;SetInput( fieldPoints );</div></div><div>  bspliner-&gt;Update();</div><div><div><div><br></div><div><br></div><div><br></div><br><div><div>On Jan 4, 2012, at 2:52 AM, Kerstin Müller wrote:</div>
<br><blockquote type="cite">Hi,<br><br>yes correctly.<br><br><div class="gmail_quote">2012/1/3 Nicholas Tustison <span dir="ltr">&lt;<a href="mailto:ntustison@gmail.com" target="_blank">ntustison@gmail.com</a>&gt;</span><br>




<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<div style="word-wrap:break-word"><div>Okay, I&#39;m assuming that these points correspond in </div><div>time, correct?  For example, point 937 in time point 0</div><div>corresponds to point 937 in time point 1, 2, 3, ...133,</div>





<div>right?  </div><div><div><div><br></div><div><br></div><br><div><div>On Jan 3, 2012, at 11:55 AM, Kerstin Müller wrote:</div><br><blockquote type="cite">Hi,<br><br>the basic problem I wanna solve:<br><br>I&#39;ll have scattered points in 3D which vary over time, they describe a motion of a surface over time. They are not ordered in x-y-and z-dimension.<br>





Now I want to represent that motion by BSplines in order to generate a dense motion vector field defined on specific 3D voxel positions.<br>
I&#39;ll have 133 timesteps and ~960 sample points in each time step. Now I want to fit the BSpline to it and resample it.<br><br>All the best,<br><br>Kerstin<br><br><div class="gmail_quote">2012/1/3 Nicholas Tustison <span dir="ltr">&lt;<a href="mailto:ntustison@gmail.com" target="_blank">ntustison@gmail.com</a>&gt;</span><br>






<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div style="word-wrap:break-word">In that case, let&#39;s start with the basic problem set-up.<div>Before, you described your problem as follows:<div>






<br><div><div><br><blockquote type="cite"><span style="border-collapse:separate;font-family:Helvetica;font-style:normal;font-variant:normal;font-weight:normal;letter-spacing:normal;line-height:normal;text-align:-webkit-auto;text-indent:0px;text-transform:none;white-space:normal;word-spacing:0px;font-size:medium"><blockquote type="cite">






<div><div>I&#39;m using the<span> </span><b>BSplineScatteredDataPointSetToImageFilter<span> </span></b>from the Insight Journal. However,<br>I cannot find the correct parameter in order to make the filter work and does not crash during the evaluation.<br>






I&#39;ll have scattered 3-D points over time and I want to fit a 4-D Bspline field to that points. Afterwards I want to evaluate the Bspline on a dense volume grid to one time step.</div></div></blockquote></span></blockquote>






</div></div></div></div><div><br></div><div>What do these scattered 3D+t points represent?  </div><div>Is it a curve, a time-varying scalar field, or something</div><div>else?</div><div><br></div></div></blockquote></div>






<br>
</blockquote></div><br></div></div></div></blockquote></div><br>
</blockquote></div><br></div></div></div></blockquote></div><br>
</blockquote></div><br></div></div></div></div></blockquote></div><br>
</blockquote></div><br></div></div></div></div></div></div></blockquote></div><br>
</blockquote></div><br></div></div></div></div></div></div></blockquote></div><br>