BICPL: create_slice.c Source File

00001 /* ---------------------------------------------------------------------------- 00002 @COPYRIGHT : 00003 Copyright 1993,1994,1995 David MacDonald, 00004 McConnell Brain Imaging Centre, 00005 Montreal Neurological Institute, McGill University. 00006 Permission to use, copy, modify, and distribute this 00007 software and its documentation for any purpose and without 00008 fee is hereby granted, provided that the above copyright 00009 notice appear in all copies. The author and McGill University 00010 make no representations about the suitability of this 00011 software for any purpose. It is provided "as is" without 00012 express or implied warranty. 00013 ---------------------------------------------------------------------------- */ 00014 00015 #include <volume_io/internal_volume_io.h> 00016 #include <bicpl/vols.h> 00017 00018 #ifndef lint 00019 static char rcsid[] = "$Header: /software/source//libraries/bicpl/Volumes/create_slice.c,v 1.45 2000/02/06 15:30:54 stever Exp $"; 00020 #endif 00021 00022 private void create_pixel_mapping( 00023 Volume volume1, 00024 int n_slices1, 00025 Real **origins1, 00026 Real x_axis1[], 00027 Real y_axis1[], 00028 Real x_translation1, 00029 Real y_translation1, 00030 Real x_scale1, 00031 Real y_scale1, 00032 Volume volume2, 00033 int n_slices2, 00034 Real **origins2, 00035 Real x_axis2[], 00036 Real y_axis2[], 00037 Real x_translation2, 00038 Real y_translation2, 00039 Real x_scale2, 00040 Real y_scale2, 00041 Real real_x_axis1[], 00042 Real real_y_axis1[], 00043 Real ***real_origins1, 00044 Real real_x_axis2[], 00045 Real real_y_axis2[], 00046 Real ***real_origins2 ) 00047 { 00048 int n_dimensions1, n_dimensions2, s; 00049 00050 n_dimensions1 = get_volume_n_dimensions( volume1 ); 00051 00052 ALLOC2D( *real_origins1, n_slices1, n_dimensions1 ); 00053 00054 for_less( s, 0, n_slices1 ) 00055 { 00056 get_mapping( volume1, origins1[s], x_axis1, y_axis1, 00057 x_translation1, y_translation1, x_scale1, y_scale1, 00058 (*real_origins1)[s], real_x_axis1, real_y_axis1 ); 00059 } 00060 00061 if( volume2 != NULL ) 00062 { 00063 n_dimensions2 = get_volume_n_dimensions( volume2 ); 00064 00065 ALLOC2D( *real_origins2, n_slices2, n_dimensions2 ); 00066 00067 for_less( s, 0, n_slices2 ) 00068 { 00069 get_mapping( volume2, origins2[s], x_axis2, y_axis2, 00070 x_translation2, y_translation2, x_scale2, y_scale2, 00071 (*real_origins2)[s], real_x_axis2, real_y_axis2 ); 00072 } 00073 } 00074 } 00075 00076 private void set_pixel_range( 00077 Volume volume1, 00078 int n_slices1, 00079 Real **real_origins1, 00080 Real real_x_axis1[], 00081 Real real_y_axis1[], 00082 Volume volume2, 00083 int n_slices2, 00084 Real **real_origins2, 00085 Real real_x_axis2[], 00086 Real real_y_axis2[], 00087 int x_viewport_size, 00088 int y_viewport_size, 00089 Pixel_types pixel_type, 00090 int *n_pixels_alloced, 00091 pixels_struct *pixels ) 00092 { 00093 int s, x_size, y_size; 00094 int x_pixel_start, x_pixel_end, y_pixel_start, y_pixel_end; 00095 00096 x_pixel_start = 0; 00097 x_pixel_end = x_viewport_size - 1; 00098 y_pixel_start = 0; 00099 y_pixel_end = y_viewport_size - 1; 00100 00101 for_less( s, 0, n_slices1 ) 00102 { 00103 clip_viewport_to_volume( volume1, real_origins1[s], 00104 real_x_axis1, real_y_axis1, 00105 &x_pixel_start, &x_pixel_end, 00106 &y_pixel_start, &y_pixel_end ); 00107 } 00108 00109 if( volume2 != NULL ) 00110 { 00111 for_less( s, 0, n_slices2 ) 00112 { 00113 clip_viewport_to_volume( volume2, real_origins2[s], 00114 real_x_axis2, real_y_axis2, 00115 &x_pixel_start, &x_pixel_end, 00116 &y_pixel_start, &y_pixel_end ); 00117 } 00118 } 00119 00120 x_size = x_pixel_end - x_pixel_start + 1; 00121 y_size = y_pixel_end - y_pixel_start + 1; 00122 if( x_size < 0 ) 00123 x_size = 0; 00124 if( y_size < 0 ) 00125 y_size = 0; 00126 00127 modify_pixels_size( n_pixels_alloced, pixels, x_size, y_size, pixel_type ); 00128 pixels->x_position = x_pixel_start; 00129 pixels->y_position = y_pixel_start; 00130 pixels->x_zoom = 1.0; 00131 pixels->y_zoom = 1.0; 00132 } 00133 00134 /* ----------------------------- MNI Header ----------------------------------- 00135 @NAME : create_weighted_volume_slices 00136 @INPUT : volume1, 00137 n_slices1 - number of parallel slices to weighted-sum. 00138 origins1 00139 x_axis1 00140 y_axis1 00141 weights1 00142 x_translation1 00143 y_translation1 00144 x_scale1 00145 y_scale1 00146 volume2, 00147 n_slices2 - number of parallel slices to weighted-sum. 00148 origins2 00149 x_axis2 00150 y_axis2 00151 weights2 00152 x_translation2 00153 y_translation2 00154 x_scale2 00155 y_scale2 00156 x_viewport_size 00157 y_viewport_size 00158 pixel_type 00159 degrees_continuity 00160 cmode_colour_map 00161 rgb_colour_map 00162 empty_colour 00163 n_pixels_alloced 00164 @OUTPUT : pixels 00165 @RETURNS : 00166 @DESCRIPTION: Renders a volume or a merge of 2 volumes to pixels. 00167 @METHOD : 00168 @GLOBALS : 00169 @CALLS : 00170 @CREATED : 1993 David MacDonald 00171 @MODIFIED : 00172 ---------------------------------------------------------------------------- */ 00173 00174 private void create_weighted_volume_slices( 00175 Volume volume1, 00176 int n_slices1, 00177 Real **origins1, 00178 Real x_axis1[], 00179 Real y_axis1[], 00180 Real weights1[], 00181 Volume volume2, 00182 int n_slices2, 00183 Real **origins2, 00184 Real x_axis2[], 00185 Real y_axis2[], 00186 Real weights2[], 00187 int x_pixel_start, 00188 int x_pixel_end, 00189 int y_pixel_start, 00190 int y_pixel_end, 00191 int degrees_continuity, 00192 unsigned short **cmode_colour_map, 00193 Colour **rgb_colour_map, 00194 Colour empty_colour, 00195 void *render_storage, 00196 pixels_struct *pixels ) 00197 { 00198 Data_types volume1_data_type, volume2_data_type; 00199 void *volume_data1, *volume_data2; 00200 int strides1[MAX_DIMENSIONS], strides2[MAX_DIMENSIONS]; 00201 int n_dimensions1, n_dimensions2; 00202 int axis, s, c; 00203 int sizes1[MAX_DIMENSIONS], sizes2[MAX_DIMENSIONS]; 00204 00205 if( pixels->x_size > 0 && pixels->y_size > 0 ) 00206 { 00207 n_dimensions1 = get_volume_n_dimensions( volume1 ); 00208 00209 for_less( s, 0, n_slices1 ) 00210 { 00211 for_less( c, 0, n_dimensions1 ) 00212 { 00213 origins1[s][c] += (Real) pixels->x_position * x_axis1[c] + 00214 (Real) pixels->y_position * y_axis1[c]; 00215 } 00216 } 00217 00218 if( !volume1->is_cached_volume ) { 00219 GET_VOXEL_PTR( volume_data1, volume1, 0, 0, 0, 0, 0 ); 00220 } 00221 volume1_data_type = get_volume_data_type( volume1 ); 00222 00223 get_volume_sizes( volume1, sizes1 ); 00224 00225 strides1[n_dimensions1-1] = 1; 00226 for_down( axis, n_dimensions1 - 2, 0 ) 00227 strides1[axis] = strides1[axis+1] * sizes1[axis+1]; 00228 00229 if( volume2 != NULL ) 00230 { 00231 n_dimensions2 = get_volume_n_dimensions( volume2 ); 00232 00233 for_less( s, 0, n_slices2 ) 00234 { 00235 for_less( c, 0, n_dimensions2 ) 00236 { 00237 origins2[s][c] += (Real) pixels->x_position * x_axis2[c] + 00238 (Real) pixels->y_position * y_axis2[c]; 00239 } 00240 } 00241 00242 if( !volume2->is_cached_volume ) { 00243 GET_VOXEL_PTR( volume_data2, volume2, 0, 0, 0, 0, 0 ); 00244 } 00245 volume2_data_type = get_volume_data_type( volume2 ); 00246 get_volume_sizes( volume2, sizes2 ); 00247 00248 strides2[n_dimensions2-1] = 1; 00249 for_down( axis, n_dimensions2 - 2, 0 ) 00250 strides2[axis] = strides2[axis+1] * sizes2[axis+1]; 00251 } 00252 else 00253 volume_data2 = NULL; 00254 00255 if( x_pixel_start > x_pixel_end || y_pixel_start > y_pixel_end ) 00256 { 00257 x_pixel_start = 0; 00258 x_pixel_end = pixels->x_size-1; 00259 y_pixel_start = 0; 00260 y_pixel_end = pixels->y_size-1; 00261 } 00262 else 00263 { 00264 if( x_pixel_start < 0 ) 00265 x_pixel_start = 0; 00266 if( x_pixel_end >= pixels->x_size ) 00267 x_pixel_end = pixels->x_size-1; 00268 00269 if( y_pixel_start < 0 ) 00270 y_pixel_start = 0; 00271 if( y_pixel_end >= pixels->y_size ) 00272 y_pixel_end = pixels->y_size-1; 00273 } 00274 00275 if( volume1->is_cached_volume || 00276 (volume2 != NULL && volume2->is_cached_volume) || 00277 (degrees_continuity != -1 && n_slices1 == 1 && 00278 (volume_data2 == NULL || n_slices2 == 1)) ) 00279 { 00280 interpolate_volume_to_slice( volume1, n_dimensions1, 00281 origins1[0], 00282 x_axis1, y_axis1, 00283 volume2, n_dimensions2, 00284 (volume_data2 == NULL) ? NULL : 00285 origins2[0], 00286 x_axis2, y_axis2, 00287 x_pixel_start, x_pixel_end, 00288 y_pixel_start, y_pixel_end, 00289 degrees_continuity, 00290 cmode_colour_map, rgb_colour_map, 00291 empty_colour, pixels ); 00292 } 00293 else 00294 { 00295 render_volume_to_slice( n_dimensions1, sizes1, volume_data1, 00296 volume1_data_type, n_slices1, weights1, 00297 strides1, origins1, x_axis1, y_axis1, 00298 n_dimensions2, sizes2, volume_data2, 00299 volume2_data_type, n_slices2, weights2, 00300 strides2, origins2, x_axis2, y_axis2, 00301 x_pixel_start, x_pixel_end, 00302 y_pixel_start, y_pixel_end, 00303 cmode_colour_map, 00304 rgb_colour_map, empty_colour, 00305 render_storage, pixels ); 00306 } 00307 } 00308 } 00309 00310 /* ----------------------------- MNI Header ----------------------------------- 00311 @NAME : get_filter_slices 00312 @INPUT : volume 00313 position 00314 x_axis 00315 y_axis 00316 filter_type 00317 filter_width 00318 @OUTPUT : n_slices 00319 origins 00320 weights 00321 @RETURNS : TRUE if successful 00322 @DESCRIPTION: Converts the filter type and width to a list of parallel slices 00323 to be used in a weighted-sum rendering. 00324 @METHOD : 00325 @GLOBALS : 00326 @CALLS : 00327 @CREATED : 1993 David MacDonald 00328 @MODIFIED : 00329 ---------------------------------------------------------------------------- */ 00330 00331 private BOOLEAN get_filter_slices( 00332 Volume volume, 00333 Real position[], 00334 Real x_axis[], 00335 Real y_axis[], 00336 Filter_types filter_type, 00337 Real filter_width, 00338 int *n_slices, 00339 Real ***origins, 00340 Real **weights ) 00341 { 00342 int c, a1, a2; 00343 Real direction[N_DIMENSIONS], separations[MAX_DIMENSIONS]; 00344 00345 if( filter_type != NEAREST_NEIGHBOUR ) 00346 { 00347 if( get_volume_n_dimensions( volume ) != N_DIMENSIONS ) 00348 { 00349 print_error( 00350 "If volume not 3D, can only do nearest neighbour filtering.\n" ); 00351 return( FALSE ); 00352 } 00353 00354 get_volume_separations( volume, separations ); 00355 for_less( c, 0, N_DIMENSIONS ) 00356 { 00357 a1 = (c + 1) % N_DIMENSIONS; 00358 a2 = (c + 2) % N_DIMENSIONS; 00359 direction[c] = x_axis[a1] * y_axis[a2] - x_axis[a2] * y_axis[a1]; 00360 direction[c] *= FABS( separations[a1] * separations[a2] / 00361 separations[c] ); 00362 } 00363 } 00364 00365 *n_slices = get_slice_weights_for_filter( volume, position, direction, 00366 filter_type, filter_width, 00367 origins, weights ); 00368 00369 return( TRUE ); 00370 } 00371 00372 /* ----------------------------- MNI Header ----------------------------------- 00373 @NAME : create_volume_slice 00374 @INPUT : volume1 - the volume to create a slice for 00375 filter_type1 - filter_type, usually NEAREST_NEIGHBOUR 00376 filter_width1 - width of filter, used for BOX, TRIANGLE, GAUSS 00377 slice_position1 - the voxel coordinate of the slice 00378 x_axis1 - the x axis in voxels 00379 y_axis1 - the y axis in voxels 00380 x_translation1 - pixel translation for viewing 00381 y_translation1 - pixel translation for viewing 00382 x_scale1 - pixel zoom for viewing 00383 y_scale1 - pixel zoom for viewing 00384 volume2 - second volume to be merged with first, or null 00385 filter_type2 - filter_type, usually NEAREST_NEIGHBOUR 00386 filter_width2 - width of filter, used for BOX, TRIANGLE, GAUSS 00387 slice_position2 - the voxel coordinate of the slice 00388 x_axis2 - the x axis in voxels 00389 y_axis2 - the y axis in voxels 00390 x_translation2 - pixel translation for viewing 00391 y_translation2 - pixel translation for viewing 00392 x_scale2 - pixel zoom for viewing 00393 y_scale2 - pixel zoom for viewing 00394 x_axis_index - X,Y, or Z 00395 y_axis_index - X,Y, or Z 00396 axis_index - X,Y, or Z 00397 x_viewport_size - will be clipped to this size 00398 y_viewport_size - will be clipped to this size 00399 pixel_type - RGB_PIXEL or COLOUR_INDEX_PIXEL for rgb/cmap 00400 interpolation_flag - ignored for now 00401 cmode_colour_map - if pixel_type == COLOUR_INDEX_PIXEL, then 00402 2d array of 16 bit colour indices for merged slices, 00403 or pointer to 1d array of colour indices for volume1 00404 rgb_colour_map - if pixel_type == RGB_PIXEL, then 00405 2d array of 24 bit colours for merged slices, 00406 or pointer to 1d array of colours for volume1 00407 @OUTPUT : n_pixels_alloced - a pointer to the size alloced. Before first 00408 call, set size alloced to zero, and all calls, 00409 pass pointer to size alloced, and pointer to pixels. 00410 pixels - 2d pixels array created, and realloced as 00411 necessary, assuming, n_pixels_alloced is a 00412 pointer to the current alloc size of pixels. 00413 @RETURNS : 00414 @DESCRIPTION: Creates a slice of one volume or merged slice of two, suitable 00415 for graphics display. 00416 @CREATED : Mar 1993 David MacDonald 00417 @MODIFIED : 00418 ---------------------------------------------------------------------------- */ 00419 00420 public void create_volume_slice( 00421 Volume volume1, 00422 Filter_types filter_type1, 00423 Real filter_width1, 00424 Real slice_position1[], 00425 Real x_axis1[], 00426 Real y_axis1[], 00427 Real x_translation1, 00428 Real y_translation1, 00429 Real x_scale1, 00430 Real y_scale1, 00431 Volume volume2, 00432 Filter_types filter_type2, 00433 Real filter_width2, 00434 Real slice_position2[], 00435 Real x_axis2[], 00436 Real y_axis2[], 00437 Real x_translation2, 00438 Real y_translation2, 00439 Real x_scale2, 00440 Real y_scale2, 00441 int x_viewport_size, 00442 int y_viewport_size, 00443 int x_pixel_start, 00444 int x_pixel_end, 00445 int y_pixel_start, 00446 int y_pixel_end, 00447 Pixel_types pixel_type, 00448 int degrees_continuity, 00449 unsigned short **cmode_colour_map, 00450 Colour **rgb_colour_map, 00451 Colour empty_colour, 00452 void *render_storage, 00453 BOOLEAN clip_pixels_flag, 00454 int *n_pixels_alloced, 00455 pixels_struct *pixels ) 00456 { 00457 int n_slices1, n_slices2; 00458 Real **positions1, **positions2, *weights1, *weights2; 00459 Real **real_origins1, **real_origins2; 00460 Real real_x_axis1[MAX_DIMENSIONS], real_y_axis1[MAX_DIMENSIONS]; 00461 Real real_x_axis2[MAX_DIMENSIONS], real_y_axis2[MAX_DIMENSIONS]; 00462 00463 if( !get_filter_slices( volume1, slice_position1, x_axis1, y_axis1, 00464 filter_type1, filter_width1, &n_slices1, 00465 &positions1, &weights1 ) ) 00466 { 00467 modify_pixels_size( n_pixels_alloced, pixels, 0, 0, pixel_type ); 00468 return; 00469 } 00470 00471 if( volume2 != NULL ) 00472 { 00473 if( !get_filter_slices( volume2, slice_position2, x_axis2, y_axis2, 00474 filter_type2, filter_width2, &n_slices2, 00475 &positions2, &weights2 ) ) 00476 { 00477 modify_pixels_size( n_pixels_alloced, pixels, 0, 0, pixel_type ); 00478 return; 00479 } 00480 } 00481 00482 create_pixel_mapping( volume1, n_slices1, positions1, x_axis1, y_axis1, 00483 x_translation1, y_translation1, x_scale1, y_scale1, 00484 volume2, n_slices2, positions2, x_axis2, y_axis2, 00485 x_translation2, y_translation2, x_scale2, y_scale2, 00486 real_x_axis1, real_y_axis1, &real_origins1, 00487 real_x_axis2, real_y_axis2, &real_origins2 ); 00488 00489 if( clip_pixels_flag ) 00490 { 00491 set_pixel_range( volume1, n_slices1, 00492 real_origins1, real_x_axis1, real_y_axis1, 00493 volume2, n_slices2, 00494 real_origins2, real_x_axis2, real_y_axis2, 00495 x_viewport_size, y_viewport_size, 00496 pixel_type, n_pixels_alloced, pixels ); 00497 } 00498 00499 create_weighted_volume_slices( volume1, n_slices1, 00500 real_origins1, real_x_axis1, real_y_axis1, 00501 weights1, 00502 volume2, n_slices2, 00503 real_origins2, real_x_axis2, real_y_axis2, 00504 weights2, 00505 x_pixel_start, x_pixel_end, 00506 y_pixel_start, y_pixel_end, 00507 degrees_continuity, 00508 cmode_colour_map, rgb_colour_map, 00509 empty_colour, render_storage, 00510 pixels ); 00511 00512 if( volume2 != NULL ) 00513 { 00514 FREE2D( positions2 ); 00515 FREE( weights2 ); 00516 FREE2D( real_origins2 ); 00517 } 00518 00519 FREE2D( positions1 ); 00520 FREE( weights1 ); 00521 FREE2D( real_origins1 ); 00522 } 00523 00524 public void set_volume_slice_pixel_range( 00525 Volume volume1, 00526 Filter_types filter_type1, 00527 Real filter_width1, 00528 Real slice_position1[], 00529 Real x_axis1[], 00530 Real y_axis1[], 00531 Real x_translation1, 00532 Real y_translation1, 00533 Real x_scale1, 00534 Real y_scale1, 00535 Volume volume2, 00536 Filter_types filter_type2, 00537 Real filter_width2, 00538 Real slice_position2[], 00539 Real x_axis2[], 00540 Real y_axis2[], 00541 Real x_translation2, 00542 Real y_translation2, 00543 Real x_scale2, 00544 Real y_scale2, 00545 int x_viewport_size, 00546 int y_viewport_size, 00547 Pixel_types pixel_type, 00548 int *n_pixels_alloced, 00549 pixels_struct *pixels ) 00550 { 00551 int n_slices1, n_slices2; 00552 Real **positions1, **positions2, *weights1, *weights2; 00553 Real **real_origins1, **real_origins2; 00554 Real real_x_axis1[MAX_DIMENSIONS], real_y_axis1[MAX_DIMENSIONS]; 00555 Real real_x_axis2[MAX_DIMENSIONS], real_y_axis2[MAX_DIMENSIONS]; 00556 00557 if( !get_filter_slices( volume1, slice_position1, x_axis1, y_axis1, 00558 filter_type1, filter_width1, &n_slices1, 00559 &positions1, &weights1 ) ) 00560 { 00561 modify_pixels_size( n_pixels_alloced, pixels, 0, 0, pixel_type ); 00562 return; 00563 } 00564 00565 if( volume2 != NULL ) 00566 { 00567 if( !get_filter_slices( volume2, slice_position2, x_axis2, y_axis2, 00568 filter_type2, filter_width2, &n_slices2, 00569 &positions2, &weights2 ) ) 00570 { 00571 modify_pixels_size( n_pixels_alloced, pixels, 0, 0, pixel_type ); 00572 return; 00573 } 00574 } 00575 00576 create_pixel_mapping( volume1, n_slices1, positions1, x_axis1, y_axis1, 00577 x_translation1, y_translation1, x_scale1, y_scale1, 00578 volume2, n_slices2, positions2, x_axis2, y_axis2, 00579 x_translation2, y_translation2, x_scale2, y_scale2, 00580 real_x_axis1, real_y_axis1, &real_origins1, 00581 real_x_axis2, real_y_axis2, &real_origins2 ); 00582 00583 set_pixel_range( volume1, n_slices1, 00584 real_origins1, real_x_axis1, real_y_axis1, 00585 volume2, n_slices2, 00586 real_origins2, real_x_axis2, real_y_axis2, 00587 x_viewport_size, y_viewport_size, 00588 pixel_type, n_pixels_alloced, pixels ); 00589 00590 if( volume2 != NULL ) 00591 { 00592 FREE2D( positions2 ); 00593 FREE( weights2 ); 00594 FREE2D( real_origins2 ); 00595 } 00596 00597 FREE2D( positions1 ); 00598 FREE( weights1 ); 00599 FREE2D( real_origins1 ); 00600 }