1 /*# ---------------------------------------------------------------------
3 # Copyright (c) CREATIS (Centre de Recherche en Acquisition et Traitement de l'Image
5 # Authors : Eduardo Davila, Frederic Cervenansky, Claire Mouton
6 # Previous Authors : Laurent Guigues, Jean-Pierre Roux
7 # CreaTools website : www.creatis.insa-lyon.fr/site/fr/creatools_accueil
9 # This software is governed by the CeCILL-B license under French law and
10 # abiding by the rules of distribution of free software. You can use,
11 # modify and/ or redistribute the software under the terms of the CeCILL-B
12 # license as circulated by CEA, CNRS and INRIA at the following URL
13 # http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
14 # or in the file LICENSE.txt.
16 # As a counterpart to the access to the source code and rights to copy,
17 # modify and redistribute granted by the license, users are provided only
18 # with a limited warranty and the software's author, the holder of the
19 # economic rights, and the successive licensors have only limited
22 # The fact that you are presently reading this means that you have had
23 # knowledge of the CeCILL-B license and that you accept its terms.
24 # ------------------------------------------------------------------------ */
26 /*=========================================================================
29 Module: $RCSfile: marAxis.cpp,v $
31 Date: $Date: 2012/11/15 14:15:31 $
32 Version: $Revision: 1.2 $
34 Copyright: (c) 2002, 2003
37 This software is distributed WITHOUT ANY WARRANTY; without even
38 the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
39 PURPOSE. See the above copyright notice for more information.
41 =========================================================================*/
42 #include <vtkPlaneSource.h>
43 #include <vtkKochanekSpline.h>
44 #include <vtkPolyDataConnectivityFilter.h>
45 #include <vtkContourFilter.h>
46 #include <vtkProbeFilter.h>
49 #include <vtkImageChangeInformation.h>
50 #include <vtkImageShiftScale.h>
51 #include <vtkDoubleArray.h>
55 #include "marVector.h"
57 #include "ExtractionAxe.h"
58 #include "FonctionsGrales.h"
60 #include <vtkCellArray.h>
62 #include <vtkCleanPolyData.h>
63 #include <vtkIdList.h>
65 #include <vtkMatrix4x4.h>
66 #include <vtkPointData.h>
67 #include <vtkPolyData.h>
69 #include <vtkSetGet.h>
70 #include <vtkStripper.h>
71 #include <vtkStructuredPoints.h>
72 #include <vtkTransform.h>
74 // ----------------------------------------------------------------------------
75 marAxis::marAxis( marParameters* p ) :
77 kCurve( 3, INDX_count - 3 ),
78 _healthySlice( -1 ), _startQuant( -1 ),
79 _finishQuant( -1 ), _actualQuant( -1 )
85 _referenceAverDiam=-1;
90 _healthySliceStart=-1;
96 // ----------------------------------------------------------------------------
97 void marAxis::addAxisPoint( double* p )
99 double pt[ INDX_count ];
101 memcpy( pt, p, sizeof( POINTAXE ) );
102 pt[ INDX_SIGNALVALUE ] = 0;
103 addControlPoint( pt, pt + 3 );
106 // ----------------------------------------------------------------------------
108 double* marAxis::getNormal( unsigned int slice ) {
110 double* ret = new double[ 3 ];
113 ret[0] = _points[ slice + 1 ][0] - _points[ slice ][0];
114 ret[1] = _points[ slice + 1 ][1] - _points[ slice ][1];
115 ret[2] = _points[ slice + 1 ][2] - _points[ slice ][2];
117 else if ( slice == _points.size( ) - 1 ) {
118 ret[0] = _points[ slice - 1 ][0] - _points[ slice ][0];
119 ret[1] = _points[ slice - 1 ][1] - _points[ slice ][1];
120 ret[2] = _points[ slice - 1 ][2] - _points[ slice ][2];
122 else if ( slice > 0 && slice < _points.size( ) - 1 ) {
123 ret[0] = _points[ slice - 1 ][0] - _points[ slice + 1 ][0];
124 ret[1] = _points[ slice - 1 ][1] - _points[ slice + 1 ][1];
125 ret[2] = _points[ slice - 1 ][2] - _points[ slice + 1 ][2];
131 // ----------------------------------------------------------------------------
132 void marAxis::changeAxisResolution( )
135 int nOut = ( int )ceil(
137 getDoubleParam( marParameters::e_axis_discret_step )
138 * ( double )( GetNumberOfPoints( ) )
141 i = GetNumberOfPoints( );
142 double d = getParameters( )->
143 getDoubleParam( marParameters::e_axis_discret_step );
145 SetResolution( nOut );
146 for( i = 0; i < _contours.size( ); i++ )
147 delete _contours[ i ];
149 for( i = 0; i < nOut; i++ )
150 _contours.push_back( new marContour( getParameters( ) ) );
154 // ----------------------------------------------------------------------------
155 void marAxis::calculateSignal( kVolume* vol )
158 * PSIGNAL_USHORT <-> unsigned short *
160 /* PSIGNAL_FLOAT PutSignalValueFloat(PSIGNAL_FLOAT signal, int pos, float value)
169 unsigned int i, mask_size = getParameters( )-> getIntParam( marParameters::e_mask_size );
171 while( _signal.size( ) > 0 ) _signal.pop_back( );
174 PSIGNAL_USHORT signal = ( PSIGNAL_USHORT ) IdSigAlloc( _points.size( ), SIG_USHORT );
175 for( i = 0; i < _points.size( ); i++ ) {
181 _points_disc, getNumberOfControlPoints( ));
183 maxlocal = vol->GetMaxIntSphere2( p, rayon );
184 PutSignalValue( signal, i, maxlocal );
187 PSIGNAL_FLOAT mask_signal = ( PSIGNAL_FLOAT )IdSigAlloc( mask_size, SIG_FLOAT );
189 for( i = 0; i < mask_size; i++ )
190 PutSignalValueFloat( mask_signal, i, 1 );
192 PSIGNAL_USHORT signalfilt = ( PSIGNAL_USHORT )IdSigAlloc( _points.size( ), SIG_USHORT );
193 signalfilt = MonSigConvolve( signal, mask_signal, signalfilt, 1.0 / ( double )mask_size, 0 );
194 for( i = 0; i < _points.size( ); i++ ){
195 _signal.push_back( GetSignalValue( signalfilt, i ) );
198 IdSigFree( mask_signal );
199 IdSigFree( signalfilt );
202 // ----------------------------------------------------------------------------
203 void marAxis::start( )
206 int swap, n = GetNumberOfPoints( );
208 swap = GSL_MIN( _startQuant, _finishQuant );
209 _finishQuant = GSL_MAX( _startQuant, _finishQuant );
212 _startQuant = ( _startQuant >= 0 )? _startQuant: 0;
213 _finishQuant = ( _finishQuant >= n )? _finishQuant: n - 1;
214 _actualQuant = _startQuant;
218 // ----------------------------------------------------------------------------
219 void marAxis::cut( int slice, bool up )
224 // ----------------------------------------------------------------------------
225 void marAxis::doSpline( )
228 unsigned int nop, nip;
230 float p1[ 3 ], p2[ 3 ];
231 double p[marAxis::INDX_count];
233 vtkKochanekSpline* spX = vtkKochanekSpline::New( );
234 vtkKochanekSpline* spY = vtkKochanekSpline::New( );
235 vtkKochanekSpline* spZ = vtkKochanekSpline::New( );
237 for( i = 0, length = 0.0; i < _controlPoints.size( ); i++ ) {
238 getControlPoint(i, p, p+3);
239 p2[ 0 ] = (float) p[marAxis::INDX_X];
240 p2[ 1 ] = (float) p[marAxis::INDX_Y];
241 p2[ 2 ] = (float) p[marAxis::INDX_Z];
242 spX->AddPoint( (float) i, p2[ 0 ] );
243 spY->AddPoint( (float) i, p2[ 1 ] );
244 spZ->AddPoint( (float) i, p2[ 2 ] );
246 length += (float) (sqrt( (p2[0]-p1[0])*(p2[0]-p1[0]) +
247 (p2[1]-p1[1])*(p2[1]-p1[1]) +
248 (p2[2]-p1[2])*(p2[2]-p1[2]) ) );
249 p1[ 0 ] = p2[ 0 ]; p1[ 1 ] = p2[ 1 ]; p1[ 2 ] = p2[ 2 ];
252 nop = ( unsigned int ) ( getParameters( )->getDoubleParam(
253 marParameters::e_axis_discret_step ) * length );
254 nip = _controlPoints.size( );
257 while( _points.size( ) > 0 ) {
258 delete _points[ _points.size( ) - 1 ];
262 while( _contours.size( ) > 0 ) {
263 delete _contours[ _contours.size( ) - 1 ];
264 _contours.pop_back( );
267 for( i = 0; i < nop; i++ ) {
268 t = (float) (( ( double ) nip - 1.0 ) / ( ( double ) nop - 1.0 ) * i);
269 double* np = new double[ 3 ];
273 np[ 0 ] = (double) (spX->Evaluate( t ));
274 np[ 1 ] = (double) (spY->Evaluate( t ));
275 np[ 2 ] = (double) (spZ->Evaluate( t ));
276 _points.push_back( np );
277 _contours.push_back( new marContour( i, getParameters( ) ) );
289 // ----------------------------------------------------------------------------
290 void marAxis::sliceVolumeAxis( kVolume* vol, bool forceCnt )
293 int nCnts = ( int ) getNumberOfSplinePoints( );
294 int sizeIma = getParameters( )->getSizeIma( );
295 double dimIma = getParameters( )->getDimIma( );
296 double voxSize = getParameters( )->getVoxelSize( );
299 for( unsigned int i = 0; i < _contours.size( ); i++ )
300 delete _contours[ i ];
303 for( unsigned int i = 0; i < _slices.size( ); i++ ) {
305 _3Dslices[ i ]->Delete( );
306 _quantificationImages[ i ]->Delete( );
311 _quantificationImages.clear( );
314 // "Cutter" initialization
315 vtkStructuredPoints* stPoints = vtkStructuredPoints::New( );
318 for( int k = 0; k < nCnts; k++ ) {
319 //s = ( double )k / ( double )( nCnts - 1 );
322 vtkPlaneSource* pSource = vtkPlaneSource::New( );
323 pSource->SetOrigin( p[ 0 ], p[ 1 ], p[ 2 ] );
324 pSource->SetPoint1( p[ 0 ] + dimIma - 1.0, p[ 1 ], p[ 2 ] );
325 pSource->SetPoint2( p[ 0 ], p[ 1 ], p[ 2 ] + dimIma - 1.0 );
326 pSource->SetResolution( sizeIma - 1, sizeIma - 1 );
328 pSource->SetCenter( p[ 0 ], p[ 1 ], p[ 2 ] );
329 pSource->SetNormal( n[ 0 ], n[ 1 ], n[ 2 ] );
332 _3Dslices.push_back( vtkProbeFilter::New( ) );
333 _3Dslices[ k ]->SetInput( ( vtkDataSet* )pSource->GetOutput( ) );
334 _3Dslices[ k ]->SetSource( vol->castVtk( ) );
335 _3Dslices[ k ]->Update( );
338 stPoints->GetPointData( )->
339 SetScalars( _3Dslices[ k ]->GetOutput( )->
340 GetPointData( )->GetScalars( ) );
341 stPoints->SetDimensions( sizeIma, sizeIma, 1 );
342 stPoints->SetScalarType( vol->castVtk( )->GetDataObjectType( ) );
346 marContour* tmp_marContour = new marContour( k, getParameters( ) ) ;
347 tmp_marContour->calculateVariables();
348 _contours.push_back( tmp_marContour );
350 _slices.push_back( new kVolume( kVolume::UCHAR, 1, 1, 1 ) );
352 *( _slices[ k ] ) = ( vtkImageData* )stPoints; //copy
353 _quantificationImages.push_back( (_slices[ k ])->castVtk( ) );
362 // ----------------------------------------------------------------------------
363 int marAxis::getNumberOfContours( )
365 return( _contours.size( ) );
368 // ----------------------------------------------------------------------------
369 int marAxis::getNumberOfSplinePoints( )
371 return( _points.size( ) );
375 // EED Esto toca limpiarlo o reimplementarlo en varias partes o algo
376 // ----------------------------------------------------------------------------
377 vtkPolyData* marAxis::setContour( int i, int x, int y, std::vector< double* >* points, marContour* c )
381 c->setS( c->getS( ) );
382 _contours[ i ]->copyFrom( *c );
386 // PS -> // int numberOfPoints, numberOfCells; // PS
387 // PS -> // int id; // PS
388 // PS -> // gslobj_vector p( 3 ); // PS
389 // PS -> // double new_p[3]; // PS
390 // PS -> // float pf[ 3 ]; // PS
391 vtkImageData* imagedata;
393 // Contour parameters
394 // PS -> int typ = getParameters( )->getIntParam( marParameters::e_algorithm_type );
395 double thr = getParameters( )->getDoubleParam( marParameters::e_threshold_isocontour );
397 // PS -> double thr = ( typ == marParameters::ISOCONTOURS )?
398 // PS -> getParameters( )->
399 // PS -> getDoubleParam( marParameters::e_threshold_isocontour ):
400 // PS -> ( typ == marParameters::SNAKE ) ?
401 // PS -> getParameters( )->
402 // PS -> getDoubleParam( marParameters::e_threshold_snake_isocontour ):
406 int sizeIma = getParameters( )->getSizeIma( );
407 double localint = getSignal( i );
408 double sigma = getParameters( )->getDoubleParam( marParameters::e_sigma );
409 double threshold = localint * thr / 100.0;
410 int dimIma = ( int ) ( ( double ) sizeIma /getParameters( )->getDoubleParam( marParameters::e_scale ) );
411 int vmin = ( int )threshold;
412 int vmax = ( int )threshold;
414 x = ( x != -1 )? x: sizeIma / 2;
415 y = ( y != -1 )? y: sizeIma / 2;
418 imagedata = (vtkImageData*) ((getSlice( i , NULL ))->castVtk( ));
419 float *range = imagedata->GetScalarRange();
421 //vtkKitwareContourFilter* cntVTK = vtkKitwareContourFilter::New( );
422 vtkContourFilter* cntVTK = vtkContourFilter::New( );
423 cntVTK->SetInput( imagedata );
424 cntVTK->SetNumberOfContours( 1 );
425 //cntVTK->SetValue( 0, vmin );
426 cntVTK->SetValue( 0, (range[1]*thr/100) );
427 //cntVTK->SetValue( 1, vmax );
430 vtkCleanPolyData* cpd = vtkCleanPolyData::New( );
431 cpd->SetInput( cntVTK->GetOutput( ) );
432 cpd->ConvertLinesToPointsOff( );
435 vtkPolyDataConnectivityFilter* conn = vtkPolyDataConnectivityFilter::New( );
436 //conn->SetMaxRecursionDepth( 3000 );
438 // PS -> //conn->SetInput( cntVTK->GetOutput( ) ); PS
439 conn->SetInput( cpd->GetOutput( ) );
441 conn->SetClosestPoint( x, y, 0 );
442 conn->SetExtractionModeToClosestPointRegion( );
445 vtkCleanPolyData* cpd2 = vtkCleanPolyData::New( );
446 cpd2->SetInput( conn->GetOutput( ) );
449 vtkPolyData* polyDataResult = cpd2->GetOutput() ;
451 double ta = _contours[ i ]->getS( );
453 _contours[ i ]->reset( );
455 // PS -> gslobj_vector o( 3 ), n( 3 ), y( 3 ), vr( 3 );
459 int numberOfPoints,numberOfPoints2;
463 // PS -> double new_p[3];
466 numberOfCells = polyDataResult->GetNumberOfCells( );
467 cell = polyDataResult->GetCell( 0 );
468 pids = cell->GetPointIds( );
469 numberOfPoints = pids->GetNumberOfIds( );
471 numberOfPoints2=polyDataResult->GetNumberOfPoints();
473 marVector o(3),n(3),y(3),vr(3);
475 y( 0 ) = y( 2 ) = 0.0;
477 getPoint( ( double* )o, ta );
478 getTangent( ( double* )n, ta );
479 vr = y.cross( n ).normalize( );
480 sint = y.cross( n ).norm2( );
483 for( int j = 0; j < numberOfPoints; j++ )
485 id = pids->GetId( j );
486 polyDataResult->GetPoint( id, pf );
492 p = ( p * cost ) + ( vr * ( vr.dot( p ) * ( 1.0 - cost ) ) ) -
493 ( p.cross( vr ) * sint ) + o;
495 // PS -> double p_cross_vr[3];
496 // PS -> vtkMath::Cross( p, vr, p_cross_vr);
497 // PS -> for(int k=0; k<3; k++)
499 // PS -> new_p[i] = (new_p[k] * cost) + ( vr(k) * ( vtkMath::Dot(vr, p) * ( 1.0 - cost ) ) ) -
500 // PS -> ( p_cross_vr[k] * sint ) + o(i);
503 //TODO: check if p and new_p are the same...
504 _contours[ i ]->addContourPoint( p );
507 // PS -> int nbPoints=polyDataResult->GetNumberOfPoints();
509 // PS -> float* vtkContourPoint;
510 // PS -> vtkContourPoint= new float[3];
511 // PS -> double* marContourPoint;
512 // PS -> marContourPoint= new double[3];
514 // PS -> for (int j=0;j< nbPoints;j++)
516 // PS -> vtkContourPoint=polyDataResult->GetPoint(j);
517 // PS -> marContourPoint[0]=(double)vtkContourPoint[0];
518 // PS -> marContourPoint[1]=(double)vtkContourPoint[1];
519 // PS -> marContourPoint[2]=(double)vtkContourPoint[2];
520 // PS -> _contours[ i ]->addContourPoint( marContourPoint );
522 marContour * tutu=_contours[ i ];
523 _contours[ i ]->getArea();
524 _contours[ i ]->getPerimeter();
525 return( polyDataResult );
526 //return( cpd->GetOutput( ) );
527 //return( cntVTK->GetOutput( ) );
530 //EED1 vtkStripper* isoStrips = vtkStripper::New( );
531 //EED1 isoStrips->SetInput( cpd2->GetOutput( ) );
532 //EED1 isoStrips->Update( );
534 //EED1 vtkPolyData* data = isoStrips->GetOutput( );
535 //EED1 data->Update( );
536 //EED1 vtkCell* cell;
537 //EED1 vtkIdList* pids;
539 //EED1 numberOfCells = data->GetNumberOfCells( );
540 //EED1 cell = data->GetCell( 0 );
541 //EED1 pids = cell->GetPointIds( );
542 //EED1 numberOfPoints = pids->GetNumberOfIds( );
545 //EED2 _contours[ i ]->reset( );
546 //EED2 if( typ == marParameters::ISOCONTOURS ) {
547 //EED2 gslobj_vector o( 3 ), n( 3 ), y( 3 ), vr( 3 );
548 //EED2 double cost, sint;
549 //EED2 y( 0 ) = y( 2 ) = 0.0;
551 //EED2 getPoint( ( double* )o, ta );
552 //EED2 getTangent( ( double* )n, ta );
553 //EED2 vr = y.cross( n ).normalize( );
554 //EED2 sint = y.cross( n ).norm2( );
555 //EED2 cost = y.dot( n );
557 //EED2 for( int j = 0; j < numberOfPoints; j++ ) {
558 //EED2 id = pids->GetId( j );
559 //EED2 data->GetPoint( id, pf );
560 //EED2 p( 0 ) = pf[ 0 ];
562 //EED2 p( 2 ) = pf[ 1 ];
564 //EED2 // 3D transformation
565 //EED2 p = ( p * cost ) + ( vr * ( vr.dot( p ) * ( 1.0 - cost ) ) ) -
566 //EED2 ( p.cross( vr ) * sint ) + o;
568 //EED2 double p_cross_vr[3];
569 //EED2 vtkMath::Cross( p, vr, p_cross_vr);
570 //EED2 for(int ii=0; ii<3; ii++) {
571 //EED2 new_p[i] = (new_p[ii] * cost) + ( vr(ii) * ( vtkMath::Dot(vr, p) * ( 1.0 - cost ) ) ) -
572 //EED2 ( p_cross_vr[ii] * sint ) + o(i);
575 //EED2 //TODO: check if p and new_p are the same...
576 //EED2 _contours[ i ]->addContourPoint( p );
578 //EED2 } else if( typ == marParameters::SNAKE || typ == marParameters::DERICHE ) {
579 //EED2 // ONLY Isocontours at the moment.
582 //EED2 cntVTK->Delete( );
583 //EED2 cpd->Delete( );
584 //EED2 conn->Delete( );
585 //EED2 cpd2->Delete( );
586 //EED2 isoStrips->Delete( );
596 // ----------------------------------------------------------------------------
597 void marAxis::createContour(int i, kVolume* vol)
601 int numberOfPoints,numberOfCells;
605 vtkPolyData* vtkPolyData_2DContour;
608 int sizeIma = getParameters( )->getSizeIma( );
609 double dimIma = getParameters( )->getDimIma( );
614 _contours[i] = new marContour( i, getParameters( ) ) ;
615 vtkPolyData_2DContour = get2Dcontour(i,vol);
617 numberOfCells = vtkPolyData_2DContour->GetNumberOfCells( );
618 cell = vtkPolyData_2DContour->GetCell( 0 );
619 pids = cell->GetPointIds( );
620 numberOfPoints = pids->GetNumberOfIds( );
622 for( j = 0; j < numberOfPoints; j++ )
624 id = pids->GetId( j );
625 vtkPolyData_2DContour->GetPoint( id, p);
628 x=x * dimIma / ( float ) sizeIma;
629 y=y * dimIma / ( float ) sizeIma;
630 _contours[ i ]->addContourPoint( x , y );
633 _contours[i]->do_spline();
634 _contours[i]->calculateVariables();
638 // ----------------------------------------------------------------------------
639 void marAxis::EraseContour(int i){
640 if (_3Dcontour[i]!=NULL){
641 _3Dcontour[i]->Delete();
645 if (_2DDiameterMax[i]!=NULL){
646 _2DDiameterMax[i]->Delete();
647 _2DDiameterMax[i]=NULL;
650 if (_2DDiameterMin[i]!=NULL){
651 _2DDiameterMin[i]->Delete();
652 _2DDiameterMin[i]=NULL;
655 if (_2Dcontours[i]!=NULL){
656 _2Dcontours[i]->Delete();
660 if (_contours[i]!=NULL){
666 // ----------------------------------------------------------------------------
667 void marAxis::replaceContour2D(int i,int size,double *vx,double *vy){
672 vtkPoints *_pts = vtkPoints::New();
673 _pts->SetNumberOfPoints(size);
676 for (j=0 ; j<size ; j++){
677 _pts->SetPoint(j, vx[j] , vy[j] , 0 );
679 // _pts->SetPoint(0, vx[0] , vy[0] , 0 );
681 vtkCellArray *lines = vtkCellArray::New();
682 lines->InsertNextCell( size );
683 for ( j=0 ; j<size+1 ; j++ ){
684 lines->InsertCellPoint(j % size );
687 vtkPolyData *_pd = vtkPolyData::New();
688 _pd->SetPoints( _pts );
689 _pd->SetLines( lines );
690 lines->Delete(); //do not delete lines ??
694 createContour(i,NULL);
696 // ----------------------------------------------------------------------------
697 void marAxis::createSlice(int i , kVolume* vol){
699 int sizeIma = getParameters( )->getSizeIma( );
701 vtkStructuredPoints* stPoints = vtkStructuredPoints::New( );
706 stPoints->GetPointData( )->SetScalars( get3DSlice(k,vol)->GetOutput()->GetPointData()->GetScalars() );
707 stPoints->SetDimensions( sizeIma, sizeIma, 1 );
710 stPoints->SetScalarType( vol->castVtk( )->GetDataObjectType( ) );
711 stPoints->SetScalarTypeToShort ();
714 vtkImageChangeInformation* change = vtkImageChangeInformation::New();
715 change->SetInput( stPoints );
716 change->Update(); //important
719 if (_slices[k]!=NULL) { delete _slices[k]; }
720 _slices[k] = new kVolume( change->GetOutput() ) ;
725 // ----------------------------------------------------------------------------
726 void marAxis::create3DSlice(int i , kVolume* vol )
730 int sizeIma = getParameters( )->getSizeIma( );
731 double dimIma = getParameters( )->getDimIma( );
736 // double voxSize = getParameters( )->getVoxelSize( );
739 // for( unsigned int i = 0; i < _contours.size( ); i++ )
740 // delete _contours[ i ];
741 // _contours.clear( );
743 // for( unsigned int i = 0; i < _slices.size( ); i++ ) {
744 // delete _slices[ i ];
745 // _3Dslices[ i ]->Delete( );
746 // _quantificationImages[ i ]->Delete( );
750 // _3Dslices.clear( );
751 // _quantificationImages.clear( );
753 // "Cutter" initialization
754 // vtkStructuredPoints* stPoints = vtkStructuredPoints::New( );
758 // for( int k = 0; k < nCnts; k++ ) {
759 //s = ( double )k / ( double )( nCnts - 1 );
763 vtkPlaneSource* pSource = vtkPlaneSource::New( );
764 pSource->SetOrigin( p[ 0 ], p[ 1 ], p[ 2 ] );
765 pSource->SetPoint1( p[ 0 ] + dimIma - 1.0, p[ 1 ], p[ 2 ] );
766 pSource->SetPoint2( p[ 0 ], p[ 1 ], p[ 2 ] + dimIma - 1.0 );
767 pSource->SetResolution( sizeIma - 1, sizeIma - 1 );
769 pSource->SetCenter( p[ 0 ], p[ 1 ], p[ 2 ] );
770 pSource->SetNormal( n[ 0 ], n[ 1 ], n[ 2 ] );
773 if (_3Dslices[ k ]!=NULL){ _3Dslices[ k ]->Delete(); }
774 _3Dslices[ k ] = vtkProbeFilter::New( ) ;
775 _3Dslices[ k ]->SetInput( ( vtkDataSet* )pSource->GetOutput( ) );
776 _3Dslices[ k ]->SetSource( vol->castVtk( ) );
777 _3Dslices[ k ]->Update( );
780 // stPoints->GetPointData( )->
781 // SetScalars( _3Dslices[ k ]->GetOutput( )->
782 // GetPointData( )->GetScalars( ) );
783 // stPoints->SetDimensions( sizeIma, sizeIma, 1 );
784 // stPoints->SetScalarType( vol->castVtk( )->GetDataObjectType( ) );
785 // stPoints->Update( );
789 // marContour* tmp_marContour = new marContour( k, getParameters( ) ) ;
790 // tmp_marContour->calculateVariables();
791 // _contours.push_back( tmp_marContour );
793 // _slices.push_back( new kVolume( kVolume::UCHAR, 1, 1, 1 ) );
794 // *( _slices[ k ] ) = ( vtkImageData* )stPoints; //copy
795 // _quantificationImages.push_back( (_slices[ k ])->castVtk( ) );
798 // stPoints->Delete( );
800 // _3Dslices[i]=NULL;
806 // ----------------------------------------------------------------------------
807 void marAxis::create3Dcontour(int i, kVolume* vol)
809 vtkPoints *_vtkPoints;
815 double c[3],p1[3],p2[3];
818 _vtkPoints = vtkPoints::New();
821 mat = vtkMatrix4x4::New();
822 trans = vtkTransform::New();
825 nn[0]=n[0]; nn[1]=n[1]; nn[2]=n[2];
826 nC=sqrt( nn[0]*nn[0] + nn[1]*nn[1] + nn[2]*nn[2] );
827 nn[0]=nn[0]/nC; nn[1]=nn[1]/nC; nn[2]=nn[2]/nC;
829 vtkPlaneSource* pSource = vtkPlaneSource::New( );
830 pSource->SetOrigin( 0, 0 , 0 );
831 pSource->SetPoint1( 1, 0 , 0 );
832 pSource->SetPoint2( 0, 0 , 1.0 );
833 pSource->SetCenter( 0,0,0 );
834 pSource->SetNormal( nn[ 0 ], nn[ 1 ], nn[ 2 ] );
836 // pSource->Update( );
837 pSource->GetOrigin(c);
838 pSource->GetPoint1(p1);
839 pSource->GetPoint2(p2);
841 p1[0]=p1[0]-c[0]; p1[1]=p1[1]-c[1]; p1[2]=p1[2]-c[2];
842 p2[0]=p2[0]-c[0]; p2[1]=p2[1]-c[1]; p2[2]=p2[2]-c[2];
843 nA=sqrt( p1[0]*p1[0] + p1[1]*p1[1] + p1[2]*p1[2] );
844 nB=sqrt( p2[0]*p2[0] + p2[1]*p2[1] + p2[2]*p2[2] );
845 p1[0]=p1[0]/nA; p1[1]=p1[1]/nA; p1[2]=p1[2]/nA;
846 p2[0]=p2[0]/nB; p2[1]=p2[1]/nB; p2[2]=p2[2]/nB;
848 mat->SetElement (0,0, nn[0]);
849 mat->SetElement (1,0, nn[1]);
850 mat->SetElement (2,0, nn[2]);
851 mat->SetElement (3,0, 0);
852 mat->SetElement (0,1, p2[0]);
853 mat->SetElement (1,1, p2[1]);
854 mat->SetElement (2,1, p2[2]);
855 mat->SetElement (3,1, 0);
856 mat->SetElement (0,2, p1[0]);
857 mat->SetElement (1,2, p1[1]);
858 mat->SetElement (2,2, p1[2]);
859 mat->SetElement (3,2, 0);
860 mat->SetElement (0,3, 0);
861 mat->SetElement (1,3, 0);
862 mat->SetElement (2,3, 0);
863 mat->SetElement (3,3, 1);
865 // double deter=mat->Determinant(mat);
867 trans->SetMatrix(mat);
870 trans->RotateWXYZ ( ang,0,1,0);
874 // trans->Scale ( scale , scale, scale);
878 vtkMatrix4x4 *m=vtkMatrix4x4::New();
879 trans->GetMatrix ( m );
882 int j,numberOfPoints;
883 marContour* contour2D = getContour(i,vol);
884 numberOfPoints = contour2D->getNumberOfPoints( );
889 for( j = 0; j <= numberOfPoints; j++ ) {
890 contour2D->getPoint( fpA,j % numberOfPoints);
892 pA[0] = fpA[0] + 0.2; // correction EED
893 pA[1] = fpA[1] + 0.2; // correction EED
897 // Why this does not works...(release version)????
898 // trans->TransformPoint(pA,ppp);
900 // trans->TransformVector(pA,ppp);
902 ppp[0] = m->GetElement(0,0)*pA[0] + m->GetElement(0,1)*pA[1] + m->GetElement(0,2)*pA[2] + m->GetElement(0,3);
903 ppp[1] = m->GetElement(1,0)*pA[0] + m->GetElement(1,1)*pA[1] + m->GetElement(1,2)*pA[2] + m->GetElement(1,3);
904 ppp[2] = m->GetElement(2,0)*pA[0] + m->GetElement(2,1)*pA[1] + m->GetElement(2,2)*pA[2] + m->GetElement(2,3);
906 // ppp[0]=ppp[0]+o[0]-c[0]; ppp[1]=ppp[1]+o[1]-c[1]; ppp[2]=ppp[2]+o[2]-c[2];
907 ppp[0]=ppp[0]+o[0]; ppp[1]=ppp[1]+o[1]; ppp[2]=ppp[2]+o[2];
908 // ppp[0]=ppp[0]+o[0]-c[0]*0.781; ppp[1]=ppp[1]+o[1]-c[1]*0.781;; ppp[2]=ppp[2]+o[2]-c[2]*0.781;;
909 _vtkPoints->InsertNextPoint( (float)ppp[0],(float)ppp[1],(float)ppp[2] );
915 if (_3Dcontour[i]!=NULL){ _3Dcontour[i]->Delete(); }
916 _3Dcontour[i]=_vtkPoints;
919 // ----------------------------------------------------------------------------
920 void marAxis::createSliceImage(int i , kVolume* vol )
922 _quantificationImages[i] = getSlice(i,vol)->castVtk( ) ;
925 // ----------------------------------------------------------------------------
926 void marAxis::create2Dcontour(int i , kVolume* vol )
928 vtkImageData* imagedata;
930 double thr = getParameters( )->getDoubleParam( marParameters::e_threshold_isocontour );
931 int sizeIma = getParameters( )->getSizeIma( );
932 double localint = getSignal( i );
933 double sigma = getParameters( )->getDoubleParam( marParameters::e_sigma );
934 double threshold = localint * thr / 100.0;
936 int dimIma = ( int ) ( ( double ) sizeIma /getParameters( )->getDoubleParam( marParameters::e_scale ) );
937 int vmin = ( int )threshold;
938 int vmax = ( int )threshold;
943 x = ( x != -1 )? x: sizeIma / 2;
944 y = ( y != -1 )? y: sizeIma / 2;
947 imagedata = (vtkImageData*) ((getSlice( i , vol ))->castVtk( ));
948 double *range = imagedata->GetScalarRange();
950 //vtkKitwareContourFilter* cntVTK = vtkKitwareContourFilter::New( );
951 vtkContourFilter* cntVTK = vtkContourFilter::New( );
952 cntVTK->SetInput( imagedata );
953 cntVTK->SetNumberOfContours( 1 );
954 //cntVTK->SetValue( 0, vmin );
955 cntVTK->SetValue( 0, (range[1]*thr/100) );
956 //cntVTK->SetValue( 1, vmax );
959 vtkCleanPolyData* cpd = vtkCleanPolyData::New( );
960 cpd->SetInput( cntVTK->GetOutput( ) );
961 cpd->ConvertLinesToPointsOff( );
964 vtkPolyDataConnectivityFilter* conn = vtkPolyDataConnectivityFilter::New( );
965 //conn->SetMaxRecursionDepth( 3000 );
967 // PS -> //conn->SetInput( cntVTK->GetOutput( ) ); PS
968 conn->SetInput( cpd->GetOutput( ) );
970 conn->SetClosestPoint( x, y, 0 );
971 conn->SetExtractionModeToClosestPointRegion( );
974 vtkCleanPolyData* cpd2 = vtkCleanPolyData::New( );
975 cpd2->SetInput( conn->GetOutput( ) );
978 vtkStripper* vtkstripper = vtkStripper::New( );
979 vtkstripper->SetInput( cpd2->GetOutput() );
980 vtkstripper->Update();
982 vtkPolyData* polyDataResult = vtkstripper->GetOutput();
983 polyDataResult->Update( );
993 int ii,size=polyDataResult->GetNumberOfPoints();
994 printf("EED marAxis::create2Dcontour size %d\n",size);
995 for (ii=0;ii<size;ii++)
997 p=polyDataResult->GetPoint(ii);
998 printf("EED marAxis::create2Dcontour %d >> %f %f %f\n",ii,p[0], p[1], p[2]);
1002 for (ii=0;ii<size;ii++)
1004 p=polyDataResult->GetPoint(ii);
1005 printf("EED marAxis::create2Dcontour xxxxx %d >> %f %f %f\n",ii,p[0], p[1], p[2]);
1009 if ( _2Dcontours[i]!=NULL ) { _2Dcontours[i]->Delete(); }
1010 _2Dcontours[i]=polyDataResult;
1012 // ----------------------------------------------------------------------------
1013 void marAxis::create2DDiameterMin( int i , kVolume* vol ){
1015 marContour *marcontour=getContour(i,vol);
1016 marcontour->getMinimumLine(p1,p2);
1017 vtkPoints *_vtkPoints = vtkPoints::New();
1018 int sizeIma = getParameters( )->getSizeIma( );
1020 double dimIma = getParameters( )->getDimIma( );
1025 double factor=( float ) sizeIma / dimIma;
1026 p1[0]=p1[0]*factor+sizeIma/2;
1027 p1[1]=p1[1]*factor+sizeIma/2;
1028 p2[0]=p2[0]*factor+sizeIma/2;
1029 p2[1]=p2[1]*factor+sizeIma/2;
1030 _vtkPoints->InsertNextPoint( p1[0] , p1[1] , 1 );
1031 _vtkPoints->InsertNextPoint( p2[0] , p2[1] , 1 );
1032 if ( _2DDiameterMin[i]!=NULL ) { _2DDiameterMin[i]->Delete(); }
1033 _2DDiameterMin[i] = _vtkPoints;
1035 // ----------------------------------------------------------------------------
1036 void marAxis::create2DDiameterMax( int i , kVolume* vol ){
1038 marContour *marcontour=getContour(i,vol);
1039 marcontour->getMaximumLine(p1,p2);
1040 vtkPoints *_vtkPoints = vtkPoints::New();
1041 int sizeIma = getParameters( )->getSizeIma( );
1043 double dimIma = getParameters( )->getDimIma( );
1047 double factor=( float ) sizeIma / dimIma;
1048 p1[0]=p1[0]*factor+sizeIma/2;
1049 p1[1]=p1[1]*factor+sizeIma/2;
1050 p2[0]=p2[0]*factor+sizeIma/2;
1051 p2[1]=p2[1]*factor+sizeIma/2;
1052 _vtkPoints->InsertNextPoint( p1[0] , p1[1] , 1 );
1053 _vtkPoints->InsertNextPoint( p2[0] , p2[1] , 1 );
1054 if ( _2DDiameterMax[i]!=NULL ) { _2DDiameterMax[i]->Delete(); }
1055 _2DDiameterMax[i] = _vtkPoints;
1057 // ----------------------------------------------------------------------------
1058 void marAxis::createEmptyVectors()
1060 int nCnts = ( int ) getNumberOfSplinePoints( );
1063 // int xxxx=getNumberOfContours();
1064 for (i=0;i<nCnts;i++){
1065 _contours.push_back( NULL );
1066 _slices.push_back( NULL );
1067 _3Dslices.push_back( NULL );
1068 _3Dcontour.push_back( NULL );
1069 _quantificationImages.push_back( NULL );
1070 _2Dcontours.push_back( NULL );
1071 _2DDiameterMin.push_back( NULL );
1072 _2DDiameterMax.push_back( NULL );
1075 // ----------------------------------------------------------------------------
1076 void marAxis::clearAllVectors()
1080 for (i=0;i<_contours.size();i++){
1081 if (_contours[i]!=NULL) { // DATA-MODEL-2D
1082 delete _contours[i];
1088 for (i=0;i<_slices.size();i++){ // DATA-MODEL-Voxel XxYx1
1089 if (_slices[i]!=NULL) {
1097 for (i=0;i<_3Dslices.size();i++){ // VISUALISATION_VTK 3D
1098 if (_3Dslices[i]!=NULL) {
1099 _3Dslices[i]->Delete();
1105 for (i=0;i<_3Dcontour.size();i++){ // VISUALISATION_VTK 3D
1106 if (_3Dcontour[i]!=NULL) {
1107 _3Dcontour[i]->Delete();
1113 for (i=0;i<_quantificationImages.size();i++){ // VISUALISATION_VTK 2D
1114 if (_quantificationImages[i]!=NULL) {
1115 // _quantificationImages[i]->Delete();
1116 // _quantificationImages[i]=NULL;
1119 _quantificationImages.clear();
1121 for (i=0;i<_2Dcontours.size();i++){ // VISUALISATION_VTK 2D
1122 if (_2Dcontours[i]!=NULL) {
1123 // EED ??? This object was allready erased bye the VTK pipeline
1124 // _2Dcontours[i]->Delete();
1125 _2Dcontours[i]=NULL;
1128 _2Dcontours.clear();
1130 for (i=0;i<_2DDiameterMin.size();i++){ // VISUALISATION_VTK 2D
1131 if (_2DDiameterMin[i]!=NULL) {
1132 _2DDiameterMin[i]->Delete() ;
1133 _2DDiameterMin[i]=NULL;
1136 _2DDiameterMin.clear();
1138 for (i=0;i<_2DDiameterMax.size();i++){ // VISUALISATION_VTK 2D
1139 if (_2DDiameterMax[i]!=NULL) {
1140 _2DDiameterMax[i]->Delete();
1141 _2DDiameterMax[i]=NULL;
1144 _2DDiameterMax.clear();
1147 // ----------------------------------------------------------------------------
1148 void marAxis::eraseContourVectorsContent()
1151 for (i=0;i<_contours.size();i++){
1152 if (_contours[i]!=NULL) { // DATA-MODEL-2D
1153 delete _contours[i];
1157 for (i=0;i<_2Dcontours.size();i++){ // VISUALISATION_VTK 2D
1158 if (_2Dcontours[i]!=NULL) {
1159 _2Dcontours[i]->Delete();
1160 _2Dcontours[i]=NULL;
1163 for (i=0;i<_2DDiameterMin.size();i++){ // VISUALISATION_VTK 2D
1164 if (_2DDiameterMin[i]!=NULL) {
1165 _2DDiameterMin[i]->Delete();
1166 _2DDiameterMin[i]=NULL;
1169 for (i=0;i<_2DDiameterMax.size();i++){ // VISUALISATION_VTK 2D
1170 if (_2DDiameterMax[i]!=NULL) {
1171 _2DDiameterMax[i]->Delete();
1172 _2DDiameterMax[i]=NULL;
1175 for (i=0;i<_3Dcontour.size();i++){ // VISUALISATION_VTK 3D
1176 if (_3Dcontour[i]!=NULL) {
1177 _3Dcontour[i]->Delete();
1183 // ----------------------------------------------------------------------------
1184 marContour* marAxis::getContour(int i, kVolume* vol) // DATA-MODEL 2D
1186 if (_contours[i]==NULL){ createContour(i, vol); }
1187 return _contours[i];
1189 // ----------------------------------------------------------------------------
1190 kVolume* marAxis::getSlice(int i, kVolume* vol) // DATA-MODEL-Voxel XxYx1
1192 if (_slices[i]==NULL){ createSlice( i , vol ); }
1195 // ----------------------------------------------------------------------------
1197 bool marAxis::if3DcontourExist(int i)
1200 if (_3Dcontour[i]==NULL)
1207 // ----------------------------------------------------------------------------
1208 void marAxis::Save3Dcontour(FILE *ff,int id)
1212 if (_3Dcontour[id]!=NULL)
1214 //Ramiro fprintf(ff,"contour_id: %d \n", id);
1215 size = _3Dcontour[id]->GetNumberOfPoints();
1216 //Ramiro fprintf(ff,"numberOfPoints: %d \n", size);
1217 for ( i=0 ; i<size ;i++ ){
1218 _3Dcontour[id]->GetPoint( i , point );
1219 fprintf(ff,"%f %f %f %d\n", point[0], point[1],point[2],id);
1223 // ----------------------------------------------------------------------------
1224 void marAxis::SaveExisting3DContours(FILE *ff)
1228 // counting existing 3Dcontours
1230 int i,size=_3Dcontour.size();
1231 for (i=0; i<size; i++)
1233 if (_3Dcontour[i]!=NULL)
1238 //Ramiro fprintf(ff, "NumberOf3DContours: %d \n", acum);
1239 // saving existing 3Dcontours
1240 for (i=0; i<size; i++)
1242 Save3Dcontour(ff,i);
1247 // ----------------------------------------------------------------------------
1248 vtkProbeFilter* marAxis::get3DSlice(int i, kVolume* vol){ // VISUALISATION-VTK 3D
1249 if (_3Dslices[i]==NULL){ create3DSlice(i,vol); }
1250 return _3Dslices[i];
1252 // ----------------------------------------------------------------------------
1253 vtkPoints* marAxis::get3Dcontour(int i, kVolume* vol){ // VISUALISATION-VTK 3D
1254 if (_3Dcontour[i]==NULL){ create3Dcontour(i, vol); }
1255 return _3Dcontour[i];
1257 // ----------------------------------------------------------------------------
1258 vtkImageData* marAxis::getSliceImage(int i, kVolume* vol){ // VISUALISATION-VTK 2D
1259 if (_quantificationImages[i]==NULL){ createSliceImage(i,vol); }
1260 return _quantificationImages[i];
1262 // ----------------------------------------------------------------------------
1263 vtkPolyData* marAxis::get2Dcontour(int i , kVolume* vol){ // VISUALISATION-VTK 2D
1264 if (_2Dcontours[i]==NULL){ create2Dcontour(i,vol); }
1265 return _2Dcontours[i];
1267 // ----------------------------------------------------------------------------
1268 vtkPoints* marAxis::get2DDiameterMin(int i , kVolume* vol){ // VISUALISATION-VTK 2D
1269 if (_2DDiameterMin[i]==NULL){ create2DDiameterMin(i,vol); }
1270 return _2DDiameterMin[i];
1272 // ----------------------------------------------------------------------------
1273 vtkPoints* marAxis::get2DDiameterMax(int i , kVolume* vol){ // VISUALISATION-VTK 2D
1274 if (_2DDiameterMax[i]==NULL){ create2DDiameterMax(i,vol); }
1275 return _2DDiameterMax[i];
1277 // ----------------------------------------------------------------------------
1278 double marAxis::getAxisLenght(int pIni,int pEnd){
1286 marVector* pO = new marVector(2);
1287 marVector* pF = new marVector(2);
1290 for( j = pIni, L = 0.0; j < pEnd-1; j++ ) {
1300 // ----------------------------------------------------------------------------
1301 void marAxis::calculateTotalAxisLenght(){
1302 _totalAxisLenght = getAxisLenght( 0 , _points.size( ) );
1304 // ----------------------------------------------------------------------------
1305 double marAxis::getTotalLength(){
1306 if (_totalAxisLenght==-1){
1307 calculateTotalAxisLenght();
1309 return _totalAxisLenght;
1312 // ----------------------------------------------------------------------------
1313 void marAxis::calculateSubAxisLength(){
1314 _subAxisLenght = getAxisLenght( _startQuant , _finishQuant );
1316 // ----------------------------------------------------------------------------
1317 double marAxis::getSubAxisLength(){
1318 if (_subAxisLenght==-1){
1319 calculateSubAxisLength();
1321 return _subAxisLenght;
1323 // ----------------------------------------------------------------------------
1324 double marAxis::getAverageArea(int pIni, int pEnd, kVolume* vol){
1325 marContour *marcontourHealthy;
1326 int ihealthySlice,itemp;
1327 double acumArea = 0;
1328 for (ihealthySlice = pIni ; ihealthySlice <= pEnd ; ihealthySlice++ ){
1329 itemp=ihealthySlice;
1330 if (itemp<0) { itemp=0; }
1331 if (itemp>=_points.size( )) { itemp=_points.size( )-1; }
1332 marcontourHealthy = getContour( itemp , vol );
1333 acumArea = acumArea + marcontourHealthy->getArea();
1335 return acumArea / (pEnd - pIni + 1);
1337 // ----------------------------------------------------------------------------
1338 void marAxis::calculateReferenceArea(kVolume* vol){
1339 _referenceArea = getAverageArea(_healthySliceStart,_healthySliceEnd, vol);
1341 // ----------------------------------------------------------------------------
1342 double marAxis::getReferenceArea(kVolume* vol){
1343 if (_referenceArea==-1){
1344 calculateReferenceArea(vol);
1346 return _referenceArea;
1348 // ----------------------------------------------------------------------------
1349 void marAxis::calculateReferenceAverDiam(kVolume* vol){
1350 marContour *marcontourHealthy;
1353 for (ihealthySlice=_healthySliceStart; ihealthySlice<=_healthySliceEnd; ihealthySlice++){
1354 marcontourHealthy = getContour( ihealthySlice , vol );
1355 acumDiam = acumDiam + marcontourHealthy->getAverageDiameter();
1357 _referenceAverDiam = acumDiam / (_healthySliceEnd-_healthySliceStart+1);
1359 // ----------------------------------------------------------------------------
1360 double marAxis::getReferenceAverDiam(kVolume* vol){
1361 if (_referenceAverDiam==-1){
1362 calculateReferenceAverDiam(vol);
1364 return _referenceAverDiam;
1367 // ----------------------------------------------------------------------------
1368 void marAxis::reset( )
1380 _subAxisLenght = -1;
1381 _referenceArea = -1;
1382 _referenceAverDiam = -1;
1388 for( i = 0; i < _contours.size( ); i++ )
1389 delete _contours[ i ];
1391 for( i = 0; i < _slices.size( ); i++ )
1392 delete _slices[ i ];
1397 // ----------------------------------------------------------------------------
1398 void marAxis::copyFrom( const marObject& from )
1402 // ----------------------------------------------------------------------------
1403 bool marAxis::save( std::ofstream& os )
1405 double p[ INDX_count ];
1408 i = _description.length( );
1409 os.write( ( const char* )&i, sizeof( int ) );
1410 os.write( ( char* )_description.c_str( ), i * sizeof( char ) );
1412 i = getNumberOfControlPoints( );
1413 os.write( ( const char* )&i, sizeof( int ) );
1414 for( i = 0; i < getNumberOfControlPoints( ); i++ ) {
1416 memcpy( p, _controlPoints[ i ], 3 * sizeof( double ) );
1417 memcpy( p + 3, _controlState[ i ],
1418 ( INDX_count - 3 ) * sizeof( double ) );
1419 os.write( ( const char* )p, INDX_count * sizeof( double ) );
1422 i = getNumberOfContours( );
1423 os.write( ( const char* )&i, sizeof( int ) );
1424 for( i = 0; i < getNumberOfContours( ); i++ )
1425 _contours[ i ]->save( os );
1430 // ----------------------------------------------------------------------------
1431 bool marAxis::load( std::ifstream& is )
1433 double p[ INDX_count ];
1438 is.read( ( char* )&n, sizeof( int ) );
1439 _description.resize( n );
1440 is.read( ( char* )_description.c_str( ), n * sizeof( char ) );
1442 is.read( ( char* )&n, sizeof( int ) );
1443 for( i = 0; i < n; i++ ) {
1445 is.read( ( char* )p, INDX_count * sizeof( double ) );
1446 addControlPoint( p, p + 3 );
1449 is.read( ( char* )&n, sizeof( int ) );
1450 for( i = 0; i < n; i++ ) {
1452 _contours.push_back( new marContour( 0, getParameters( ) ) );
1453 _contours[ i ]->load( is );
1458 // ----------------------------------------------------------------------------
1459 vtkPolyData* marAxis::Draw( )
1464 vtkPoints* allPoints = vtkPoints::New( );
1465 vtkCellArray* allTopology = vtkCellArray::New( );
1467 allTopology->InsertNextCell( _points.size( ) );
1468 for( i = 0, j=0; i < _points.size( ); i++, j++ ) {
1470 allPoints->InsertNextPoint( p[ 0 ], p[ 1 ], p[ 2 ] );
1471 allTopology->InsertCellPoint( j );
1474 _allData = vtkPolyData::New( );
1475 _allData->SetPoints( allPoints );
1476 _allData->SetLines( allTopology );
1477 allPoints->Delete();
1478 allTopology->Delete();
1480 return ( _allData );
1482 // ----------------------------------------------------------------------------
1483 vtkPolyData *marAxis::GetAxisData( )
1488 double p[marAxis::INDX_count];
1491 vtkPoints *allPoints = vtkPoints::New( );
1492 vtkCellArray *allTopology = vtkCellArray::New( );
1493 vtkDoubleArray *allRadios = vtkDoubleArray::New( );
1494 allRadios->SetName("radio");
1496 for( i = 0; i < _controlPoints.size( ); i++ ) {
1497 getControlPoint(i, p, p+3);
1498 point[ 0 ] = p[marAxis::INDX_X];
1499 point[ 1 ] = p[marAxis::INDX_Y];
1500 point[ 2 ] = p[marAxis::INDX_Z];
1501 radio = p[marAxis::INDX_RAYON]*2.0;
1502 allPoints -> InsertPoint( i, point[ 0 ] , point[ 1 ] , point[ 2 ] ); //para saber exactamante el indice que se le asigno
1504 allRadios -> InsertValue(i,radio);
1507 allTopology->InsertNextCell(2);
1508 allTopology->InsertCellPoint(i);
1509 allTopology->InsertCellPoint(i-1);
1514 _allData = vtkPolyData::New( );
1515 _allData -> SetPoints( allPoints );
1516 _allData -> SetLines( allTopology );
1517 _allData -> GetPointData()->SetScalars(allRadios);
1520 return ( _allData );
1522 // ----------------------------------------------------------------------------
1523 void marAxis::Delete( ){
1527 for (i=0;i<_points.size();i++){
1528 if (_points[i]!=NULL) {
1544 // ----------------------------------------------------------------------------
1545 void marAxis::setHealthySlice( int hsS, int hs, int hsE ){
1546 _healthySliceStart = hsS;
1548 _healthySliceEnd = hsE;
1550 _referenceAverDiam=-1;
1552 // ----------------------------------------------------------------------------
1553 void marAxis::setStartQuant( int sq ){
1555 _subAxisLenght = -1;
1557 // ----------------------------------------------------------------------------
1558 void marAxis::setFinishQuant( int fq ){
1560 _subAxisLenght = -1;
1563 // ----------------------------------------------------------------------------
1564 void marAxis::AddPointToList(double x, double y, double z, int signal)
1567 p=(double *)malloc(sizeof(double)*3);
1571 _points.push_back( p );
1572 _signal.push_back( signal );