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: volume.cxx,v $
31 Date: $Date: 2012/11/15 14:16:13 $
32 Version: $Revision: 1.10 $
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 =========================================================================*/
43 // PS -> //#include <gsl/gsl_math.h>//PS
49 #include <vtkCharArray.h>
50 #include <vtkDataArray.h>
51 #include <vtkDoubleArray.h>
52 #include <vtkFloatArray.h>
53 #include <vtkIntArray.h>
54 #include <vtkPointData.h>
55 #include <vtkShortArray.h>
56 #include <vtkUnsignedCharArray.h>
57 #include <vtkUnsignedIntArray.h>
58 #include <vtkUnsignedShortArray.h>
60 // -------------------------------------------------------------------------
61 const vtkIdType kVolume::VTKTypes[] = { VTK_CHAR, VTK_FLOAT, VTK_DOUBLE,
62 VTK_INT, VTK_SHORT, VTK_UNSIGNED_CHAR,
63 VTK_UNSIGNED_INT, VTK_UNSIGNED_SHORT };
65 #endif // KGFO_USE_VTK
67 // -------------------------------------------------------------------------
68 const void* kVolume::BLANK = ( void* ) 0;
69 const void* kVolume::NOALLOC = ( void* )-1;
70 const int kVolume::SIZETypes[] = { sizeof( char ), sizeof( float ),
71 sizeof( double ), sizeof( int ),
72 sizeof( short ), sizeof( uint8_t ),
73 sizeof( uint32_t ), sizeof( uint16_t ) };
75 // ---------------------------------------------------------------------------
76 template< class FROM, class TO >
77 static void convertCastT( FROM* src, TO* dest, ulong size )
79 FROM* end = src + size;
80 for( ; src < end; src++, dest++ ) *dest = ( TO )*src;
84 // ---------------------------------------------------------------------------
85 template< class FROM, class TO >
86 static void convertScaleT( FROM *src, TO *dest, ulong size,
87 double smin, double tmin, double slope )
89 FROM* end = src + size;
90 for( ; src < end; src++, dest++ )
91 *dest = ( TO )( ( double( *src ) - smin ) * slope + tmin );
94 // ---------------------------------------------------------------------------
95 template< class TYPE >
96 static void getMinMaxT( TYPE *src, ulong size,
97 double& min, double& max )
99 TYPE* end = src + size;
105 for( st = true; src < end; src++, st = false ) {
107 if( *src < m || st ) m = *src;
108 if( *src > M || st ) M = *src;
116 // -------------------------------------------------------------------------
126 #endif // KGFO_USE_VTK
128 _dims[ CX ] = 1; _dims[ CY ] = 1; _dims[ CZ ] = 1;
129 _sizes[ CX ] = 1; _sizes[ CY ] = 1; _sizes[ CZ ] = 1;
134 // -------------------------------------------------------------------------
135 kVolume::kVolume( Type type,
136 uint32_t xdim, uint32_t ydim, uint32_t zdim,
137 double xsize, double ysize, double zsize,
147 #endif // KGFO_USE_VTK
149 _dims[ CX ] = xdim; _dims[ CY ] = ydim; _dims[ CZ ] = zdim;
150 _sizes[ CX ] = xsize; _sizes[ CY ] = ysize; _sizes[ CZ ] = zsize;
151 if( data != NOALLOC ) {
162 // -------------------------------------------------------------------------
163 kVolume::kVolume( Type type,
164 const uint32_t *dims,
169 _raw( NULL ), _columns( NULL ),
174 #endif // KGFO_USE_VTK
177 memcpy( _dims, dims, 3 * sizeof( uint32_t ) );
178 memcpy( _sizes, sizes, 3 * sizeof( double ) );
179 if( data != NOALLOC ) {
190 // -------------------------------------------------------------------------
191 kVolume::kVolume( const kVolume& org )
200 #endif // KGFO_USE_VTK
206 // -------------------------------------------------------------------------
207 kVolume& kVolume::operator=( const kVolume& org )
213 // -------------------------------------------------------------------------
214 void kVolume::copyFrom( const kVolume& org )
222 memcpy( _dims, org._dims, 3 * sizeof( uint32_t ) );
223 memcpy( _sizes, org._sizes, 3 * sizeof( double ) );
228 memcpy( _raw, org._raw, getRawSizeInBytes( ) );
233 // -------------------------------------------------------------------------
234 bool kVolume::sameDimension( const kVolume& comp )
236 return( ( _dims[ CX ] == comp._dims[ CX ] &&
237 _dims[ CY ] == comp._dims[ CY ] &&
238 _dims[ CZ ] == comp._dims[ CZ ] ) );
241 // -------------------------------------------------------------------------
242 double kVolume::getPixel( uint32_t x, uint32_t y, uint32_t z ) const
248 case CHAR: p = ( double )( ( int8_t*** )_images )[ z ][ y ][ x ]; break;
249 case FLOAT: p = ( double )( ( float*** )_images )[ z ][ y ][ x ]; break;
250 case DOUBLE: p = ( double )( ( double*** )_images )[ z ][ y ][ x ]; break;
251 case INT: p = ( double )( ( int32_t*** )_images )[ z ][ y ][ x ]; break;
252 case SHORT: p = ( double )( ( int16_t*** )_images )[ z ][ y ][ x ]; break;
253 case UCHAR: p = ( double )( ( uint8_t*** )_images )[ z ][ y ][ x ]; break;
254 case UINT: p = ( double )( ( uint32_t*** )_images )[ z ][ y ][ x ]; break;
255 case USHORT: p = ( double )( ( uint16_t*** )_images )[ z ][ y ][ x ]; break;
256 default: p = 0.0; break;
263 // -------------------------------------------------------------------------
264 void kVolume::setPixel( double v, uint32_t x, uint32_t y, uint32_t z )
269 case CHAR: ( ( int8_t*** )_images )[ z ][ y ][ x ] = ( int8_t )v; break;
270 case FLOAT: ( ( float*** )_images )[ z ][ y ][ x ] = ( float )v; break;
271 case DOUBLE: ( ( double*** )_images )[ z ][ y ][ x ] = ( double )v; break;
272 case INT: ( ( int32_t*** )_images )[ z ][ y ][ x ] = ( int32_t )v; break;
273 case SHORT: ( ( int16_t*** )_images )[ z ][ y ][ x ] = ( int16_t )v; break;
274 case UCHAR: ( ( uint8_t*** )_images )[ z ][ y ][ x ] = ( uint8_t )v; break;
275 case UINT: ( ( uint32_t*** )_images )[ z ][ y ][ x ] = ( uint32_t )v; break;
276 case USHORT: ( ( uint16_t*** )_images )[ z ][ y ][ x ] = ( uint16_t )v; break;
282 // -------------------------------------------------------------------------
283 void kVolume::convertCast( Type type )
285 if( _type != type ) {
288 ulong size = getRawSize( );
290 buffer = ( void* )new uint8_t[ size * SIZETypes[ type ] ];
297 case SHORT: convertCastT( ( char* )_raw, ( int16_t* )buffer, size ); break;
298 case INT: convertCastT( ( char* )_raw, ( int32_t* )buffer, size ); break;
299 case USHORT: convertCastT( ( char* )_raw, ( uint16_t* )buffer, size ); break;
300 case UINT: convertCastT( ( char* )_raw, ( uint32_t* )buffer, size ); break;
301 case FLOAT: convertCastT( ( char* )_raw, ( float* )buffer, size ); break;
302 case DOUBLE: convertCastT( ( char* )_raw, ( double* )buffer, size ); break;
303 case UCHAR: convertCastT( ( char* )_raw, ( uint8_t* )buffer, size ); break;
312 case CHAR: convertCastT( ( int16_t* )_raw, ( uint8_t* )buffer, size ); break;
313 case INT: convertCastT( ( int16_t* )_raw, ( int32_t* )buffer, size ); break;
314 case USHORT: convertCastT( ( int16_t* )_raw, ( uint16_t* )buffer, size ); break;
315 case UINT: convertCastT( ( int16_t* )_raw, ( uint32_t* )buffer, size ); break;
316 case FLOAT: convertCastT( ( int16_t* )_raw, ( float* )buffer, size ); break;
317 case DOUBLE: convertCastT( ( int16_t* )_raw, ( double* )buffer, size ); break;
318 case UCHAR: convertCastT( ( int16_t* )_raw, ( uint8_t* )buffer, size ); break;
326 case CHAR: convertCastT( ( int32_t* )_raw, ( int8_t* )buffer, size ); break;
327 case SHORT: convertCastT( ( int32_t* )_raw, ( int16_t* )buffer, size ); break;
328 case USHORT: convertCastT( ( int32_t* )_raw, ( uint16_t* )buffer, size ); break;
329 case UINT: convertCastT( ( int32_t* )_raw, ( uint32_t* )buffer, size ); break;
330 case FLOAT: convertCastT( ( int32_t* )_raw, ( float* )buffer, size ); break;
331 case DOUBLE: convertCastT( ( int32_t* )_raw, ( double* )buffer, size ); break;
332 case UCHAR: convertCastT( ( int32_t* )_raw, ( uint8_t* )buffer, size ); break;
340 case CHAR: convertCastT( ( uint16_t* )_raw, ( int8_t* )buffer, size ); break;
341 case SHORT: convertCastT( ( uint16_t* )_raw, ( int16_t* )buffer, size ); break;
342 case INT: convertCastT( ( uint16_t* )_raw, ( int32_t* )buffer, size ); break;
343 case UINT: convertCastT( ( uint16_t* )_raw, ( uint32_t* )buffer, size ); break;
344 case FLOAT: convertCastT( ( uint16_t* )_raw, ( float* )buffer, size ); break;
345 case DOUBLE: convertCastT( ( uint16_t* )_raw, ( double* )buffer, size ); break;
346 case UCHAR: convertCastT( ( uint16_t* )_raw, ( uint8_t* )buffer, size ); break;
354 case CHAR: convertCastT( ( uint32_t* )_raw, ( int8_t* )buffer, size ); break;
355 case SHORT: convertCastT( ( uint32_t* )_raw, ( int16_t* )buffer, size ); break;
356 case INT: convertCastT( ( uint32_t* )_raw, ( int32_t* )buffer, size ); break;
357 case USHORT: convertCastT( ( uint32_t* )_raw, ( uint16_t* )buffer, size ); break;
358 case FLOAT: convertCastT( ( uint32_t* )_raw, ( float* )buffer, size ); break;
359 case DOUBLE: convertCastT( ( uint32_t* )_raw, ( double* )buffer, size ); break;
360 case UCHAR: convertCastT( ( uint32_t* )_raw, ( uint8_t* )buffer, size ); break;
368 case CHAR: convertCastT( ( float* )_raw, ( int8_t* )buffer, size ); break;
369 case SHORT: convertCastT( ( float* )_raw, ( int16_t* )buffer, size ); break;
370 case INT: convertCastT( ( float* )_raw, ( int32_t* )buffer, size ); break;
371 case USHORT: convertCastT( ( float* )_raw, ( uint16_t* )buffer, size ); break;
372 case UINT: convertCastT( ( float* )_raw, ( uint32_t* )buffer, size ); break;
373 case DOUBLE: convertCastT( ( float* )_raw, ( double* )buffer, size ); break;
374 case UCHAR: convertCastT( ( float* )_raw, ( uint8_t* )buffer, size ); break;
382 case CHAR: convertCastT( ( double* )_raw, ( int8_t* )buffer, size ); break;
383 case SHORT: convertCastT( ( double* )_raw, ( int16_t* )buffer, size ); break;
384 case INT: convertCastT( ( double* )_raw, ( int32_t* )buffer, size ); break;
385 case USHORT: convertCastT( ( double* )_raw, ( uint16_t* )buffer, size ); break;
386 case UINT: convertCastT( ( double* )_raw, ( uint32_t* )buffer, size ); break;
387 case FLOAT: convertCastT( ( double* )_raw, ( float* )buffer, size ); break;
388 case UCHAR: convertCastT( ( double* )_raw, ( uint8_t* )buffer, size ); break;
396 case CHAR: convertCastT( ( uint8_t* )_raw, ( int8_t* )buffer, size ); break;
397 case SHORT: convertCastT( ( uint8_t* )_raw, ( int16_t* )buffer, size ); break;
398 case INT: convertCastT( ( uint8_t* )_raw, ( int32_t* )buffer, size ); break;
399 case USHORT: convertCastT( ( uint8_t* )_raw, ( uint16_t* )buffer, size ); break;
400 case UINT: convertCastT( ( uint8_t* )_raw, ( uint32_t* )buffer, size ); break;
401 case FLOAT: convertCastT( ( uint8_t* )_raw, ( float* )buffer, size ); break;
402 case DOUBLE: convertCastT( ( uint8_t* )_raw, ( double* )buffer, size ); break;
418 // -------------------------------------------------------------------------
419 void kVolume::convertScale( Type type, double min, double max )
421 double tmin, tmax, smin, smax;
423 // PS -> //tmin = GSL_MIN( min, max ); //PS
424 tmin= ((min<max)?min:max);
425 // PS -> //tmax = GSL_MAX( min, max ); //PS
426 tmax= ((min>max)?min:max);
428 getMinMax( smin, smax );
430 // compute scaling slope
431 double a = ( tmax - tmin ) / ( smax - smin );
434 ulong size = getRawSize( );
436 buffer = ( void* )new uint8_t[ size * SIZETypes[ type ] ];
443 case CHAR: convertScaleT( ( int8_t* )_raw, ( int8_t* )buffer, size, smin, tmin, a ); break;
444 case SHORT: convertScaleT( ( int8_t* )_raw, ( int16_t* )buffer, size, smin, tmin, a ); break;
445 case INT: convertScaleT( ( int8_t* )_raw, ( int32_t* )buffer, size, smin, tmin, a ); break;
446 case USHORT: convertScaleT( ( int8_t* )_raw, ( uint16_t* )buffer, size, smin, tmin, a ); break;
447 case UINT: convertScaleT( ( int8_t* )_raw, ( uint32_t* )buffer, size, smin, tmin, a ); break;
448 case FLOAT: convertScaleT( ( int8_t* )_raw, ( float* )buffer, size, smin, tmin, a ); break;
449 case DOUBLE: convertScaleT( ( int8_t* )_raw, ( double* )buffer, size, smin, tmin, a ); break;
450 case UCHAR: convertScaleT( ( int8_t* )_raw, ( uint8_t* )buffer, size, smin, tmin, a ); break;
457 case CHAR: convertScaleT( ( int16_t* )_raw, ( int8_t* )buffer, size, smin, tmin, a ); break;
458 case SHORT: convertScaleT( ( int16_t* )_raw, ( int16_t* )buffer, size, smin, tmin, a ); break;
459 case INT: convertScaleT( ( int16_t* )_raw, ( int32_t* )buffer, size, smin, tmin, a ); break;
460 case USHORT: convertScaleT( ( int16_t* )_raw, ( uint16_t* )buffer, size, smin, tmin, a ); break;
461 case UINT: convertScaleT( ( int16_t* )_raw, ( uint32_t* )buffer, size, smin, tmin, a ); break;
462 case FLOAT: convertScaleT( ( int16_t* )_raw, ( float* )buffer, size, smin, tmin, a ); break;
463 case DOUBLE: convertScaleT( ( int16_t* )_raw, ( double* )buffer, size, smin, tmin, a ); break;
464 case UCHAR: convertScaleT( ( int16_t* )_raw, ( uint8_t* )buffer, size, smin, tmin, a ); break;
471 case CHAR: convertScaleT( ( int32_t* )_raw, ( int8_t* )buffer, size, smin, tmin, a ); break;
472 case SHORT: convertScaleT( ( int32_t* )_raw, ( int16_t* )buffer, size, smin, tmin, a ); break;
473 case INT: convertScaleT( ( int32_t* )_raw, ( int32_t* )buffer, size, smin, tmin, a ); break;
474 case USHORT: convertScaleT( ( int32_t* )_raw, ( uint16_t* )buffer, size, smin, tmin, a ); break;
475 case UINT: convertScaleT( ( int32_t* )_raw, ( uint32_t* )buffer, size, smin, tmin, a ); break;
476 case FLOAT: convertScaleT( ( int32_t* )_raw, ( float* )buffer, size, smin, tmin, a ); break;
477 case DOUBLE: convertScaleT( ( int32_t* )_raw, ( double* )buffer, size, smin, tmin, a ); break;
478 case UCHAR: convertScaleT( ( int32_t* )_raw, ( uint8_t* )buffer, size, smin, tmin, a ); break;
485 case CHAR: convertScaleT( ( uint16_t* )_raw, ( int8_t* )buffer, size, smin, tmin, a ); break;
486 case SHORT: convertScaleT( ( uint16_t* )_raw, ( int16_t* )buffer, size, smin, tmin, a ); break;
487 case INT: convertScaleT( ( uint16_t* )_raw, ( int32_t* )buffer, size, smin, tmin, a ); break;
488 case USHORT: convertScaleT( ( uint16_t* )_raw, ( uint16_t* )buffer, size, smin, tmin, a ); break;
489 case UINT: convertScaleT( ( uint16_t* )_raw, ( uint32_t* )buffer, size, smin, tmin, a ); break;
490 case FLOAT: convertScaleT( ( uint16_t* )_raw, ( float* )buffer, size, smin, tmin, a ); break;
491 case DOUBLE: convertScaleT( ( uint16_t* )_raw, ( double* )buffer, size, smin, tmin, a ); break;
492 case UCHAR: convertScaleT( ( uint16_t* )_raw, ( uint8_t* )buffer, size, smin, tmin, a ); break;
499 case CHAR: convertScaleT( ( uint32_t* )_raw, ( int8_t* )buffer, size, smin, tmin, a ); break;
500 case SHORT: convertScaleT( ( uint32_t* )_raw, ( int16_t* )buffer, size, smin, tmin, a ); break;
501 case INT: convertScaleT( ( uint32_t* )_raw, ( int32_t* )buffer, size, smin, tmin, a ); break;
502 case USHORT: convertScaleT( ( uint32_t* )_raw, ( uint16_t* )buffer, size, smin, tmin, a ); break;
503 case UINT: convertScaleT( ( uint32_t* )_raw, ( uint32_t* )buffer, size, smin, tmin, a ); break;
504 case FLOAT: convertScaleT( ( uint32_t* )_raw, ( float* )buffer, size, smin, tmin, a ); break;
505 case DOUBLE: convertScaleT( ( uint32_t* )_raw, ( double* )buffer, size, smin, tmin, a ); break;
506 case UCHAR: convertScaleT( ( uint32_t* )_raw, ( uint8_t* )buffer, size, smin, tmin, a ); break;
513 case CHAR: convertScaleT( ( uint8_t* )_raw, ( int8_t* )buffer, size, smin, tmin, a ); break;
514 case SHORT: convertScaleT( ( uint8_t* )_raw, ( int16_t* )buffer, size, smin, tmin, a ); break;
515 case INT: convertScaleT( ( uint8_t* )_raw, ( int32_t* )buffer, size, smin, tmin, a ); break;
516 case USHORT: convertScaleT( ( uint8_t* )_raw, ( uint16_t* )buffer, size, smin, tmin, a ); break;
517 case UINT: convertScaleT( ( uint8_t* )_raw, ( uint32_t* )buffer, size, smin, tmin, a ); break;
518 case FLOAT: convertScaleT( ( uint8_t* )_raw, ( float* )buffer, size, smin, tmin, a ); break;
519 case DOUBLE: convertScaleT( ( uint8_t* )_raw, ( double* )buffer, size, smin, tmin, a ); break;
520 case UCHAR: convertScaleT( ( uint8_t* )_raw, ( uint8_t* )buffer, size, smin, tmin, a ); break;
527 case CHAR: convertScaleT( ( double* )_raw, ( int8_t* )buffer, size, smin, tmin, a ); break;
528 case SHORT: convertScaleT( ( double* )_raw, ( int16_t* )buffer, size, smin, tmin, a ); break;
529 case INT: convertScaleT( ( double* )_raw, ( int32_t* )buffer, size, smin, tmin, a ); break;
530 case USHORT: convertScaleT( ( double* )_raw, ( uint16_t* )buffer, size, smin, tmin, a ); break;
531 case UINT: convertScaleT( ( double* )_raw, ( uint32_t* )buffer, size, smin, tmin, a ); break;
532 case FLOAT: convertScaleT( ( double* )_raw, ( float* )buffer, size, smin, tmin, a ); break;
533 case DOUBLE: convertScaleT( ( double* )_raw, ( double* )buffer, size, smin, tmin, a ); break;
534 case UCHAR: convertScaleT( ( double* )_raw, ( uint8_t* )buffer, size, smin, tmin, a ); break;
541 case CHAR: convertScaleT( ( float* )_raw, ( int8_t* )buffer, size, smin, tmin, a ); break;
542 case SHORT: convertScaleT( ( float* )_raw, ( int16_t* )buffer, size, smin, tmin, a ); break;
543 case INT: convertScaleT( ( float* )_raw, ( int32_t* )buffer, size, smin, tmin, a ); break;
544 case USHORT: convertScaleT( ( float* )_raw, ( uint16_t* )buffer, size, smin, tmin, a ); break;
545 case UINT: convertScaleT( ( float* )_raw, ( uint32_t* )buffer, size, smin, tmin, a ); break;
546 case FLOAT: convertScaleT( ( float* )_raw, ( float* )buffer, size, smin, tmin, a ); break;
547 case DOUBLE: convertScaleT( ( float* )_raw, ( double* )buffer, size, smin, tmin, a ); break;
548 case UCHAR: convertScaleT( ( float* )_raw, ( uint8_t* )buffer, size, smin, tmin, a ); break;
562 // -------------------------------------------------------------------------
563 void kVolume::getMinMax( double& min, double& max ) const
565 ulong size = getRawSize( );
569 case CHAR: getMinMaxT( ( int8_t* )_raw, size, min, max ); break;
570 case FLOAT: getMinMaxT( ( float* )_raw, size, min, max ); break;
571 case DOUBLE: getMinMaxT( ( double* )_raw, size, min, max ); break;
572 case INT: getMinMaxT( ( int32_t* )_raw, size, min, max ); break;
573 case SHORT: getMinMaxT( ( int16_t* )_raw, size, min, max ); break;
574 case UCHAR: getMinMaxT( ( uint8_t* )_raw, size, min, max ); break;
575 case UINT: getMinMaxT( ( uint32_t* )_raw, size, min, max ); break;
576 case USHORT: getMinMaxT( ( uint16_t* )_raw, size, min, max ); break;
581 // -------------------------------------------------------------------------
582 double kVolume::getMin( ) const
590 // -------------------------------------------------------------------------
591 double kVolume::getMax( ) const
599 // -------------------------------------------------------------------------
600 /*double kVolume::GetMaxIntSphere( double* p, double r )
602 int minX, minY, minZ, maxX, maxY, maxZ;
603 gslobj_vector vP( p, 3 ), v( 3 );
607 minX = int( floor( p[ 0 ] - r ) );
608 minY = int( floor( p[ 1 ] - r ) );
609 minZ = int( floor( p[ 2 ] - r ) );
610 maxX = int( ceil( p[ 0 ] + r ) );
611 maxY = int( ceil( p[ 1 ] + r ) );
612 maxZ = int( ceil( p[ 2 ] + r ) );
614 minX = GSL_MAX( minX, 0 );
615 minY = GSL_MAX( minY, 0 );
616 minZ = GSL_MAX( minZ, 0 );
618 maxX = GSL_MIN( maxX, this->getXdim( ) - 1 );
619 maxY = GSL_MIN( maxY, this->getYdim( ) - 1 );
620 maxZ = GSL_MIN( maxZ, this->getZdim( ) - 1 );
624 for( v( 0 ) = minX; v( 0 ) <= maxX; v( 0 )++ )
625 for( v( 1 ) = minY; v( 1 ) <= maxY; v( 1 )++ )
626 for( v( 2 ) = minZ; v( 2 ) <= maxZ; v( 2 )++ )
627 if( ( v - vP ).norm2( ) <= r ) {
628 tmp = this->getPixel( ( uint32_t )v(0), ( uint32_t )v(1), ( uint32_t )v(2));
629 maxint = ( tmp > maxint || start )? tmp: maxint;
635 // -------------------------------------------------------------------------
636 unsigned short kVolume::GetMaxIntSphere2( double* p, double r )
639 * unsigned short range : 0 -> 65535
640 * unsigned int range : 0 -> 2147483647 // 2^31 - 1
642 int minX, minY, minZ, maxX, maxY, maxZ;
644 unsigned short maxint = 0, tmp;
646 minX = int( floor( p[ 0 ] - r ) );
647 minY = int( floor( p[ 1 ] - r ) );
648 minZ = int( floor( p[ 2 ] - r ) );
649 maxX = int( ceil( p[ 0 ] + r ) );
650 maxY = int( ceil( p[ 1 ] + r ) );
651 maxZ = int( ceil( p[ 2 ] + r ) );
653 // PS -> //minX = GSL_MAX( minX, 0 );//PS
654 // PS -> //minY = GSL_MAX( minY, 0 );//PS
655 // PS -> //minZ = GSL_MAX( minZ, 0 );//PS
656 minX=((minX>0)?minX:0);
657 minY=((minY>0)?minY:0);
658 minZ=((minZ>0)?minZ:0);
660 // PS -> //maxX = GSL_MIN( maxX, this->getXdim( ) - 1 );//PS
661 // PS -> //maxY = GSL_MIN( maxY, this->getYdim( ) - 1 );//PS
662 // PS -> //maxZ = GSL_MIN( maxZ, this->getZdim( ) - 1 );//PS
663 int xDim=this->getXdim( ) - 1;
664 maxX= (maxX<xDim?maxX:xDim);
666 int yDim=this->getYdim( ) - 1;
667 maxY= (maxY<yDim?maxY:yDim);
669 int zDim=this->getZdim( ) - 1;
670 maxZ= (maxZ<zDim?maxZ:zDim);
672 double r2 = r*r; //need to do comparison in double ... don't know why...
673 for( v[0] = minX; v[0] <= maxX; v[0]++ )
674 for( v[1] = minY; v[1] <= maxY; v[1]++ )
675 for( v[2] = minZ; v[2] <= maxZ; v[2]++ )
676 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 )
678 tmp = (unsigned short)this->getPixel( v[0], v[1], v[2] );
679 maxint = ( tmp > maxint ) ? tmp : maxint;
685 // -------------------------------------------------------------------------
686 void kVolume::allocate( )
688 ulong size = getRawSizeInBytes( );
690 if( _creator == SELF ) {
692 _raw = ( void* )new uint8_t[ size ];
693 memset( _raw, 0, size );
698 // -------------------------------------------------------------------------
699 void kVolume::buildIndex( )
703 size = ( _dims[ CZ ] * sizeof( uint8_t** ) ) +
704 ( _dims[ CZ ] * _dims[ CY ] * sizeof( void* ) );
706 _images = ( void*** )new uint8_t[ size ];
708 _columns = ( void** )( _images + _dims[ CZ ] );
709 void** plane = _columns;
710 for( uint32_t z = 0; z < _dims[ CZ ]; z++ ) {
712 _images[ z ] = plane;
713 plane += _dims[ CY ];
717 for( uint32_t y = 0; y < _dims[ CZ ] * _dims[ CY ]; y++ ) {
719 _columns[ y ] = line;
720 line = ( void* )( ( uint8_t* ) line +
721 _dims[ CX ] * SIZETypes[ _type ] );
727 vtkCharArray *carray;
728 vtkDoubleArray *darray;
729 vtkFloatArray *farray;
731 vtkShortArray *sarray;
732 vtkUnsignedCharArray *ucarray;
733 vtkUnsignedIntArray *uiarray;
734 vtkUnsignedShortArray *usarray;
736 size = _dims[ CX ] * _dims[ CY ] * _dims[ CZ ];
738 if( _creator == SELF || _creator == IDO ) {
740 _vtk = vtkImageData::New( );
741 _vtk->SetDimensions( _dims[ CX ], _dims[ CY ], _dims[ CZ ] );
742 _vtk->SetSpacing( _sizes[ CX ], _sizes[ CY ], _sizes[ CZ ] );
746 carray = vtkCharArray::New( );
747 carray->SetArray( ( char* )( _raw ), size, 1 );
748 //EED 2017-01-01 Migration VTK7
749 #if VTK_MAJOR_VERSION <= 5
750 carray->SetNumberOfComponents( 1 );
751 _vtk->SetScalarType( VTK_CHAR );
753 vtkInformation* infoC=_vtk->GetInformation();
754 vtkDataObject::SetPointDataActiveScalarInfo(infoC, VTK_CHAR, 1);
756 _vtk->GetPointData( )->SetScalars( carray );
762 ucarray = vtkUnsignedCharArray::New( );
763 ucarray->SetArray( ( uint8_t* )( _raw ), size, 1 );
764 //EED 2017-01-01 Migration VTK7
765 #if VTK_MAJOR_VERSION <= 5
766 ucarray->SetNumberOfComponents( 1 );
767 _vtk->SetScalarType( VTK_UNSIGNED_CHAR );
769 vtkInformation* infoUC=_vtk->GetInformation();
770 vtkDataObject::SetPointDataActiveScalarInfo(infoUC, VTK_UNSIGNED_CHAR, 1);
772 _vtk->GetPointData( )->SetScalars( ucarray );
782 sarray = vtkShortArray::New( );
783 sarray->SetArray( ( int16_t* )( _raw ), size, 1 );
784 //EED 2017-01-01 Migration VTK7
785 #if VTK_MAJOR_VERSION <= 5
786 sarray->SetNumberOfComponents( 1 );
787 _vtk->SetScalarType( VTK_SHORT );
789 vtkInformation* infoS=_vtk->GetInformation();
790 vtkDataObject::SetPointDataActiveScalarInfo(infoS, VTK_SHORT, 1);
792 _vtk->GetPointData( )->SetScalars( sarray );
798 iarray = vtkIntArray::New( );
799 iarray->SetArray( ( int32_t* )( _raw ), size, 1 );
800 //EED 2017-01-01 Migration VTK7
801 #if VTK_MAJOR_VERSION <= 5
802 iarray->SetNumberOfComponents( 1 );
803 _vtk->SetScalarType( VTK_INT );
805 vtkInformation* infoI=_vtk->GetInformation();
806 vtkDataObject::SetPointDataActiveScalarInfo(infoI, VTK_INT, 1);
808 _vtk->GetPointData( )->SetScalars( iarray );
815 usarray = vtkUnsignedShortArray::New( );
816 usarray->SetArray( ( uint16_t* )( _raw ), size, 1 );
817 //EED 2017-01-01 Migration VTK7
818 #if VTK_MAJOR_VERSION <= 5
819 usarray->SetNumberOfComponents( 1 );
820 _vtk->SetScalarType( VTK_UNSIGNED_SHORT );
822 vtkInformation* infoUS=_vtk->GetInformation();
823 vtkDataObject::SetPointDataActiveScalarInfo(infoUS, VTK_UNSIGNED_SHORT, 1);
825 _vtk->GetPointData( )->SetScalars( usarray );
832 uiarray = vtkUnsignedIntArray::New( );
833 uiarray->SetArray( ( uint32_t* )( _raw ), size, 1 );
834 //EED 2017-01-01 Migration VTK7
835 #if VTK_MAJOR_VERSION <= 5
836 uiarray->SetNumberOfComponents( 1 );
837 _vtk->SetScalarType( VTK_UNSIGNED_INT );
839 vtkInformation* infoUI=_vtk->GetInformation();
840 vtkDataObject::SetPointDataActiveScalarInfo(infoUI, VTK_UNSIGNED_INT, 1);
842 _vtk->GetPointData( )->SetScalars( uiarray );
849 farray = vtkFloatArray::New( );
850 farray->SetArray( ( float* )( _raw ), size, 1 );
851 //EED 2017-01-01 Migration VTK7
852 #if VTK_MAJOR_VERSION <= 5
853 farray->SetNumberOfComponents( 1 );
854 _vtk->SetScalarType( VTK_FLOAT );
856 vtkInformation* infoF=_vtk->GetInformation();
857 vtkDataObject::SetPointDataActiveScalarInfo(infoF, VTK_FLOAT, 1);
859 _vtk->GetPointData( )->SetScalars( farray );
866 darray = vtkDoubleArray::New( );
867 darray->SetArray( ( double* )( _raw ), size, 1 );
868 //EED 2017-01-01 Migration VTK7
869 #if VTK_MAJOR_VERSION <= 5
870 darray->SetNumberOfComponents( 1 );
871 _vtk->SetScalarType( VTK_DOUBLE );
873 vtkInformation* infoD=_vtk->GetInformation();
874 vtkDataObject::SetPointDataActiveScalarInfo(infoD, VTK_DOUBLE, 1);
878 _vtk->GetPointData( )->SetScalars( darray );
891 carray = vtkCharArray::New( );
892 carray->SetArray( ( char* )( _raw ), size, 1 );
894 //EED 2017-01-01 Migration VTK7
895 //#if VTK_MAJOR_VERSION <= 5
896 carray->SetNumberOfComponents( 1 );
897 _vtk->SetScalarType( VTK_CHAR );
899 vtkInformation* infoC=_vtk->GetInformation();
900 vtkDataObject::SetPointDataActiveScalarInfo(infoC, VTK_CHAR, 1);
903 _vtk->GetPointData( )->SetScalars( carray );
908 ucarray = vtkUnsignedCharArray::New( );
909 ucarray->SetArray( ( uint8_t* )( _raw ), size, 1 );
910 //EED 2017-01-01 Migration VTK7
911 #if VTK_MAJOR_VERSION <= 5
912 ucarray->SetNumberOfComponents( 1 );
913 _vtk->SetScalarType( VTK_UNSIGNED_CHAR );
915 vtkInformation* infoUC=_vtk->GetInformation();
916 vtkDataObject::SetPointDataActiveScalarInfo(infoUC, VTK_UNSIGNED_CHAR, 1);
918 _vtk->GetPointData( )->SetScalars( ucarray );
923 sarray = vtkShortArray::New( );
924 sarray->SetArray( ( int16_t* )( _raw ), size, 1 );
925 //EED 2017-01-01 Migration VTK7
926 #if VTK_MAJOR_VERSION <= 5
927 sarray->SetNumberOfComponents( 1 );
928 _vtk->SetScalarType( VTK_SHORT );
930 vtkInformation* infoS=_vtk->GetInformation();
931 vtkDataObject::SetPointDataActiveScalarInfo(infoS, VTK_SHORT, 1);
933 _vtk->GetPointData( )->SetScalars( sarray );
938 iarray = vtkIntArray::New( );
939 iarray->SetArray( ( int32_t* )( _raw ), size, 1 );
940 //EED 2017-01-01 Migration VTK7
941 #if VTK_MAJOR_VERSION <= 5
942 iarray->SetNumberOfComponents( 1 );
943 _vtk->SetScalarType( VTK_INT );
945 vtkInformation* infoI=_vtk->GetInformation();
946 vtkDataObject::SetPointDataActiveScalarInfo(infoI, VTK_INT, 1);
948 _vtk->GetPointData( )->SetScalars( iarray );
953 usarray = vtkUnsignedShortArray::New( );
954 usarray->SetArray( ( uint16_t* )( _raw ), size, 1 );
955 //EED 2017-01-01 Migration VTK7
956 #if VTK_MAJOR_VERSION <= 5
957 usarray->SetNumberOfComponents( 1 );
958 _vtk->SetScalarType( VTK_UNSIGNED_SHORT );
960 vtkInformation* infoUS=_vtk->GetInformation();
961 vtkDataObject::SetPointDataActiveScalarInfo(infoUS, VTK_UNSIGNED_SHORT, 1);
963 _vtk->GetPointData( )->SetScalars( usarray );
968 uiarray = vtkUnsignedIntArray::New( );
969 uiarray->SetArray( ( uint32_t* )( _raw ), size, 1 );
970 //EED 2017-01-01 Migration VTK7
971 #if VTK_MAJOR_VERSION <= 5
972 uiarray->SetNumberOfComponents( 1 );
973 _vtk->SetScalarType( VTK_UNSIGNED_INT );
975 vtkInformation* infoUI=_vtk->GetInformation();
976 vtkDataObject::SetPointDataActiveScalarInfo(infoUI, VTK_UNSIGNED_INT, 1);
978 _vtk->GetPointData( )->SetScalars( uiarray );
983 farray = vtkFloatArray::New( );
984 farray->SetArray( ( float* )( _raw ), size, 1 );
985 //EED 2017-01-01 Migration VTK7
986 #if VTK_MAJOR_VERSION <= 5
987 farray->SetNumberOfComponents( 1 );
988 _vtk->SetScalarType( VTK_FLOAT );
990 vtkInformation* infoF=_vtk->GetInformation();
991 vtkDataObject::SetPointDataActiveScalarInfo(infoF, VTK_FLOAT, 1);
993 _vtk->GetPointData( )->SetScalars( farray );
998 darray = vtkDoubleArray::New( );
999 darray->SetArray( ( double* )( _raw ), size, 1 );
1000 //EED 2017-01-01 Migration VTK7
1001 #if VTK_MAJOR_VERSION <= 5
1002 darray->SetNumberOfComponents( 1 );
1003 _vtk->SetScalarType( VTK_DOUBLE );
1005 vtkInformation* infoD=_vtk->GetInformation();
1006 vtkDataObject::SetPointDataActiveScalarInfo(infoD, VTK_DOUBLE, 1);
1010 _vtk->GetPointData( )->SetScalars( darray );
1019 #endif // KGFO_USE_VTK
1022 // -------------------------------------------------------------------------
1023 void kVolume::deallocate( )
1026 if( _vtk ) _vtk->Delete();
1028 #endif // KGFO_USE_VTK
1029 delete[] ( uint8_t* )_images;
1030 if( _raw && _creator == SELF )
1032 //EED purify 12/sept/2006
1033 // delete[] ( uint8_t* )_raw;
1045 // -------------------------------------------------------------------------
1046 kVolume::kVolume( vtkImageData* org )
1059 switch( org->GetScalarType( ) ) {
1061 case VTK_CHAR: _type = CHAR; break;
1062 case VTK_UNSIGNED_CHAR: _type = UCHAR; break;
1063 case VTK_SHORT: _type = SHORT; break;
1064 case VTK_INT: _type = INT; break;
1065 case VTK_UNSIGNED_SHORT: _type = USHORT; break;
1066 case VTK_UNSIGNED_INT: _type = UINT; break;
1067 case VTK_FLOAT: _type = FLOAT; break;
1068 case VTK_DOUBLE: _type = DOUBLE; break;
1073 org->GetDimensions( itmp );
1074 _dims[ CX ] = ( uint32_t )itmp[ 0 ];
1075 _dims[ CY ] = ( uint32_t )itmp[ 1 ];
1076 _dims[ CZ ] = ( uint32_t )itmp[ 2 ];
1077 org->GetSpacing( ftmp );
1078 _sizes[ CX ] = ( double )ftmp[ 0 ];
1079 _sizes[ CY ] = ( double )ftmp[ 1 ];
1080 _sizes[ CZ ] = ( double )ftmp[ 2 ];
1082 _raw = org->GetPointData( )->GetScalars( )->GetVoidPointer( 0 );
1084 _vtk = vtkImageData::New();
1085 _vtk->ShallowCopy( org );
1089 // -------------------------------------------------------------------------
1090 kVolume& kVolume::operator=( vtkImageData* org )
1096 // -------------------------------------------------------------------------
1097 void kVolume::copyFrom( vtkImageData* org )
1107 //#ifdef KGFO_USE_IDO
1109 // _privateIdo = NULL;
1111 //#endif // KGFO_USE_IDO
1119 switch( org->GetScalarType( ) ) {
1121 case VTK_CHAR: _type = CHAR; break;
1122 case VTK_UNSIGNED_CHAR: _type = UCHAR; break;
1123 case VTK_SHORT: _type = SHORT; break;
1124 case VTK_INT: _type = INT; break;
1125 case VTK_UNSIGNED_SHORT: _type = USHORT; break;
1126 case VTK_UNSIGNED_INT: _type = UINT; break;
1127 case VTK_FLOAT: _type = FLOAT; break;
1128 case VTK_DOUBLE: _type = DOUBLE; break;
1133 org->GetDimensions( itmp );
1134 _dims[ CX ] = ( uint32_t )itmp[ 0 ];
1135 _dims[ CY ] = ( uint32_t )itmp[ 1 ];
1136 _dims[ CZ ] = ( uint32_t )itmp[ 2 ];
1137 org->GetSpacing( ftmp );
1138 _sizes[ CX ] = ( double )ftmp[ 0 ];
1139 _sizes[ CY ] = ( double )ftmp[ 1 ];
1140 _sizes[ CZ ] = ( double )ftmp[ 2 ];
1145 // This avoids vtk extent crap conversion...
1146 org->GetExtent( ext );
1147 for( i = ext[ 0 ]; i <= ext[ 1 ]; i++ ) {
1148 for( j = ext[ 2 ]; j <= ext[ 3 ]; j++ ) {
1149 for( k = ext[ 4 ]; k <= ext[ 5 ]; k++ ) {
1150 v = org->GetScalarComponentAsDouble( i, j, k, 0 );
1151 setPixel( v, i - ext[ 0 ], j - ext[ 2 ], k - ext[ 4 ] );
1157 #endif // KGFO_USE_VTK