1 /*=========================================================================
4 Module: $RCSfile: volume.cxx,v $
6 Date: $Date: 2009/05/14 13:54:43 $
7 Version: $Revision: 1.5 $
9 Copyright: (c) 2002, 2003
12 This software is distributed WITHOUT ANY WARRANTY; without even
13 the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
14 PURPOSE. See the above copyright notice for more information.
16 =========================================================================*/
18 // PS -> //#include <gsl/gsl_math.h>//PS
24 #include <vtkCharArray.h>
25 #include <vtkDataArray.h>
26 #include <vtkDoubleArray.h>
27 #include <vtkFloatArray.h>
28 #include <vtkIntArray.h>
29 #include <vtkPointData.h>
30 #include <vtkShortArray.h>
31 #include <vtkUnsignedCharArray.h>
32 #include <vtkUnsignedIntArray.h>
33 #include <vtkUnsignedShortArray.h>
35 // -------------------------------------------------------------------------
36 const vtkIdType kVolume::VTKTypes[] = { VTK_CHAR, VTK_FLOAT, VTK_DOUBLE,
37 VTK_INT, VTK_SHORT, VTK_UNSIGNED_CHAR,
38 VTK_UNSIGNED_INT, VTK_UNSIGNED_SHORT };
40 #endif // KGFO_USE_VTK
44 // -------------------------------------------------------------------------
45 //const int kVolume::IDOTypes[] = { VOL_CHAR, VOL_FLOAT, VOL_DOUBLE, VOL_LONG,
46 // VOL_SHORT, VOL_UCHAR, VOL_ULONG,
49 //#endif // KGFO_USE_IDO
51 // -------------------------------------------------------------------------
52 const void* kVolume::BLANK = ( void* ) 0;
53 const void* kVolume::NOALLOC = ( void* )-1;
54 const int kVolume::SIZETypes[] = { sizeof( char ), sizeof( float ),
55 sizeof( double ), sizeof( int ),
56 sizeof( short ), sizeof( uchar ),
57 sizeof( uint ), sizeof( ushort ) };
59 // ---------------------------------------------------------------------------
60 template< class FROM, class TO >
61 static void convertCastT( FROM* src, TO* dest, ulong size )
63 FROM* end = src + size;
64 for( ; src < end; src++, dest++ ) *dest = ( TO )*src;
68 // ---------------------------------------------------------------------------
69 template< class FROM, class TO >
70 static void convertScaleT( FROM *src, TO *dest, ulong size,
71 double smin, double tmin, double slope )
73 FROM* end = src + size;
74 for( ; src < end; src++, dest++ )
75 *dest = ( TO )( ( double( *src ) - smin ) * slope + tmin );
78 // ---------------------------------------------------------------------------
79 template< class TYPE >
80 static void getMinMaxT( TYPE *src, ulong size,
81 double& min, double& max )
83 TYPE* end = src + size;
89 for( st = true; src < end; src++, st = false ) {
91 if( *src < m || st ) m = *src;
92 if( *src > M || st ) M = *src;
100 // -------------------------------------------------------------------------
102 : _type( UCHAR ), _creator( SELF ),
105 #endif // KGFO_USE_VTK
106 //#ifdef KGFO_USE_IDO
107 // _privateIdo( NULL ),
108 //#endif // KGFO_USE_IDO
109 _raw( NULL ), _columns( NULL ),
112 _dims[ CX ] = 1; _dims[ CY ] = 1; _dims[ CZ ] = 1;
113 _sizes[ CX ] = 1; _sizes[ CY ] = 1; _sizes[ CZ ] = 1;
118 // -------------------------------------------------------------------------
119 kVolume::kVolume( Type type,
120 uint xdim, uint ydim, uint zdim,
121 double xsize, double ysize, double zsize,
123 : _type( type ), _creator( SELF ),
126 #endif // KGFO_USE_VTK
127 //#ifdef KGFO_USE_IDO
128 // _privateIdo( NULL ),
129 //#endif // KGFO_USE_IDO
130 _raw( NULL ), _columns( NULL ),
133 _dims[ CX ] = xdim; _dims[ CY ] = ydim; _dims[ CZ ] = zdim;
134 _sizes[ CX ] = xsize; _sizes[ CY ] = ysize; _sizes[ CZ ] = zsize;
135 if( data != NOALLOC ) {
146 // -------------------------------------------------------------------------
147 kVolume::kVolume( Type type,
151 : _type( type ), _creator( SELF ),
154 #endif // KGFO_USE_VTK
155 //#ifdef KGFO_USE_IDO
156 // _privateIdo( NULL ),
157 //#endif // KGFO_USE_IDO
158 _raw( NULL ), _columns( NULL ),
161 memcpy( _dims, dims, 3 * sizeof( uint ) );
162 memcpy( _sizes, sizes, 3 * sizeof( double ) );
163 if( data != NOALLOC ) {
174 // -------------------------------------------------------------------------
175 kVolume::kVolume( const kVolume& org )
176 : _type( UCHAR ), _creator( SELF ),
179 #endif // KGFO_USE_VTK
180 //#ifdef KGFO_USE_IDO
181 // _privateIdo( NULL ),
182 //#endif // KGFO_USE_IDO
183 _raw( NULL ), _columns( NULL ),
189 // -------------------------------------------------------------------------
190 kVolume& kVolume::operator=( const kVolume& org )
196 // -------------------------------------------------------------------------
197 void kVolume::copyFrom( const kVolume& org )
205 memcpy( _dims, org._dims, 3 * sizeof( uint ) );
206 memcpy( _sizes, org._sizes, 3 * sizeof( double ) );
211 memcpy( _raw, org._raw, getRawSizeInBytes( ) );
216 // -------------------------------------------------------------------------
217 bool kVolume::sameDimension( const kVolume& comp )
219 return( ( _dims[ CX ] == comp._dims[ CX ] &&
220 _dims[ CY ] == comp._dims[ CY ] &&
221 _dims[ CZ ] == comp._dims[ CZ ] ) );
224 // -------------------------------------------------------------------------
225 double kVolume::getPixel( uint x, uint y, uint z ) const
231 case CHAR: p = ( double )( ( char*** )_images )[ z ][ y ][ x ]; break;
232 case FLOAT: p = ( double )( ( float*** )_images )[ z ][ y ][ x ]; break;
233 case DOUBLE: p = ( double )( ( double*** )_images )[ z ][ y ][ x ]; break;
234 case INT: p = ( double )( ( int*** )_images )[ z ][ y ][ x ]; break;
235 case SHORT: p = ( double )( ( short*** )_images )[ z ][ y ][ x ]; break;
236 case UCHAR: p = ( double )( ( uchar*** )_images )[ z ][ y ][ x ]; break;
237 case UINT: p = ( double )( ( uint*** )_images )[ z ][ y ][ x ]; break;
238 case USHORT: p = ( double )( ( ushort*** )_images )[ z ][ y ][ x ]; break;
239 default: p = 0.0; break;
246 // -------------------------------------------------------------------------
247 void kVolume::setPixel( double v, uint x, uint y, uint z )
252 case CHAR: ( ( char*** )_images )[ z ][ y ][ x ] = ( char )v; break;
253 case FLOAT: ( ( float*** )_images )[ z ][ y ][ x ] = ( float )v; break;
254 case DOUBLE: ( ( double*** )_images )[ z ][ y ][ x ] = ( double )v; break;
255 case INT: ( ( int*** )_images )[ z ][ y ][ x ] = ( int )v; break;
256 case SHORT: ( ( short*** )_images )[ z ][ y ][ x ] = ( short )v; break;
257 case UCHAR: ( ( uchar*** )_images )[ z ][ y ][ x ] = ( uchar )v; break;
258 case UINT: ( ( uint*** )_images )[ z ][ y ][ x ] = ( uint )v; break;
259 case USHORT: ( ( ushort*** )_images )[ z ][ y ][ x ] = ( ushort )v; break;
265 // -------------------------------------------------------------------------
266 void kVolume::convertCast( Type type )
268 if( _type != type ) {
271 ulong size = getRawSize( );
273 buffer = ( void* )new uchar[ size * SIZETypes[ type ] ];
280 case SHORT: convertCastT( ( char* )_raw, ( short* )buffer, size ); break;
281 case INT: convertCastT( ( char* )_raw, ( int* )buffer, size ); break;
282 case USHORT: convertCastT( ( char* )_raw, ( ushort* )buffer, size ); break;
283 case UINT: convertCastT( ( char* )_raw, ( uint* )buffer, size ); break;
284 case FLOAT: convertCastT( ( char* )_raw, ( float* )buffer, size ); break;
285 case DOUBLE: convertCastT( ( char* )_raw, ( double* )buffer, size ); break;
286 case UCHAR: convertCastT( ( char* )_raw, ( uchar* )buffer, size ); break;
295 case CHAR: convertCastT( ( short* )_raw, ( char* )buffer, size ); break;
296 case INT: convertCastT( ( short* )_raw, ( int* )buffer, size ); break;
297 case USHORT: convertCastT( ( short* )_raw, ( ushort* )buffer, size ); break;
298 case UINT: convertCastT( ( short* )_raw, ( uint* )buffer, size ); break;
299 case FLOAT: convertCastT( ( short* )_raw, ( float* )buffer, size ); break;
300 case DOUBLE: convertCastT( ( short* )_raw, ( double* )buffer, size ); break;
301 case UCHAR: convertCastT( ( short* )_raw, ( uchar* )buffer, size ); break;
309 case CHAR: convertCastT( ( int* )_raw, ( char* )buffer, size ); break;
310 case SHORT: convertCastT( ( int* )_raw, ( short* )buffer, size ); break;
311 case USHORT: convertCastT( ( int* )_raw, ( ushort* )buffer, size ); break;
312 case UINT: convertCastT( ( int* )_raw, ( uint* )buffer, size ); break;
313 case FLOAT: convertCastT( ( int* )_raw, ( float* )buffer, size ); break;
314 case DOUBLE: convertCastT( ( int* )_raw, ( double* )buffer, size ); break;
315 case UCHAR: convertCastT( ( int* )_raw, ( uchar* )buffer, size ); break;
323 case CHAR: convertCastT( ( ushort* )_raw, ( char* )buffer, size ); break;
324 case SHORT: convertCastT( ( ushort* )_raw, ( short* )buffer, size ); break;
325 case INT: convertCastT( ( ushort* )_raw, ( int* )buffer, size ); break;
326 case UINT: convertCastT( ( ushort* )_raw, ( uint* )buffer, size ); break;
327 case FLOAT: convertCastT( ( ushort* )_raw, ( float* )buffer, size ); break;
328 case DOUBLE: convertCastT( ( ushort* )_raw, ( double* )buffer, size ); break;
329 case UCHAR: convertCastT( ( ushort* )_raw, ( uchar* )buffer, size ); break;
337 case CHAR: convertCastT( ( uint* )_raw, ( char* )buffer, size ); break;
338 case SHORT: convertCastT( ( uint* )_raw, ( short* )buffer, size ); break;
339 case INT: convertCastT( ( uint* )_raw, ( int* )buffer, size ); break;
340 case USHORT: convertCastT( ( uint* )_raw, ( ushort* )buffer, size ); break;
341 case FLOAT: convertCastT( ( uint* )_raw, ( float* )buffer, size ); break;
342 case DOUBLE: convertCastT( ( uint* )_raw, ( double* )buffer, size ); break;
343 case UCHAR: convertCastT( ( uint* )_raw, ( uchar* )buffer, size ); break;
351 case CHAR: convertCastT( ( float* )_raw, ( char* )buffer, size ); break;
352 case SHORT: convertCastT( ( float* )_raw, ( short* )buffer, size ); break;
353 case INT: convertCastT( ( float* )_raw, ( int* )buffer, size ); break;
354 case USHORT: convertCastT( ( float* )_raw, ( ushort* )buffer, size ); break;
355 case UINT: convertCastT( ( float* )_raw, ( uint* )buffer, size ); break;
356 case DOUBLE: convertCastT( ( float* )_raw, ( double* )buffer, size ); break;
357 case UCHAR: convertCastT( ( float* )_raw, ( uchar* )buffer, size ); break;
365 case CHAR: convertCastT( ( double* )_raw, ( char* )buffer, size ); break;
366 case SHORT: convertCastT( ( double* )_raw, ( short* )buffer, size ); break;
367 case INT: convertCastT( ( double* )_raw, ( int* )buffer, size ); break;
368 case USHORT: convertCastT( ( double* )_raw, ( ushort* )buffer, size ); break;
369 case UINT: convertCastT( ( double* )_raw, ( uint* )buffer, size ); break;
370 case FLOAT: convertCastT( ( double* )_raw, ( float* )buffer, size ); break;
371 case UCHAR: convertCastT( ( double* )_raw, ( uchar* )buffer, size ); break;
379 case CHAR: convertCastT( ( uchar* )_raw, ( char* )buffer, size ); break;
380 case SHORT: convertCastT( ( uchar* )_raw, ( short* )buffer, size ); break;
381 case INT: convertCastT( ( uchar* )_raw, ( int* )buffer, size ); break;
382 case USHORT: convertCastT( ( uchar* )_raw, ( ushort* )buffer, size ); break;
383 case UINT: convertCastT( ( uchar* )_raw, ( uint* )buffer, size ); break;
384 case FLOAT: convertCastT( ( uchar* )_raw, ( float* )buffer, size ); break;
385 case DOUBLE: convertCastT( ( uchar* )_raw, ( double* )buffer, size ); break;
401 // -------------------------------------------------------------------------
402 void kVolume::convertScale( Type type, double min, double max )
404 double tmin, tmax, smin, smax;
406 // PS -> //tmin = GSL_MIN( min, max ); //PS
407 tmin= ((min<max)?min:max);
408 // PS -> //tmax = GSL_MAX( min, max ); //PS
409 tmax= ((min>max)?min:max);
411 getMinMax( smin, smax );
413 // compute scaling slope
414 double a = ( tmax - tmin ) / ( smax - smin );
417 ulong size = getRawSize( );
419 buffer = ( void* )new uchar[ size * SIZETypes[ type ] ];
426 case CHAR: convertScaleT( ( char* )_raw, ( char* )buffer, size, smin, tmin, a ); break;
427 case SHORT: convertScaleT( ( char* )_raw, ( short* )buffer, size, smin, tmin, a ); break;
428 case INT: convertScaleT( ( char* )_raw, ( int* )buffer, size, smin, tmin, a ); break;
429 case USHORT: convertScaleT( ( char* )_raw, ( ushort* )buffer, size, smin, tmin, a ); break;
430 case UINT: convertScaleT( ( char* )_raw, ( uint* )buffer, size, smin, tmin, a ); break;
431 case FLOAT: convertScaleT( ( char* )_raw, ( float* )buffer, size, smin, tmin, a ); break;
432 case DOUBLE: convertScaleT( ( char* )_raw, ( double* )buffer, size, smin, tmin, a ); break;
433 case UCHAR: convertScaleT( ( char* )_raw, ( uchar* )buffer, size, smin, tmin, a ); break;
440 case CHAR: convertScaleT( ( short* )_raw, ( char* )buffer, size, smin, tmin, a ); break;
441 case SHORT: convertScaleT( ( short* )_raw, ( short* )buffer, size, smin, tmin, a ); break;
442 case INT: convertScaleT( ( short* )_raw, ( int* )buffer, size, smin, tmin, a ); break;
443 case USHORT: convertScaleT( ( short* )_raw, ( ushort* )buffer, size, smin, tmin, a ); break;
444 case UINT: convertScaleT( ( short* )_raw, ( uint* )buffer, size, smin, tmin, a ); break;
445 case FLOAT: convertScaleT( ( short* )_raw, ( float* )buffer, size, smin, tmin, a ); break;
446 case DOUBLE: convertScaleT( ( short* )_raw, ( double* )buffer, size, smin, tmin, a ); break;
447 case UCHAR: convertScaleT( ( short* )_raw, ( uchar* )buffer, size, smin, tmin, a ); break;
454 case CHAR: convertScaleT( ( int* )_raw, ( char* )buffer, size, smin, tmin, a ); break;
455 case SHORT: convertScaleT( ( int* )_raw, ( short* )buffer, size, smin, tmin, a ); break;
456 case INT: convertScaleT( ( int* )_raw, ( int* )buffer, size, smin, tmin, a ); break;
457 case USHORT: convertScaleT( ( int* )_raw, ( ushort* )buffer, size, smin, tmin, a ); break;
458 case UINT: convertScaleT( ( int* )_raw, ( uint* )buffer, size, smin, tmin, a ); break;
459 case FLOAT: convertScaleT( ( int* )_raw, ( float* )buffer, size, smin, tmin, a ); break;
460 case DOUBLE: convertScaleT( ( int* )_raw, ( double* )buffer, size, smin, tmin, a ); break;
461 case UCHAR: convertScaleT( ( int* )_raw, ( uchar* )buffer, size, smin, tmin, a ); break;
468 case CHAR: convertScaleT( ( ushort* )_raw, ( char* )buffer, size, smin, tmin, a ); break;
469 case SHORT: convertScaleT( ( ushort* )_raw, ( short* )buffer, size, smin, tmin, a ); break;
470 case INT: convertScaleT( ( ushort* )_raw, ( int* )buffer, size, smin, tmin, a ); break;
471 case USHORT: convertScaleT( ( ushort* )_raw, ( ushort* )buffer, size, smin, tmin, a ); break;
472 case UINT: convertScaleT( ( ushort* )_raw, ( uint* )buffer, size, smin, tmin, a ); break;
473 case FLOAT: convertScaleT( ( ushort* )_raw, ( float* )buffer, size, smin, tmin, a ); break;
474 case DOUBLE: convertScaleT( ( ushort* )_raw, ( double* )buffer, size, smin, tmin, a ); break;
475 case UCHAR: convertScaleT( ( ushort* )_raw, ( uchar* )buffer, size, smin, tmin, a ); break;
482 case CHAR: convertScaleT( ( uint* )_raw, ( char* )buffer, size, smin, tmin, a ); break;
483 case SHORT: convertScaleT( ( uint* )_raw, ( short* )buffer, size, smin, tmin, a ); break;
484 case INT: convertScaleT( ( uint* )_raw, ( int* )buffer, size, smin, tmin, a ); break;
485 case USHORT: convertScaleT( ( uint* )_raw, ( ushort* )buffer, size, smin, tmin, a ); break;
486 case UINT: convertScaleT( ( uint* )_raw, ( uint* )buffer, size, smin, tmin, a ); break;
487 case FLOAT: convertScaleT( ( uint* )_raw, ( float* )buffer, size, smin, tmin, a ); break;
488 case DOUBLE: convertScaleT( ( uint* )_raw, ( double* )buffer, size, smin, tmin, a ); break;
489 case UCHAR: convertScaleT( ( uint* )_raw, ( uchar* )buffer, size, smin, tmin, a ); break;
496 case CHAR: convertScaleT( ( uchar* )_raw, ( char* )buffer, size, smin, tmin, a ); break;
497 case SHORT: convertScaleT( ( uchar* )_raw, ( short* )buffer, size, smin, tmin, a ); break;
498 case INT: convertScaleT( ( uchar* )_raw, ( int* )buffer, size, smin, tmin, a ); break;
499 case USHORT: convertScaleT( ( uchar* )_raw, ( ushort* )buffer, size, smin, tmin, a ); break;
500 case UINT: convertScaleT( ( uchar* )_raw, ( uint* )buffer, size, smin, tmin, a ); break;
501 case FLOAT: convertScaleT( ( uchar* )_raw, ( float* )buffer, size, smin, tmin, a ); break;
502 case DOUBLE: convertScaleT( ( uchar* )_raw, ( double* )buffer, size, smin, tmin, a ); break;
503 case UCHAR: convertScaleT( ( uchar* )_raw, ( uchar* )buffer, size, smin, tmin, a ); break;
510 case CHAR: convertScaleT( ( double* )_raw, ( char* )buffer, size, smin, tmin, a ); break;
511 case SHORT: convertScaleT( ( double* )_raw, ( short* )buffer, size, smin, tmin, a ); break;
512 case INT: convertScaleT( ( double* )_raw, ( int* )buffer, size, smin, tmin, a ); break;
513 case USHORT: convertScaleT( ( double* )_raw, ( ushort* )buffer, size, smin, tmin, a ); break;
514 case UINT: convertScaleT( ( double* )_raw, ( uint* )buffer, size, smin, tmin, a ); break;
515 case FLOAT: convertScaleT( ( double* )_raw, ( float* )buffer, size, smin, tmin, a ); break;
516 case DOUBLE: convertScaleT( ( double* )_raw, ( double* )buffer, size, smin, tmin, a ); break;
517 case UCHAR: convertScaleT( ( uchar* )_raw, ( double* )buffer, size, smin, tmin, a ); break;
524 case CHAR: convertScaleT( ( float* )_raw, ( char* )buffer, size, smin, tmin, a ); break;
525 case SHORT: convertScaleT( ( float* )_raw, ( short* )buffer, size, smin, tmin, a ); break;
526 case INT: convertScaleT( ( float* )_raw, ( int* )buffer, size, smin, tmin, a ); break;
527 case USHORT: convertScaleT( ( float* )_raw, ( ushort* )buffer, size, smin, tmin, a ); break;
528 case UINT: convertScaleT( ( float* )_raw, ( uint* )buffer, size, smin, tmin, a ); break;
529 case FLOAT: convertScaleT( ( float* )_raw, ( float* )buffer, size, smin, tmin, a ); break;
530 case DOUBLE: convertScaleT( ( float* )_raw, ( double* )buffer, size, smin, tmin, a ); break;
531 case UCHAR: convertScaleT( ( float* )_raw, ( double* )buffer, size, smin, tmin, a ); break;
545 // -------------------------------------------------------------------------
546 void kVolume::getMinMax( double& min, double& max ) const
548 ulong size = getRawSize( );
552 case CHAR: getMinMaxT( ( char* )_raw, size, min, max ); break;
553 case FLOAT: getMinMaxT( ( float* )_raw, size, min, max ); break;
554 case DOUBLE: getMinMaxT( ( double* )_raw, size, min, max ); break;
555 case INT: getMinMaxT( ( int* )_raw, size, min, max ); break;
556 case SHORT: getMinMaxT( ( short* )_raw, size, min, max ); break;
557 case UCHAR: getMinMaxT( ( uchar* )_raw, size, min, max ); break;
558 case UINT: getMinMaxT( ( uint* )_raw, size, min, max ); break;
559 case USHORT: getMinMaxT( ( ushort* )_raw, size, min, max ); break;
564 // -------------------------------------------------------------------------
565 double kVolume::getMin( ) const
573 // -------------------------------------------------------------------------
574 double kVolume::getMax( ) const
582 // -------------------------------------------------------------------------
583 /*double kVolume::GetMaxIntSphere( double* p, double r )
585 int minX, minY, minZ, maxX, maxY, maxZ;
586 gslobj_vector vP( p, 3 ), v( 3 );
590 minX = int( floor( p[ 0 ] - r ) );
591 minY = int( floor( p[ 1 ] - r ) );
592 minZ = int( floor( p[ 2 ] - r ) );
593 maxX = int( ceil( p[ 0 ] + r ) );
594 maxY = int( ceil( p[ 1 ] + r ) );
595 maxZ = int( ceil( p[ 2 ] + r ) );
597 minX = GSL_MAX( minX, 0 );
598 minY = GSL_MAX( minY, 0 );
599 minZ = GSL_MAX( minZ, 0 );
601 maxX = GSL_MIN( maxX, this->getXdim( ) - 1 );
602 maxY = GSL_MIN( maxY, this->getYdim( ) - 1 );
603 maxZ = GSL_MIN( maxZ, this->getZdim( ) - 1 );
607 for( v( 0 ) = minX; v( 0 ) <= maxX; v( 0 )++ )
608 for( v( 1 ) = minY; v( 1 ) <= maxY; v( 1 )++ )
609 for( v( 2 ) = minZ; v( 2 ) <= maxZ; v( 2 )++ )
610 if( ( v - vP ).norm2( ) <= r ) {
611 tmp = this->getPixel( ( uint )v(0), ( uint )v(1), ( uint )v(2));
612 maxint = ( tmp > maxint || start )? tmp: maxint;
618 // -------------------------------------------------------------------------
619 unsigned short kVolume::GetMaxIntSphere2( double* p, double r )
622 * unsigned short range : 0 -> 65535
623 * unsigned int range : 0 -> 2147483647 // 2^31 - 1
625 int minX, minY, minZ, maxX, maxY, maxZ;
627 unsigned short maxint = 0, tmp;
629 minX = int( floor( p[ 0 ] - r ) );
630 minY = int( floor( p[ 1 ] - r ) );
631 minZ = int( floor( p[ 2 ] - r ) );
632 maxX = int( ceil( p[ 0 ] + r ) );
633 maxY = int( ceil( p[ 1 ] + r ) );
634 maxZ = int( ceil( p[ 2 ] + r ) );
636 // PS -> //minX = GSL_MAX( minX, 0 );//PS
637 // PS -> //minY = GSL_MAX( minY, 0 );//PS
638 // PS -> //minZ = GSL_MAX( minZ, 0 );//PS
639 minX=((minX>0)?minX:0);
640 minY=((minY>0)?minY:0);
641 minZ=((minZ>0)?minZ:0);
643 // PS -> //maxX = GSL_MIN( maxX, this->getXdim( ) - 1 );//PS
644 // PS -> //maxY = GSL_MIN( maxY, this->getYdim( ) - 1 );//PS
645 // PS -> //maxZ = GSL_MIN( maxZ, this->getZdim( ) - 1 );//PS
646 int xDim=this->getXdim( ) - 1;
647 maxX= (maxX<xDim?maxX:xDim);
649 int yDim=this->getYdim( ) - 1;
650 maxY= (maxY<yDim?maxY:yDim);
652 int zDim=this->getZdim( ) - 1;
653 maxZ= (maxZ<zDim?maxZ:zDim);
655 double r2 = r*r; //need to do comparison in double ... don't know why...
656 for( v[0] = minX; v[0] <= maxX; v[0]++ )
657 for( v[1] = minY; v[1] <= maxY; v[1]++ )
658 for( v[2] = minZ; v[2] <= maxZ; v[2]++ )
659 if( ((v[0]-p[0])*(v[0]-p[0])+(v[1]-p[1])*(v[1]-p[1])+(v[2]-p[2])*(v[2]-p[2])) <= r2 )
661 tmp = (unsigned short)this->getPixel( v[0], v[1], v[2] );
662 maxint = ( tmp > maxint ) ? tmp : maxint;
668 // -------------------------------------------------------------------------
669 void kVolume::allocate( )
671 ulong size = getRawSizeInBytes( );
673 if( _creator == SELF ) {
675 _raw = ( void* )new uchar[ size ];
676 memset( _raw, 0, size );
681 // -------------------------------------------------------------------------
682 void kVolume::buildIndex( )
686 size = ( _dims[ CZ ] * sizeof( uchar** ) ) +
687 ( _dims[ CZ ] * _dims[ CY ] * sizeof( void* ) );
689 _images = ( void*** )new uchar[ size ];
691 _columns = ( void** )( _images + _dims[ CZ ] );
692 void** plane = _columns;
693 for( uint z = 0; z < _dims[ CZ ]; z++ ) {
695 _images[ z ] = plane;
696 plane += _dims[ CY ];
700 for( uint y = 0; y < _dims[ CZ ] * _dims[ CY ]; y++ ) {
702 _columns[ y ] = line;
703 line = ( void* )( ( uchar* ) line +
704 _dims[ CX ] * SIZETypes[ _type ] );
710 vtkCharArray* carray;
711 vtkDoubleArray* darray;
712 vtkFloatArray* farray;
714 vtkShortArray* sarray;
715 vtkUnsignedCharArray* ucarray;
716 vtkUnsignedIntArray* uiarray;
717 vtkUnsignedShortArray* usarray;
719 size = _dims[ CX ] * _dims[ CY ] * _dims[ CZ ];
721 if( _creator == SELF || _creator == IDO ) {
723 _vtk = vtkImageData::New( );
724 _vtk->SetDimensions( _dims[ CX ], _dims[ CY ], _dims[ CZ ] );
725 _vtk->SetSpacing( _sizes[ CX ], _sizes[ CY ], _sizes[ CZ ] );
730 carray = vtkCharArray::New( );
731 carray->SetNumberOfComponents( 1 );
732 carray->SetArray( ( char* )( _raw ), size, 1 );
733 _vtk->SetScalarType( VTK_CHAR );
734 _vtk->GetPointData( )->SetScalars( carray );
739 ucarray = vtkUnsignedCharArray::New( );
740 ucarray->SetNumberOfComponents( 1 );
741 ucarray->SetArray( ( uchar* )( _raw ), size, 1 );
742 _vtk->SetScalarType( VTK_UNSIGNED_CHAR );
743 _vtk->GetPointData( )->SetScalars( ucarray );
748 sarray = vtkShortArray::New( );
749 sarray->SetNumberOfComponents( 1 );
750 sarray->SetArray( ( short* )( _raw ), size, 1 );
751 _vtk->SetScalarType( VTK_SHORT );
752 _vtk->GetPointData( )->SetScalars( sarray );
757 iarray = vtkIntArray::New( );
758 iarray->SetNumberOfComponents( 1 );
759 iarray->SetArray( ( int* )( _raw ), size, 1 );
760 _vtk->SetScalarType( VTK_INT );
761 _vtk->GetPointData( )->SetScalars( iarray );
766 usarray = vtkUnsignedShortArray::New( );
767 usarray->SetNumberOfComponents( 1 );
768 usarray->SetArray( ( ushort* )( _raw ), size, 1 );
769 _vtk->SetScalarType( VTK_UNSIGNED_SHORT );
770 _vtk->GetPointData( )->SetScalars( usarray );
775 uiarray = vtkUnsignedIntArray::New( );
776 uiarray->SetNumberOfComponents( 1 );
777 uiarray->SetArray( ( uint* )( _raw ), size, 1 );
778 _vtk->SetScalarType( VTK_UNSIGNED_INT );
779 _vtk->GetPointData( )->SetScalars( uiarray );
784 farray = vtkFloatArray::New( );
785 farray->SetNumberOfComponents( 1 );
786 farray->SetArray( ( float* )( _raw ), size, 1 );
787 _vtk->SetScalarType( VTK_FLOAT );
788 _vtk->GetPointData( )->SetScalars( farray );
793 darray = vtkDoubleArray::New( );
794 darray->SetNumberOfComponents( 1 );
795 darray->SetArray( ( double* )( _raw ), size, 1 );
796 _vtk->SetScalarType( VTK_DOUBLE );
797 _vtk->GetPointData( )->SetScalars( darray );
805 #endif // KGFO_USE_VTK
808 // -------------------------------------------------------------------------
809 void kVolume::deallocate( )
813 if( _vtk ) _vtk->Delete();
816 #endif // KGFO_USE_VTK
818 //#ifdef KGFO_USE_IDO
820 // if( _creator == SELF || _creator == VTK ) {
822 // delete[] _privateIdo;
823 // _privateIdo = NULL;
829 delete[] ( uchar* )_images;
831 //#endif // KGFO_USE_IDO
833 if( _raw && _creator == SELF )
835 //EED purify 12/sept/2006
836 // delete[] ( uchar* )_raw;
848 // -------------------------------------------------------------------------
849 kVolume::kVolume( vtkImageData* org )
850 : _raw( 0 ), _columns( 0 ), _images( 0 ),
852 //#ifdef KGFO_USE_IDO
853 // _privateIdo( NULL ),
854 //#endif // KGFO_USE_IDO
862 switch( org->GetScalarType( ) ) {
864 case VTK_CHAR: _type = CHAR; break;
865 case VTK_UNSIGNED_CHAR: _type = UCHAR; break;
866 case VTK_SHORT: _type = SHORT; break;
867 case VTK_INT: _type = INT; break;
868 case VTK_UNSIGNED_SHORT: _type = USHORT; break;
869 case VTK_UNSIGNED_INT: _type = UINT; break;
870 case VTK_FLOAT: _type = FLOAT; break;
871 case VTK_DOUBLE: _type = DOUBLE; break;
876 org->GetDimensions( itmp );
877 _dims[ CX ] = ( uint )itmp[ 0 ];
878 _dims[ CY ] = ( uint )itmp[ 1 ];
879 _dims[ CZ ] = ( uint )itmp[ 2 ];
880 org->GetSpacing( ftmp );
881 _sizes[ CX ] = ( double )ftmp[ 0 ];
882 _sizes[ CY ] = ( double )ftmp[ 1 ];
883 _sizes[ CZ ] = ( double )ftmp[ 2 ];
885 _raw = org->GetPointData( )->GetScalars( )->GetVoidPointer( 0 );
887 _vtk = vtkImageData::New();
888 _vtk->ShallowCopy( org );
892 // -------------------------------------------------------------------------
893 kVolume& kVolume::operator=( vtkImageData* org )
899 // -------------------------------------------------------------------------
900 void kVolume::copyFrom( vtkImageData* org )
910 //#ifdef KGFO_USE_IDO
912 // _privateIdo = NULL;
914 //#endif // KGFO_USE_IDO
922 switch( org->GetScalarType( ) ) {
924 case VTK_CHAR: _type = CHAR; break;
925 case VTK_UNSIGNED_CHAR: _type = UCHAR; break;
926 case VTK_SHORT: _type = SHORT; break;
927 case VTK_INT: _type = INT; break;
928 case VTK_UNSIGNED_SHORT: _type = USHORT; break;
929 case VTK_UNSIGNED_INT: _type = UINT; break;
930 case VTK_FLOAT: _type = FLOAT; break;
931 case VTK_DOUBLE: _type = DOUBLE; break;
936 org->GetDimensions( itmp );
937 _dims[ CX ] = ( uint )itmp[ 0 ];
938 _dims[ CY ] = ( uint )itmp[ 1 ];
939 _dims[ CZ ] = ( uint )itmp[ 2 ];
940 org->GetSpacing( ftmp );
941 _sizes[ CX ] = ( double )ftmp[ 0 ];
942 _sizes[ CY ] = ( double )ftmp[ 1 ];
943 _sizes[ CZ ] = ( double )ftmp[ 2 ];
948 // This avoids vtk extent crap conversion...
949 org->GetExtent( ext );
950 for( i = ext[ 0 ]; i <= ext[ 1 ]; i++ ) {
951 for( j = ext[ 2 ]; j <= ext[ 3 ]; j++ ) {
952 for( k = ext[ 4 ]; k <= ext[ 5 ]; k++ ) {
953 v = org->GetScalarComponentAsDouble( i, j, k, 0 );
954 setPixel( v, i - ext[ 0 ], j - ext[ 2 ], k - ext[ 4 ] );
960 #endif // KGFO_USE_VTK