[vtkusers] vtkTransform problem: really roundoff error?

Wes Turner wes.turner at kitware.com
Mon Feb 8 11:48:52 EST 2010


Michele,

Matrix inversion is very sensitive, and the successive inversion operations
can easily result in the accumulation of round-off error.  You matrix does
not seem to be badly ill-conditioned, which would worsen the issue --
admittedly I just eyeballed it instead of calculating it -- and yet after 2
million matrix inversions, I find it believable that you have lost all of
your significant digits.  Why would you do this?  If you need to use both
the inverse and forward transforms, make a copy and invert it once.  You
could also use the GetInverse() function to get the 4x4 matrix that
implements the inverse of the current transform.  It will always be
calculated from the original matrix and will not exhibit this problem.
 Constantly going back and forth between the forward and inverse transforms
is a bad idea from both a mathematical standpoint and from a processing
standpoint.

- Wes

On Mon, Feb 8, 2010 at 10:45 AM, Michele Conconi
<michele.conconi at unibo.it>wrote:

> Hi again.
> I'm sorry, I probably did not explain my problem in the best way.
> I do not think that it is a round off error. Indeed, what I'm doing is not
> reading the output after a bing number of inversion. Instead, I set the
> transformation matrix to some values, invert it many times, than I reset the
> matrix into my check configuration and see if the output still corresponds
> with what I impose. So the problem seems to be that by inverting the
> vktTransofrm enough times, it "obeys" no more, i.e. it goes no more where it
> should.
>
> As both David and Murat suggested, I also thought at the beginning that the
> problem was caused by some roundoff error. So, to compensateit, I  reset the
> matrix of my transformation to the identity before any change. I also
> performed my checks introducing some tolerance as you may see in the new
> code below. Unfortunately, with epsilon up to 1.5 I  obtain the same
> problem. The filter change its behavior always after the same number of
> cycle (214 for your knowledge). Also, the error I obtain is not of the
> dimension of a roundoff error. Indeed, at the failure the center of the
> polydata is
>
> x = -0.0438237 , y =0.0876474 , z = -4.98359  instead of  x = 9.95618,  y =
> -23.462, z = 26.4141
>
> and the matrix is the identity instead of
>
> 1   0   0   10
> 0    0.707107   0.707107   -20
> 0   -0.707107   0.707107   30
> 0   0   0   1
>
> Where could be the problem?
> Michele
>
> #include "vtkTransformPolyDataFilter.h"
> #include "vtkSTLReader.h"
> #include "vtkTransform.h"
> #include <iostream>
>
>
> using namespace std;
>
>
> int main(int argc, char *argv[])
> {
>   double center[]={0,0,0};
>   double c[]={0,0,0};
>
>   vtkSTLReader *volumeReader = vtkSTLReader::New();
>   vtkTransformPolyDataFilter *TF = vtkTransformPolyDataFilter::New();
>   vtkTransform *ActualPosition = vtkTransform::New();
>   vtkPolyData *poly = vtkPolyData::New();
>
>   volumeReader->SetFileName("cube.stl");
>   TF->SetInput(volumeReader->GetOutput());
>   volumeReader->Delete();
>
>   TF->SetTransform(ActualPosition);
>   poly=TF->GetOutput();
>
>   poly->GetCenter(center);
>   cout<<center[0]<<" "<<center[1]<<" "<<center[2]<<endl;//the center of the
> cube containing the polydata in its ZERO configuration, i.e. no rotation or
> translation
>
>   ActualPosition->Identity();
>   ActualPosition->PostMultiply();
>   ActualPosition->RotateX(-45);
>   ActualPosition->Translate(10,-20,30);
>   poly->Update();
>
>   poly->GetCenter(center);
>   cout<<center[0]<<" "<<center[1]<<" "<<center[2]<<endl;//the center of the
> cube containing the polydata in the control configuration  configuration
>   ActualPosition->Print(cout);
>     int i=0;
>   while(i<200000)
>   {
>       //modify the vtkTransform
>       ActualPosition->Identity();
>       ActualPosition->PostMultiply();
>       ActualPosition->RotateX(0);
>       ActualPosition->Translate(0,0,0);
>       poly->Update();
>
>       for (int j=0; j<1e7; j++)
>       {
>           ActualPosition->Inverse();
>           //my operations...
>           ActualPosition->Inverse();
>       }
>
>
>       if(i%10000==0)
>           cout<<i<<endl;
>
>       //bring it back to the control configuration
>       ActualPosition->Identity();
>       ActualPosition->PostMultiply();
>       ActualPosition->RotateX(-45);
>       ActualPosition->Translate(10,-20,30);
>       poly->Update();
>
>
>       poly->GetCenter(c);
>       //check if the center of the cube containing the polydata is still in
> the right position
>       double epsilon=1.5;
>       if(abs(c[0]-center[0])<epsilon && abs(c[1]-center[1])<epsilon &&
> abs(c[2]-center[2])<epsilon)
>           i++;
>       else
>       {
>           cout<<c[0]<<" "<<c[1]<<" "<<c[2]<<endl;
>           cout<<i<<endl;
>           ActualPosition->Print(cout);
>           i=2000000;
>       }
>   }
>
>   ActualPosition->Delete();
>   TF->Delete();
> }
> --
> * Eng. Michele Conconi*
> Ph.D. student
> DIEM - Dept. of Mechanical Engineering
> Alma Mater Studiorum - University of Bologna
> Viale Risorgimento 2, 40136, Bologna, Italy
> Email: michele.conconi at unibo.it
> Website: http://www.diem.ing.unibo.it/grab
> Office: (+39) 051 20 93451
> Mobile: (+39) 329 0287996
>
>
>
> * INFORMAZIONE RISERVATA E CONFIDENZIALE *
> Questo messaggio contiene informazioni riservate e confidenziali.
> Se il lettore non fosse il destinatario del messaggio, inoltrato e ricevuto
> per errore,
> il testo dovrà essere immediatamente cancellato dal computer del ricevente.
> E' assolutamente proibita qualunque circolazione, disseminazione o
> copia del messaggio spedito e/o ricevuto per errore.
>
> * CONFIDENTIALITY and WARNING NOTICE *
> This message may contain legally privileged or confidential information.
> If the reader is not the intended recipient, you received this message in
> error
> and it should therefore be deleted from your computer at once.
> Any dissemination, distribution and copying of a message
> mistakenly sent or received is strictly forbidden.
> _______________________________________________
> Powered by www.kitware.com
>
> Visit other Kitware open-source projects at
> http://www.kitware.com/opensource/opensource.html
>
> Please keep messages on-topic and check the VTK FAQ at:
> http://www.vtk.org/Wiki/VTK_FAQ
>
> Follow this link to subscribe/unsubscribe:
> http://www.vtk.org/mailman/listinfo/vtkusers
>



-- 
Wesley D. Turner, Ph.D.
Kitware, Inc.
Technical Leader
28 Corporate Drive
Clifton Park, NY 12065-8662
Phone: 518-881-4920
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://www.vtk.org/pipermail/vtkusers/attachments/20100208/6d0b0625/attachment.htm>


More information about the vtkusers mailing list