[Insight-users] itkDeformableMesh3DFilter : the violation access
Thomas Boettger
t.boettger at dkfz-heidelberg.de
Wed, 17 Mar 2004 12:03:13 +0100
Hi Julien,
at least I can answer a few of your questions.
I played around alot with this filter and also got bad results.
GetPoint() and GetPointData() are defined in itk::PointSet (the
superclass of mesh).
The stiffness parameters did not change my results either.
Julien Mercenier wrote:
> Hi Luis,
>
> this e-mail will be quite long....
> I haven't yet tried the exemple2 (registarion model2) you created with
> my example because I must format my PC...
>
> (1) First, I found how, when using the itkDeformableMesh3DFilter with my own initialization, not to have a crash.
>
> When initializing my own triangle mesh, I must associate (here by giving it a value) a
> data to each triangle cell.
> I found this because this was the main difference between my initializing and a meshSource (binary or sphere).
> I put an simple example in attachment :
> This is a tetraedra composed of triangle cells.
> The image (100*100*100) is a sphere (I've tried with a radius of 40.0) created by the user-guide
> example. Sorry, I don't push the shere.img in attachment (otherwise the e-mail will be refused).
>
> But, even this example doesn't work well :
> (a) variating the stiffness doesn't seem to variate the results (the nodes' displacements)
> (b) if I increase the cell data, the displacements seem to decrease.
> !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
> !!!!!!!!!!! I've tried with a data cell value = 100 and the program crashed !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
> !!!!!!!!!!! with the same access violation than before !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
> Is this the consequence of itkDeformableMesh3DFilter.txx/line 279 ?
> while (celldata != myCellData->End()){
> x = celldata.Value();
> m_K[j] = m_StiffnessMatrix+((int) x);
> while casting the double to an int ???????????????????
> (c) when I increase the force scale (that must be big to obtain a deformation) or the time step,
> the deformations seem to increase
> (d) the step number doesn' t seem to influence the results
What step number do you mean?
> (e) the step threshold have effects
> Although these, the results are bad !... : )
> I didn't manage to find good parameters to this example.
>
> (2) Can you tell me on which article you based your work to construct this filter...?
> This could be very helpful to me to understand the stiffnessmatrix computing and maybe other things in the code.
> I must understand this one or create another one (this is a part of my "ingeneering thesis"...I don't know if this is called like that
> in English).
Ask the programmer (Ting Chen chenting at graphics.cis.upenn.edu) on what
article it is based .
>
> (3) Now, I have questions and remarks about the source code itself...(itkDeformableMesh3DFilter.txx)
>
> (a) itkDeformableMesh3DFilter.txx/line 702
> InputPointDataContainerPointer myForceData = m_Forces->GetPointData();
> where is defined (in the documentation) the function GetPointData ? I didn't find it in the itk::mesh class.
itk::PointSet
>
> (b) where is called the function
> DeformableMesh3DFilter<TInputMesh, TOutputMesh>
> ::SetStiffnessMatrix ( vnl_matrix_fixed<double, 4, 4> *stiff, int i )
> defined at itkDeformableMesh3DFilter.txx/line 238 ??
>
> (c) I don't understand this command
> x = celldata.Value();
> m_K[j] = m_StiffnessMatrix+((int) x);
> defined at itkDeformableMesh3DFilter.txx/line 280
this is pointer arithmetics where a certain stiffness matrix depending
on the celldata.Value() (which provides the index) is assigned to each cell
>
> and why do you use data associated with cell ?
>
> (d) itkDeformableMesh3DFilter.txx/line 319
> m_Displacements->GetPoint (tp[0], v1_pt);
> m_Displacements->GetPoint (tp[1], v2_pt);
> m_Displacements->GetPoint (tp[2], v3_pt);
>
> where is defined the function GetPoint (same reason as in (a)) ???
itk::PointSet
>
> (e) itkDeformableMesh3DFilter.txx/line 489
> output->SetCellData(i, (PixelType)x);
>
> why do you modify the cell data at the output ?
>
> (f) itkDeformableMesh3DFilter.txx/line 551
>
> coord[0] = (int) vec_p[0];
> coord[1] = (int) vec_p[1];
> coord[2] = (int) vec_p[2];
>
> are you sure there is no problem if the vec_p[0] type is double ???
>
> (g) itkDeformableMesh3DFilter.txx/line 555
>
> if ( m_Potential->GetPixel(coord) != m_ObjectLabel )
>
> what is m_ObjectLabel ?
>
> (h)itkDeformableMesh3DFilter.txx/line 521
>
> template <typename TInputMesh, typename TOutputMesh>
> void
> DeformableMesh3DFilter<TInputMesh, TOutputMesh>
> ::PotentialFit()
>
> what is the use of this function ? I don't see it is called anywhere.
>
> (i) itkDeformableMesh3DFilter.txx/line 819
>
> coa = -(v1[1]*(v2[2]-v3[2]) +
> v2[1]*(v3[2]-v1[2]) +
> v3[1]*(v1[2]-v2[2])) ;
> cob = -(v1[2] * (v2[0]-v3[0]) +
> v2[2]*(v3[0]-v1[0]) +
> v3[2]*(v1[0]-v2[0])) ;
> coc = -(v1[0] * (v2[1]-v3[1]) +
> v2[0]*(v3[1]-v1[1]) +
> v3[0]*(v1[1]-v2[1])) ;
>
> I don't understand these steps in the normals computing. What is it based on ?
>
> (j) itkDeformableMesh3DFilter.txx/line 783
>
> Remark : in this function, you compute several times the normal at one point.
> The value of the normal saved is the last one computed. The other computings are lost.
> Is it right ?
> So the normal at one node seems to be the consequence of only one cell (triangle)....
I think the normals are computed for each triangle...
>
>
> (k) itkDeformableMesh3DFilter.txx/line 749
>
> vec_for[0] = vec_for[0] + (vec_loc[0]-coord[0])*tmp_vec_1[0]
> + (vec_loc[1]-coord[1])*tmp_vec_2[0] + (vec_loc[2]-coord[2])*tmp_vec_3[0];
> vec_for[1] = vec_for[1] + (vec_loc[1]-coord[1])*tmp_vec_2[1]
> + (vec_loc[0]-coord[0])*tmp_vec_1[1] + (vec_loc[2]-coord[2])*tmp_vec_3[1];
> vec_for[2] = vec_for[2] + (vec_loc[2]-coord[2])*tmp_vec_3[2]
> + (vec_loc[1]-coord[1])*tmp_vec_2[2] + (vec_loc[0]-coord[0])*tmp_vec_1[2];
>
> remark : here, in the GradientFit function, using generalized Taylor theorem, I think you
> assumed the spacing to be equal to 1 in each direction.
> The voxel spacing can although be strongly different.
>
> (l) The last but not least.
> itkDeformableMesh3DFilter.txx/line 761
>
>
> mag = vec_for[0]*vec_nor[0] + vec_for[1]*vec_nor[1]+ vec_for[2]*vec_nor[2];
>
> vec_for[0] = m_GradientMagnitude*mag*vec_nor[0]/*num_for*/;
> vec_for[1] = m_GradientMagnitude*mag*vec_nor[1]/*num_for*/;
> vec_for[2] = m_GradientMagnitude*mag*vec_nor[2]/*num_for*/;
>
> mag = sqrt (vec_for[0]*vec_for[0] + vec_for[1]*vec_for[1]+ vec_for[2]*vec_for[2]);
> if (mag > 0.5)
> for (int i=0; i<3; i++) vec_for[i] = (0.5 * vec_for[i])/mag;
> forces.Value() = vec_for;
>
> You compute the force vector on one point by projection of the gradient vector fitted
> on the normal vector of this point.
> So, the force direction is defined by the normal and not by the gradient field.
> Am I right ?
>
> So, considering I'm right, using GVF field in place of the simple gradient field
> would not conclude to better results using this itkDeformableMesh3DFilter...
> This GVF filed would although be so helpful to recover concave regions....
>
> What do you think ?
>
>
>
> I hope you will not be angry against me for writing such a long email, but this is
> very important for my thesis...
>
> Thanks in advance,
>
> Julien,BELGIUM
--
Dipl.-Inform. Thomas Boettger
Deutsches Krebsforschungszentrum (German Cancer Research Center)
Div. Medical and Biological Informatics B010 Tel: (+49) 6221-42 2328
Im Neuenheimer Feld 280 Fax: (+49) 6221-42 2345
D-69120 Heidelberg e-mail: t.boettger at dkfz.de
Germany http://www.dkfz.de/mbi/people/thomasb.shtml