00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013 #include <stdio.h>
00014 #include <ctype.h>
00015 #include <string.h>
00016 #include <sys/stat.h>
00017
00018
00019
00020 #ifndef SEEK_END
00021 #define SEEK_END 2
00022 #endif
00023
00024 #ifndef SEEK_SET
00025 #define SEEK_SET 0
00026 #endif
00027
00028 #include "mrilib.h"
00029 #include "ge4_header.h"
00030
00031
00032 short check_dicom_magic_num( char * );
00033
00034
00035
00036 char MRILIB_orients[8] = "\0\0\0\0\0\0\0\0" ;
00037
00038
00039
00040 float MRILIB_zoff = 0.0 ;
00041
00042
00043
00044 float MRILIB_tr = 0.0 ;
00045
00046
00047
00048 float MRILIB_xoff = 0.0 ;
00049
00050
00051
00052 float MRILIB_yoff = 0.0 ;
00053
00054
00055
00056 int use_MRILIB_xoff = 0 ;
00057
00058
00059
00060 int use_MRILIB_yoff = 0 ;
00061
00062
00063
00064 int use_MRILIB_zoff = 0 ;
00065
00066
00067
00068 int use_MRILIB_xcos = 0 ;
00069
00070
00071
00072 float MRILIB_xcos[3] = { 1.0 , 0.0 , 0.0 } ;
00073
00074
00075
00076 int use_MRILIB_ycos = 0 ;
00077
00078
00079
00080 float MRILIB_ycos[3] = { 0.0 , 1.0 , 0.0 } ;
00081
00082
00083
00084 int use_MRILIB_zcos = 0 ;
00085
00086
00087
00088 float MRILIB_zcos[3] = { 0.0 , 0.0 , 1.0 } ;
00089
00090
00091
00092 int use_MRILIB_slicespacing = 0 ;
00093
00094
00095
00096 float MRILIB_slicespacing = 0.0 ;
00097
00098
00099
00100 MRI_IMAGE *mri_try_mri( FILE * , int * ) ;
00101 MRI_IMAGE *mri_try_7D ( FILE * , int * ) ;
00102 MRI_IMAGE *mri_try_pgm( FILE * , int * ) ;
00103
00104
00105
00106
00107
00108
00109
00110
00111
00112
00113
00114 typedef struct {
00115 int size ;
00116 int head ;
00117 char * prefix ;
00118 } MCW_imsize ;
00119
00120
00121
00122 #define MAX_MCW_IMSIZE 99
00123
00124
00125
00126 static MCW_imsize imsize[MAX_MCW_IMSIZE] ;
00127
00128
00129
00130 static int MCW_imsize_good = -1 ;
00131
00132
00133
00134 #undef swap_4
00135 #undef swap_2
00136
00137
00138
00139 static void swap_4(void *ppp)
00140 {
00141 unsigned char *pntr = (unsigned char *) ppp ;
00142 unsigned char b0, b1, b2, b3;
00143
00144 b0 = *pntr; b1 = *(pntr+1); b2 = *(pntr+2); b3 = *(pntr+3);
00145 *pntr = b3; *(pntr+1) = b2; *(pntr+2) = b1; *(pntr+3) = b0;
00146 }
00147
00148
00149
00150
00151
00152 static void swap_8(void *ppp)
00153 {
00154 unsigned char *pntr = (unsigned char *) ppp ;
00155 unsigned char b0, b1, b2, b3;
00156 unsigned char b4, b5, b6, b7;
00157
00158 b0 = *pntr ; b1 = *(pntr+1); b2 = *(pntr+2); b3 = *(pntr+3);
00159 b4 = *(pntr+4); b5 = *(pntr+5); b6 = *(pntr+6); b7 = *(pntr+7);
00160
00161 *pntr = b7; *(pntr+1) = b6; *(pntr+2) = b5; *(pntr+3) = b4;
00162 *(pntr+4) = b3; *(pntr+5) = b2; *(pntr+6) = b1; *(pntr+7) = b0;
00163 }
00164
00165
00166
00167
00168
00169 static void swap_2(void *ppp)
00170 {
00171 unsigned char *pntr = (unsigned char *) ppp ;
00172 unsigned char b0, b1;
00173
00174 b0 = *pntr; b1 = *(pntr+1);
00175 *pntr = b1; *(pntr+1) = b0;
00176 }
00177
00178
00179
00180
00181
00182
00183
00184
00185
00186
00187
00188 MRI_IMAGE *mri_read( char *fname )
00189 {
00190 FILE *imfile ;
00191 MRI_IMAGE *im ;
00192 int length , skip=0 , swap=0 ;
00193 void *data ;
00194
00195 ENTRY("mri_read") ;
00196
00197 if( fname == NULL || *fname == '\0' ) RETURN(NULL) ;
00198
00199
00200
00201 if( strstr(fname,".jpg" ) != NULL ||
00202 strstr(fname,".JPG" ) != NULL ||
00203 strstr(fname,".jpeg") != NULL ||
00204 strstr(fname,".JPEG") != NULL ||
00205 strstr(fname,".gif" ) != NULL ||
00206 strstr(fname,".GIF" ) != NULL ||
00207 strstr(fname,".tif" ) != NULL ||
00208 strstr(fname,".TIF" ) != NULL ||
00209 strstr(fname,".tiff") != NULL ||
00210 strstr(fname,".TIFF") != NULL ||
00211 strstr(fname,".bmp" ) != NULL ||
00212 strstr(fname,".BMP" ) != NULL ||
00213 strstr(fname,".pbm" ) != NULL ||
00214 strstr(fname,".PBM" ) != NULL ||
00215 strstr(fname,".pgm" ) != NULL ||
00216 strstr(fname,".PGM" ) != NULL ||
00217 strstr(fname,".ppm" ) != NULL ||
00218 strstr(fname,".PPM" ) != NULL ||
00219 strstr(fname,".png" ) != NULL ||
00220 strstr(fname,".PNG" ) != NULL ){
00221
00222 im = mri_read_stuff(fname) ; if( im != NULL ) RETURN(im) ;
00223 }
00224
00225
00226
00227 imfile = fopen( fname , "r" ) ;
00228 if( imfile == NULL ){
00229 fprintf( stderr , "couldn't open image file %s\n" , fname ) ;
00230 RETURN( NULL );
00231 }
00232
00233 fseek( imfile , 0L , SEEK_END ) ;
00234 length = ftell( imfile ) ;
00235
00236
00237
00238
00239 { char str[5]="AFNI" ;
00240 int nx , ny , bpp , cflag , hdroff , extraskip=0 ;
00241 rewind(imfile) ; fread(str,1,4,imfile) ;
00242
00243 if( str[0]=='G' && str[1]=='E' && str[2]=='M' && str[3]=='S' ){
00244 char buf[4096]; int bb,cc;
00245 rewind(imfile); cc=fread(buf,1,4096,imfile); cc-=4 ;
00246 for( bb=4; bb < cc ; bb++ )
00247 if( buf[bb]=='I' && buf[bb+1]=='M' && buf[bb+2]=='G' && buf[bb+3]=='F' ) break ;
00248 if( bb < cc ){
00249 fseek( imfile , (long)bb , SEEK_SET ) ; extraskip = bb ;
00250 fread(str,1,4,imfile) ;
00251 }
00252 }
00253
00254
00255
00256 if( str[0]=='I' && str[1]=='M' && str[2]=='G' && str[3]=='F' ){
00257
00258 fread( &skip , 4,1, imfile ) ;
00259 fread( &nx , 4,1, imfile ) ;
00260 fread( &ny , 4,1, imfile ) ;
00261 fread( &bpp , 4,1, imfile ) ;
00262 fread( &cflag, 4,1, imfile ) ;
00263
00264 if( nx < 0 || nx > 8192 ){
00265 swap = 1 ;
00266 swap_4(&skip); swap_4(&nx) ;
00267 swap_4(&ny) ; swap_4(&bpp); swap_4(&cflag);
00268 }
00269 skip += extraskip ;
00270 if( nx < 0 || nx > 8192 || ny < 0 || ny > 8192 ) goto The_Old_Way ;
00271 if( skip < 0 || skip >= length ) goto The_Old_Way ;
00272 if( bpp != 16 || cflag != 1 ) goto The_Old_Way ;
00273
00274
00275
00276 im = mri_new( nx , ny , MRI_short ) ;
00277
00278
00279
00280 length = fseek( imfile , 148L+extraskip , SEEK_SET ) ;
00281 if( length == 0 ){
00282 fread( &hdroff , 4,1 , imfile ) ;
00283 if( swap ) swap_4(&hdroff) ;
00284 if( hdroff > 0 ){
00285 float dx,dy,dz, dxx,dyy,dzz, xyz[9], zz ; int itr, ii,jj,kk, qq ;
00286 static int nzoff=0 ;
00287 static float zoff ;
00288
00289
00290
00291 fseek( imfile , hdroff+26+extraskip , SEEK_SET ) ;
00292 fread( &dzz , 4,1 , imfile ) ;
00293
00294 fseek( imfile , hdroff+50+extraskip , SEEK_SET ) ;
00295 fread( &dxx , 4,1 , imfile ) ;
00296 fread( &dyy , 4,1 , imfile ) ;
00297
00298 if( swap ){ swap_4(&dxx); swap_4(&dyy); swap_4(&dzz); }
00299
00300
00301
00302 if( dxx > 0.01 && dyy > 0.01 && dzz > 0.01 ){
00303 im->dx = dxx; im->dy = dyy; im->dz = dzz; im->dw = 1.0;
00304 }
00305
00306
00307
00308
00309
00310
00311
00312
00313
00314 fseek( imfile , hdroff+154+extraskip , SEEK_SET ) ;
00315 fread( xyz , 4,9 , imfile ) ;
00316 if( swap ) swap_fourbytes(9,xyz) ;
00317
00318
00319
00320
00321
00322
00323 dx = fabs(xyz[3]-xyz[0]) ; ii = 1 ;
00324 dy = fabs(xyz[4]-xyz[1]) ; if( dy > dx ){ ii=2; dx=dy; }
00325 dz = fabs(xyz[5]-xyz[2]) ; if( dz > dx ){ ii=3; }
00326 dx = xyz[ii+2]-xyz[ii-1] ; if( dx < 0. ){ ii = -ii; }
00327 switch( ii ){
00328 case 1: MRILIB_orients[0] = 'L'; MRILIB_orients[1] = 'R'; break;
00329 case -1: MRILIB_orients[0] = 'R'; MRILIB_orients[1] = 'L'; break;
00330 case 2: MRILIB_orients[0] = 'P'; MRILIB_orients[1] = 'A'; break;
00331 case -2: MRILIB_orients[0] = 'A'; MRILIB_orients[1] = 'P'; break;
00332 case 3: MRILIB_orients[0] = 'I'; MRILIB_orients[1] = 'S'; break;
00333 case -3: MRILIB_orients[0] = 'S'; MRILIB_orients[1] = 'I'; break;
00334 default: MRILIB_orients[0] ='\0'; MRILIB_orients[1] ='\0'; break;
00335 }
00336
00337
00338
00339
00340
00341
00342 dx = fabs(xyz[6]-xyz[3]) ; jj = 1 ;
00343 dy = fabs(xyz[7]-xyz[4]) ; if( dy > dx ){ jj=2; dx=dy; }
00344 dz = fabs(xyz[8]-xyz[5]) ; if( dz > dx ){ jj=3; }
00345 dx = xyz[jj+5]-xyz[jj+2] ; if( dx < 0. ){ jj = -jj; }
00346 switch( jj ){
00347 case 1: MRILIB_orients[2] = 'L'; MRILIB_orients[3] = 'R'; break;
00348 case -1: MRILIB_orients[2] = 'R'; MRILIB_orients[3] = 'L'; break;
00349 case 2: MRILIB_orients[2] = 'P'; MRILIB_orients[3] = 'A'; break;
00350 case -2: MRILIB_orients[2] = 'A'; MRILIB_orients[3] = 'P'; break;
00351 case 3: MRILIB_orients[2] = 'I'; MRILIB_orients[3] = 'S'; break;
00352 case -3: MRILIB_orients[2] = 'S'; MRILIB_orients[3] = 'I'; break;
00353 default: MRILIB_orients[2] ='\0'; MRILIB_orients[3] ='\0'; break;
00354 }
00355
00356 MRILIB_orients[6] = '\0' ;
00357
00358 kk = 6 - abs(ii)-abs(jj) ;
00359
00360 zz = xyz[kk-1] ;
00361
00362 im->zo = zz ;
00363
00364
00365
00366 if( nzoff == 0 ){
00367
00368 zoff = zz ;
00369 switch( kk ){
00370 case 1: MRILIB_orients[4] = 'L'; MRILIB_orients[5] = 'R'; break;
00371 case 2: MRILIB_orients[4] = 'P'; MRILIB_orients[5] = 'A'; break;
00372 case 3: MRILIB_orients[4] = 'I'; MRILIB_orients[5] = 'S'; break;
00373 default: MRILIB_orients[4] ='\0'; MRILIB_orients[5] ='\0'; break;
00374 }
00375
00376 } else if( nzoff == 1 ){
00377
00378 float qoff = zz - zoff ;
00379 if( qoff < 0 ) kk = -kk ;
00380
00381 if( !use_MRILIB_slicespacing && qoff != 0.0 ){
00382 use_MRILIB_slicespacing = 1 ;
00383 MRILIB_slicespacing = fabs(qoff) ;
00384 }
00385
00386 switch( kk ){
00387 case 1: MRILIB_orients[4] = 'L'; MRILIB_orients[5] = 'R'; break;
00388 case -1: MRILIB_orients[4] = 'R'; MRILIB_orients[5] = 'L'; break;
00389 case 2: MRILIB_orients[4] = 'P'; MRILIB_orients[5] = 'A'; break;
00390 case -2: MRILIB_orients[4] = 'A'; MRILIB_orients[5] = 'P'; break;
00391 case 3: MRILIB_orients[4] = 'I'; MRILIB_orients[5] = 'S'; break;
00392 case -3: MRILIB_orients[4] = 'S'; MRILIB_orients[5] = 'I'; break;
00393 default: MRILIB_orients[4] ='\0'; MRILIB_orients[5] ='\0'; break;
00394 }
00395
00396
00397
00398
00399
00400 MRILIB_zoff = zoff ; use_MRILIB_zoff = 1 ;
00401 if( kk == 1 || kk == 2 || kk == 3 ) MRILIB_zoff = -MRILIB_zoff ;
00402
00403
00404
00405
00406
00407
00408
00409
00410 qq = abs(ii) ;
00411 MRILIB_xoff = ( xyz[qq-1]*(nx-0.5) + xyz[qq+2]*0.5 ) / nx ;
00412 if( ii == 1 || ii == 2 || ii == 3 ) MRILIB_xoff = -MRILIB_xoff ;
00413 use_MRILIB_xoff = ( fabs(xyz[qq+2]-xyz[qq+5]) < 0.01*dxx ) ;
00414
00415
00416
00417
00418
00419
00420
00421 qq = abs(jj) ;
00422 MRILIB_yoff = ( xyz[qq-1]*(ny-0.5) + xyz[qq+5]*0.5 ) / ny ;
00423 if( jj == 1 || jj == 2 || jj == 3 ) MRILIB_yoff = -MRILIB_yoff ;
00424 use_MRILIB_yoff = ( fabs(xyz[qq-1]-xyz[qq+2]) < 0.01*dyy ) ;
00425 }
00426 nzoff++ ;
00427
00428
00429
00430 if( MRILIB_tr <= 0.0 ){
00431 fseek( imfile , hdroff+194+extraskip , SEEK_SET ) ;
00432 fread( &itr , 4,1 , imfile ) ;
00433 if( swap ) swap_4(&itr) ;
00434 MRILIB_tr = im->dt = 1.0e-6 * itr ;
00435 }
00436 }
00437 }
00438
00439 goto Ready_To_Roll ;
00440 }
00441 }
00442
00443
00444
00445 The_Old_Way:
00446
00447 #if 0
00448 MRILIB_orients[0] = '\0' ; MRILIB_zoff = MRILIB_tr = 0.0 ;
00449 #endif
00450
00451 switch( length ){
00452
00453 case 512:
00454 im = mri_new( 16 , 16 , MRI_short ) ;
00455 break ;
00456
00457 case 2048:
00458 im = mri_new( 32 , 32 , MRI_short ) ;
00459 break ;
00460
00461 case 4096:
00462 im = mri_new( 64 , 64 , MRI_byte ) ;
00463 break ;
00464
00465 case 8192:
00466 case 16096:
00467 im = mri_new( 64 , 64 , MRI_short ) ;
00468 skip = length - 8192 ;
00469 break ;
00470
00471 #if 0
00472 case 18432:
00473 im = mri_new( 96 , 96 , MRI_short ) ;
00474 break ;
00475 #endif
00476
00477 case 16384:
00478 im = mri_new( 128 , 128 , MRI_byte ) ;
00479 break ;
00480
00481 case 32768:
00482 case 40672:
00483 im = mri_new( 128 , 128 , MRI_short ) ;
00484 skip = length - 32768 ;
00485 break ;
00486
00487 case 65536:
00488 im = mri_new( 256 , 256 , MRI_byte ) ;
00489 break ;
00490
00491 case 131072:
00492 case 138976:
00493 case 145408:
00494
00495 im = mri_new( 256 , 256 , MRI_short ) ;
00496 skip = length - 131072 ;
00497 break ;
00498
00499 #if 0
00500 case 262144:
00501 im = mri_new( 256 , 256 , MRI_float ) ;
00502 break ;
00503 #else
00504 case 262144:
00505 im = mri_new( 512 , 512 , MRI_byte ) ;
00506 break ;
00507
00508 case 524288:
00509 im = mri_new( 512 , 512 , MRI_short ) ;
00510 break ;
00511
00512 case 1048576:
00513 im = mri_new( 1024 , 1024 , MRI_byte ) ;
00514 break ;
00515
00516 case 2097152:
00517 im = mri_new( 1024 , 1024 , MRI_short ) ;
00518 break ;
00519 #endif
00520
00521
00522
00523 default:
00524 im = mri_try_mri( imfile , &skip ) ;
00525 if( im == NULL ) im = mri_try_7D ( imfile , &skip ) ;
00526 if( im == NULL ) im = mri_try_pgm( imfile , &skip ) ;
00527 if( im != NULL ) break ;
00528
00529 fclose( imfile ) ;
00530
00531 im = mri_read_ascii( fname ) ;
00532 if( im != NULL ) RETURN( im );
00533
00534 im = mri_read_ppm( fname ) ;
00535 if( im != NULL ) RETURN( im );
00536
00537 im = mri_read_stuff( fname ) ;
00538 if( im != NULL ) RETURN( im );
00539
00540 fprintf( stderr , "do not recognize image file %s\n" , fname );
00541 fprintf( stderr , "length seen as %d\n" , length ) ;
00542 RETURN( NULL );
00543 }
00544
00545
00546
00547 Ready_To_Roll:
00548
00549 data = mri_data_pointer( im ) ;
00550
00551 length = fseek( imfile , skip , SEEK_SET ) ;
00552 if( length != 0 ){
00553 fprintf( stderr , "mri_read error in skipping in file %s\n" , fname ) ;
00554 mri_free( im ) ;
00555 RETURN( NULL );
00556 }
00557
00558 length = fread( data , im->pixel_size , im->nvox , imfile ) ;
00559 fclose( imfile ) ;
00560
00561 if( length != im->nvox ){
00562 mri_free( im ) ;
00563 fprintf( stderr , "couldn't read image data from file %s\n" , fname ) ;
00564 RETURN( NULL );
00565 }
00566
00567 mri_add_name( fname , im ) ;
00568
00569
00570
00571 if( swap ){
00572 switch( im->pixel_size ){
00573 default: break ;
00574 case 2: swap_twobytes ( im->nvox, data ) ; break ;
00575 case 4: swap_fourbytes( im->nvox, data ) ; break ;
00576 case 8: swap_fourbytes( 2*im->nvox, data ) ; break ;
00577 }
00578
00579 im->was_swapped = 1 ;
00580 }
00581
00582 RETURN( im );
00583 }
00584
00585
00586
00587
00588
00589
00590
00591
00592
00593
00594
00595
00596
00597 MRI_IMAGE * mri_read_ge4( char * filename )
00598 {
00599 MRI_IMAGE * im;
00600 ge4_header H;
00601
00602 ENTRY( "mri_read_ge4" );
00603
00604 if ( filename == NULL )
00605 {
00606 fprintf( stderr, "** mri_read_ge4 - missing filename\n" );
00607 RETURN( NULL );
00608 }
00609
00610
00611 if ( ge4_read_header( &H, filename, True ) != 0 )
00612 RETURN( NULL );
00613
00614
00615 if ( (im = mri_new(256, 256, MRI_short)) == NULL )
00616 {
00617 free(H.image);
00618 RETURN( NULL );
00619 }
00620
00621
00622 im->zo = H.im_h.im_loc;
00623 im->dt = H.im_h.tr;
00624 im->was_swapped = H.swap;
00625
00626 if ( ( H.ser_h.fov > 1.0 ) &&
00627 ( H.ser_h.fov < 1000.0 ) &&
00628 ( H.ser_h.scan_mat_x > 0 ) &&
00629 ( H.ser_h.scan_mat_x < 1000 ) &&
00630 ( H.ser_h.scan_mat_y > 0 ) &&
00631 ( H.ser_h.scan_mat_y < 1000 ) )
00632 {
00633
00634
00635 im->dx = 2 * H.ser_h.fov / H.ser_h.scan_mat_x;
00636 im->dy = im->dx;
00637 im->dz = 2 * H.ser_h.fov / H.ser_h.scan_mat_y;
00638 im->dw = 1;
00639 }
00640
00641 memcpy( mri_data_pointer(im), H.image, H.im_bytes );
00642
00643 mri_add_name( filename, im );
00644
00645 free(H.image);
00646
00647 RETURN( im );
00648 }
00649
00650
00651
00652
00653 #if 0
00654 #define NUMSCAN(var) \
00655 { while( isspace(ch) ) {ch = getc(imfile) ;} \
00656 if(ch == '#') do{ch = getc(imfile) ;}while(ch != '\n' && ch != EOF) ; \
00657 if(ch =='\n') ch = getc(imfile) ; \
00658 for( nch=0 ; isdigit(ch) ; nch++,ch=getc(imfile) ) {buf[nch] = ch ;} \
00659 buf[nch]='\0'; \
00660 var = strtol( buf , NULL , 10 ) ; }
00661 #else
00662
00663
00664
00665 #define SKIPCOM \
00666 {if(ch == '#') do{ch = getc(imfile) ;}while(ch != '\n' && ch != EOF);}
00667
00668
00669
00670 #define NUMSCAN(var) \
00671 { SKIPCOM ; \
00672 while( ch!=EOF && !isdigit(ch) ){ch = getc(imfile); SKIPCOM; } \
00673 for( nch=0 ; isdigit(ch) ; nch++,ch=getc(imfile) ) {buf[nch] = ch ;} \
00674 buf[nch]='\0'; \
00675 var = strtol( buf , NULL , 10 ) ; }
00676 #endif
00677
00678
00679
00680
00681
00682
00683
00684
00685
00686
00687
00688 MRI_IMAGE *mri_try_mri( FILE *imfile , int *skip )
00689 {
00690 int ch , nch , nx,ny,imcode ;
00691 char buf[64] ;
00692 MRI_IMAGE *im ;
00693
00694 ENTRY("mri_try_mri") ;
00695
00696 fseek( imfile , 0 , SEEK_SET ) ;
00697
00698 ch = getc( imfile ) ; if( ch != 'M' ) RETURN( NULL );
00699 ch = getc( imfile ) ; if( ch != 'R' ) RETURN( NULL );
00700 ch = getc( imfile ) ; if( ch != 'I' ) RETURN( NULL );
00701
00702
00703
00704 ch = getc(imfile) ;
00705
00706 NUMSCAN(imcode) ;
00707 NUMSCAN(nx) ; if( nx <= 0 ) RETURN( NULL );
00708 NUMSCAN(ny) ; if( ny <= 0 ) RETURN( NULL );
00709
00710 *skip = ftell(imfile) ;
00711 im = mri_new( nx , ny , imcode ) ;
00712 RETURN( im );
00713 }
00714
00715
00716
00717
00718
00719
00720
00721 MRI_IMAGE *mri_try_7D( FILE *imfile , int *skip )
00722 {
00723 int ch , nch , nx,ny,nz,nt,nu,nv,nw , imcode , ndim ;
00724 char buf[64] ;
00725 MRI_IMAGE *im ;
00726
00727 ENTRY("mri_try_7D") ;
00728
00729 fseek( imfile , 0 , SEEK_SET ) ;
00730
00731 ch = getc( imfile ) ; if( ch != 'M' ) RETURN( NULL );
00732 ch = getc( imfile ) ; if( ch != 'R' ) RETURN( NULL );
00733 ch = getc( imfile ) ;
00734 switch( ch ){
00735 default: RETURN( NULL );
00736
00737 case '1': ndim = 1 ; break ;
00738 case '2': ndim = 2 ; break ;
00739 case '3': ndim = 3 ; break ;
00740 case '4': ndim = 4 ; break ;
00741 case '5': ndim = 5 ; break ;
00742 case '6': ndim = 6 ; break ;
00743 case '7': ndim = 7 ; break ;
00744 }
00745
00746
00747 ch = getc(imfile) ;
00748 NUMSCAN(imcode) ;
00749
00750 nx = ny = nz = nt = nu = nv = nw = 1 ;
00751
00752 NUMSCAN(nx) ; if( nx <= 0 ) RETURN( NULL );
00753 if( ndim > 1 ){ NUMSCAN(ny) ; if( ny <= 0 ) RETURN( NULL ); }
00754 if( ndim > 2 ){ NUMSCAN(nz) ; if( nz <= 0 ) RETURN( NULL ); }
00755 if( ndim > 3 ){ NUMSCAN(nt) ; if( nt <= 0 ) RETURN( NULL ); }
00756 if( ndim > 4 ){ NUMSCAN(nu) ; if( nu <= 0 ) RETURN( NULL ); }
00757 if( ndim > 5 ){ NUMSCAN(nv) ; if( nv <= 0 ) RETURN( NULL ); }
00758 if( ndim > 6 ){ NUMSCAN(nw) ; if( nw <= 0 ) RETURN( NULL ); }
00759
00760 *skip = ftell(imfile) ;
00761 im = mri_new_7D_generic( nx,ny,nz,nt,nu,nv,nw , imcode , TRUE ) ;
00762 RETURN( im );
00763 }
00764
00765
00766
00767
00768
00769
00770
00771
00772
00773
00774
00775
00776
00777 MRI_IMAGE *mri_try_pgm( FILE *imfile , int *skip )
00778 {
00779 int ch , nch , nx,ny,maxval ;
00780 char buf[64] ;
00781 MRI_IMAGE *im ;
00782
00783 ENTRY("mri_try_pgm") ;
00784
00785 fseek( imfile , 0 , SEEK_SET ) ;
00786
00787 ch = getc( imfile ) ; if( ch != 'P' ) RETURN(NULL);
00788 ch = getc( imfile ) ; if( ch != '5' ) RETURN(NULL);
00789
00790
00791
00792 ch = getc(imfile) ;
00793
00794 NUMSCAN(nx) ; if( nx <= 0 ) RETURN(NULL);
00795 NUMSCAN(ny) ; if( ny <= 0 ) RETURN(NULL);
00796 NUMSCAN(maxval) ; if( maxval <= 0 || maxval > 255 ) RETURN(NULL);
00797
00798 *skip = ftell(imfile) ;
00799 im = mri_new( nx , ny , MRI_byte ) ;
00800 RETURN(im);
00801 }
00802
00803
00804
00805
00806
00807
00808
00809
00810
00811
00812
00813
00814
00815
00816
00817 MRI_IMARR * mri_read_3D( char * tname )
00818 {
00819 int hglobal , himage , nx , ny , nz ;
00820 char fname[256] , buf[512] ;
00821 int ngood , length , kim , koff , datum_type , datum_len , swap ;
00822 MRI_IMARR * newar ;
00823 MRI_IMAGE * newim ;
00824 void * imar ;
00825 FILE * imfile ;
00826
00827 ENTRY("mri_read_3D") ;
00828
00829
00830
00831 if( tname == NULL || strlen(tname) < 10 ) RETURN(NULL) ;
00832
00833 switch( tname[2] ){
00834
00835 default:
00836 case ':':
00837 ngood = sscanf( tname , "3D:%d:%d:%d:%d:%d:%s" ,
00838 &hglobal , &himage , &nx , &ny , &nz , fname ) ;
00839
00840 swap = 0 ;
00841 datum_type = MRI_short ;
00842 datum_len = sizeof(short) ;
00843 break ;
00844
00845 case 's':
00846 ngood = sscanf( tname , "3Ds:%d:%d:%d:%d:%d:%s" ,
00847 &hglobal , &himage , &nx , &ny , &nz , fname ) ;
00848
00849 swap = 1 ;
00850 datum_type = MRI_short ;
00851 datum_len = sizeof(short) ;
00852 break ;
00853
00854 case 'b':
00855 ngood = sscanf( tname , "3Db:%d:%d:%d:%d:%d:%s" ,
00856 &hglobal , &himage , &nx , &ny , &nz , fname ) ;
00857
00858 swap = 0 ;
00859 datum_type = MRI_byte ;
00860 datum_len = sizeof(byte) ;
00861 break ;
00862
00863 case 'f':
00864 ngood = sscanf( tname , "3Df:%d:%d:%d:%d:%d:%s" ,
00865 &hglobal , &himage , &nx , &ny , &nz , fname ) ;
00866
00867 swap = 0 ;
00868 datum_type = MRI_float ;
00869 datum_len = sizeof(float) ;
00870 break ;
00871
00872 case 'd':
00873 ngood = sscanf( tname , "3Dd:%d:%d:%d:%d:%d:%s" ,
00874 &hglobal , &himage , &nx , &ny , &nz , fname ) ;
00875
00876 swap = 0 ;
00877 datum_type = MRI_double ;
00878 datum_len = sizeof(double) ;
00879 break ;
00880
00881 case 'i':
00882 ngood = sscanf( tname , "3Di:%d:%d:%d:%d:%d:%s" ,
00883 &hglobal , &himage , &nx , &ny , &nz , fname ) ;
00884
00885 swap = 0 ;
00886 datum_type = MRI_int ;
00887 datum_len = sizeof(int) ;
00888 break ;
00889
00890 case 'c':
00891 ngood = sscanf( tname , "3Dc:%d:%d:%d:%d:%d:%s" ,
00892 &hglobal , &himage , &nx , &ny , &nz , fname ) ;
00893
00894 swap = 0 ;
00895 datum_type = MRI_complex ;
00896 datum_len = sizeof(complex) ;
00897 break ;
00898
00899 case 'r':
00900 ngood = sscanf( tname , "3Dr:%d:%d:%d:%d:%d:%s" ,
00901 &hglobal , &himage , &nx , &ny , &nz , fname ) ;
00902
00903 swap = 0 ;
00904 datum_type = MRI_rgb ;
00905 datum_len = 3*sizeof(byte) ;
00906 break ;
00907 }
00908
00909 if( ngood < 6 || himage < 0 ||
00910 nx <= 0 || ny <= 0 || nz <= 0 ||
00911 strlen(fname) <= 0 ) RETURN(NULL);
00912
00913
00914
00915 if( strcmp(fname,"ALLZERO") == 0 ){
00916 INIT_IMARR(newar) ;
00917 for( kim=0 ; kim < nz ; kim++ ){
00918 newim = mri_new( nx , ny , datum_type ) ;
00919 imar = mri_data_pointer( newim ) ;
00920 memset( imar , 0 , newim->nvox * newim->pixel_size ) ;
00921 sprintf( buf , "%s#%d" , fname,kim ) ;
00922 mri_add_name( buf , newim ) ;
00923 ADDTO_IMARR(newar,newim) ;
00924 }
00925 RETURN(newar);
00926 }
00927
00928
00929
00930 imfile = fopen( fname , "r" ) ;
00931 if( imfile == NULL ){
00932 fprintf( stderr , "couldn't open image file %s\n" , fname ) ;
00933 RETURN(NULL);
00934 }
00935
00936 fseek( imfile , 0L , SEEK_END ) ;
00937 length = ftell( imfile ) ;
00938
00939
00940
00941
00942 #if 0
00943 if( hglobal < 0 ){
00944 hglobal = length - nz*(datum_len*nx*ny+himage) ;
00945 if( hglobal < 0 ) hglobal = 0 ;
00946 }
00947 #else
00948 if( hglobal == -1 || hglobal+himage < 0 ){
00949 hglobal = length - nz*(datum_len*nx*ny+himage) ;
00950 if( hglobal < 0 ) hglobal = 0 ;
00951 }
00952 #endif
00953
00954 ngood = hglobal + nz*(datum_len*nx*ny+himage) ;
00955 if( length < ngood ){
00956 fprintf( stderr,
00957 "image file %s is %d bytes long but must be at least %d bytes long\n"
00958 "for hglobal=%d himage=%d nx=%d ny=%d nz=%d and voxel=%d bytes\n",
00959 fname,length,ngood,hglobal,himage,nx,ny,nz,datum_len ) ;
00960 fclose( imfile ) ;
00961 RETURN(NULL);
00962 }
00963
00964
00965
00966 INIT_IMARR(newar) ;
00967
00968 for( kim=0 ; kim < nz ; kim++ ){
00969 koff = hglobal + (kim+1)*himage + datum_len*nx*ny*kim ;
00970 fseek( imfile , koff , SEEK_SET ) ;
00971
00972 newim = mri_new( nx , ny , datum_type ) ;
00973 imar = mri_data_pointer( newim ) ;
00974 length = fread( imar , datum_len , nx * ny , imfile ) ;
00975 if( swap ){
00976 mri_swapbytes( newim ) ;
00977 newim->was_swapped = 1 ;
00978 }
00979
00980 if( nz == 1 ) mri_add_name( fname , newim ) ;
00981 else {
00982 sprintf( buf , "%s#%d" , fname,kim ) ;
00983 mri_add_name( buf , newim ) ;
00984 }
00985
00986 ADDTO_IMARR(newar,newim) ;
00987 }
00988
00989 fclose(imfile) ;
00990 RETURN(newar);
00991 }
00992
00993
00994
00995
00996
00997
00998
00999
01000
01001
01002
01003
01004
01005
01006
01007
01008
01009
01010
01011
01012
01013
01014
01015
01016 MRI_IMARR * mri_read_file( char * fname )
01017 {
01018 MRI_IMARR * newar=NULL ;
01019 MRI_IMAGE * newim ;
01020 char * new_fname ;
01021
01022 ENTRY("mri_read_file") ;
01023
01024
01025
01026 new_fname = imsized_fname( fname ) ;
01027 if( new_fname == NULL ) RETURN( NULL );
01028
01029
01030
01031 if( strlen(new_fname) > 9 &&
01032 new_fname[0] == '3' &&
01033 new_fname[1] == 'D' &&
01034 (new_fname[2] == ':' || new_fname[3] == ':') ){
01035
01036 newar = mri_read_3D( new_fname ) ;
01037
01038 } else if( strlen(new_fname) > 9 &&
01039 new_fname[0] == '3' && new_fname[1] == 'A' && new_fname[3] == ':' ){
01040
01041 newar = mri_read_3A( new_fname ) ;
01042
01043 } else if( check_dicom_magic_num( new_fname ) ) {
01044
01045 newar = mri_read_dicom( new_fname );
01046
01047 } else if( strstr(new_fname,".hdr") != NULL ||
01048 strstr(new_fname,".HDR") != NULL ){
01049
01050 newar = mri_read_analyze75( new_fname ) ;
01051
01052 } else if( strstr(new_fname,".ima") != NULL ||
01053 strstr(new_fname,".IMA") != NULL ){
01054
01055 newar = mri_read_siemens( new_fname ) ;
01056
01057 } else if( strncmp(new_fname,"I.",2) == 0 ||
01058 strstr(new_fname,"/I.") != NULL ||
01059 strstr(new_fname,".ppm") != NULL ||
01060 strstr(new_fname,".pgm") != NULL ||
01061 strstr(new_fname,".pnm") != NULL ||
01062 strstr(new_fname,".PPM") != NULL ||
01063 strstr(new_fname,".PNM") != NULL ||
01064 strstr(new_fname,".PGM") != NULL ){
01065
01066 newim = mri_read( new_fname ) ;
01067
01068 if ( newim == NULL )
01069 newim = mri_read_ge4( new_fname ) ;
01070
01071 if( newim != NULL ){
01072 INIT_IMARR(newar) ;
01073 ADDTO_IMARR(newar,newim) ;
01074 }
01075
01076 } else if( strncmp(new_fname,"i.",2) == 0 ||
01077 strstr(new_fname,"/i.") != NULL ){
01078
01079 newim = mri_read_ge4( new_fname ) ;
01080
01081 if( newim != NULL ){
01082 INIT_IMARR(newar) ;
01083 ADDTO_IMARR(newar,newim) ;
01084 }
01085
01086 } else if( strstr(new_fname,".jpg" ) != NULL ||
01087 strstr(new_fname,".JPG" ) != NULL ||
01088 strstr(new_fname,".jpeg") != NULL ||
01089 strstr(new_fname,".JPEG") != NULL ||
01090 strstr(new_fname,".gif" ) != NULL ||
01091 strstr(new_fname,".GIF" ) != NULL ||
01092 strstr(new_fname,".tif" ) != NULL ||
01093 strstr(new_fname,".TIF" ) != NULL ||
01094 strstr(new_fname,".tiff") != NULL ||
01095 strstr(new_fname,".TIFF") != NULL ||
01096 strstr(new_fname,".bmp" ) != NULL ||
01097 strstr(new_fname,".BMP" ) != NULL ||
01098 strstr(new_fname,".pbm" ) != NULL ||
01099 strstr(new_fname,".PBM" ) != NULL ||
01100 strstr(new_fname,".png" ) != NULL ||
01101 strstr(new_fname,".PNG" ) != NULL ){
01102
01103 newim = mri_read_stuff( new_fname ) ;
01104 if( newim != NULL ){
01105 INIT_IMARR(newar) ;
01106 ADDTO_IMARR(newar,newim) ;
01107 }
01108
01109 } else if( strstr(new_fname,".mpg" ) != NULL ||
01110 strstr(new_fname,".MPG" ) != NULL ||
01111 strstr(new_fname,".mpeg") != NULL ||
01112 strstr(new_fname,".MPEG") != NULL ){
01113
01114 newar = mri_read_mpeg( new_fname ) ;
01115 }
01116
01117
01118
01119
01120 if( newar == NULL ){
01121
01122 if ( !AFNI_yesenv("AFNI_TRY_DICOM_LAST")) {
01123 newar = mri_read_dicom( new_fname ) ;
01124 }
01125
01126
01127
01128 if( newar == NULL ){
01129 newim = mri_read( new_fname ) ;
01130 if( newim == NULL ){ free(new_fname); RETURN( NULL ); }
01131 INIT_IMARR(newar) ;
01132 ADDTO_IMARR(newar,newim) ;
01133 }
01134
01135 if ( (newar == NULL) && AFNI_yesenv("AFNI_TRY_DICOM_LAST")) {
01136 newar = mri_read_dicom( new_fname ) ;
01137 }
01138 }
01139
01140 free(new_fname) ;
01141
01142
01143
01144 if( newar != NULL && newar->num > 0 ){
01145 int ii ;
01146 for( ii=0 ; ii < newar->num ; ii++ ){
01147 newim = IMARR_SUBIM(newar,ii) ;
01148 if( newim != NULL && newim->fname == NULL )
01149 newim->fname = strdup(fname) ;
01150 }
01151 }
01152
01153 RETURN( newar );
01154 }
01155
01156
01157
01158
01159
01160
01161
01162
01163
01164 MRI_IMAGE * mri_read_just_one( char * fname )
01165 {
01166 MRI_IMARR * imar ;
01167 MRI_IMAGE * im ;
01168 char * new_fname ;
01169
01170 ENTRY("mri_read_just_one") ;
01171
01172 new_fname = imsized_fname( fname ) ;
01173 if( new_fname == NULL ) RETURN( NULL );
01174
01175 imar = mri_read_file( new_fname ) ; free(new_fname) ;
01176 if( imar == NULL ) RETURN( NULL );
01177 if( imar->num != 1 ){ DESTROY_IMARR(imar) ; RETURN( NULL ); }
01178 im = IMAGE_IN_IMARR(imar,0) ;
01179 FREE_IMARR(imar) ;
01180 RETURN( im );
01181 }
01182
01183
01184
01185
01186
01187
01188
01189
01190
01191
01192
01193
01194
01195 static int mri_imcount_analyze75( char * ) ;
01196 static int mri_imcount_siemens( char * ) ;
01197
01198 int mri_imcount( char * tname )
01199 {
01200 int hglobal , himage , nx , ny , nz , ngood ;
01201 char fname[256]="\0" ;
01202 char * new_fname ;
01203
01204 ENTRY("mri_imcount") ;
01205
01206 if( tname == NULL ) RETURN( 0 );
01207 new_fname = imsized_fname( tname ) ;
01208 if( new_fname == NULL ) RETURN( 0 );
01209
01210
01211
01212 if( strlen(new_fname) > 9 && new_fname[0] == '3' && new_fname[1] == 'D' &&
01213 (new_fname[2] == ':' || new_fname[3] == ':') ){
01214
01215 switch( new_fname[2] ){
01216
01217 default:
01218 case ':':
01219 ngood = sscanf( new_fname , "3D:%d:%d:%d:%d:%d:%s" ,
01220 &hglobal , &himage , &nx , &ny , &nz , fname ) ;
01221 break ;
01222
01223 case 's':
01224 ngood = sscanf( new_fname , "3Ds:%d:%d:%d:%d:%d:%s" ,
01225 &hglobal , &himage , &nx , &ny , &nz , fname ) ;
01226 break ;
01227
01228 case 'b':
01229 ngood = sscanf( new_fname , "3Db:%d:%d:%d:%d:%d:%s" ,
01230 &hglobal , &himage , &nx , &ny , &nz , fname ) ;
01231 break ;
01232
01233 case 'f':
01234 ngood = sscanf( new_fname , "3Df:%d:%d:%d:%d:%d:%s" ,
01235 &hglobal , &himage , &nx , &ny , &nz , fname ) ;
01236 break ;
01237
01238 case 'd':
01239 ngood = sscanf( new_fname , "3Dd:%d:%d:%d:%d:%d:%s" ,
01240 &hglobal , &himage , &nx , &ny , &nz , fname ) ;
01241 break ;
01242
01243 case 'i':
01244 ngood = sscanf( new_fname , "3Di:%d:%d:%d:%d:%d:%s" ,
01245 &hglobal , &himage , &nx , &ny , &nz , fname ) ;
01246 break ;
01247
01248 case 'c':
01249 ngood = sscanf( new_fname , "3Dc:%d:%d:%d:%d:%d:%s" ,
01250 &hglobal , &himage , &nx , &ny , &nz , fname ) ;
01251 break ;
01252
01253 case 'r':
01254 ngood = sscanf( new_fname , "3Dr:%d:%d:%d:%d:%d:%s" ,
01255 &hglobal , &himage , &nx , &ny , &nz , fname ) ;
01256 break ;
01257 }
01258
01259 free( new_fname ) ;
01260 if( ngood < 6 || himage < 0 ||
01261 nx <= 0 || ny <= 0 || nz <= 0 ||
01262 strlen(fname) <= 0 ) RETURN( 0 );
01263 else RETURN( nz );
01264 }
01265
01266
01267
01268 if( strlen(new_fname) > 9 &&
01269 new_fname[0] == '3' && new_fname[1] == 'A' && new_fname[3] == ':' ){
01270
01271 switch( new_fname[2] ){
01272
01273 default: ngood = 0 ; break ;
01274
01275 case 's':
01276 ngood = sscanf( new_fname, "3As:%d:%d:%d:%s", &nx, &ny, &nz, fname ) ;
01277 break ;
01278
01279 case 'b':
01280 ngood = sscanf( new_fname, "3Ab:%d:%d:%d:%s", &nx, &ny, &nz, fname ) ;
01281 break ;
01282
01283 case 'f':
01284 ngood = sscanf( new_fname, "3Af:%d:%d:%d:%s", &nx, &ny, &nz, fname ) ;
01285 break ;
01286 }
01287
01288 free( new_fname ) ;
01289 if( ngood < 4 || nx <= 0 || ny <= 0 || nz <= 0 || strlen(fname) <= 0 ) RETURN( 0 );
01290 else RETURN( nz );
01291 }
01292
01293
01294
01295 if( strstr(new_fname,".hdr") != NULL ||
01296 strstr(new_fname,".HDR") != NULL ){
01297
01298 nz = mri_imcount_analyze75( new_fname ) ;
01299 if( nz > 0 ){ free(new_fname); RETURN(nz); }
01300 }
01301
01302 if( strstr(new_fname,".ima") != NULL ||
01303 strstr(new_fname,".IMA") != NULL ){
01304
01305 nz = mri_imcount_siemens( new_fname ) ;
01306 if( nz > 0 ){ free(new_fname); RETURN(nz); }
01307 }
01308
01309 if( strstr(new_fname,".mpg" ) != NULL ||
01310 strstr(new_fname,".MPG" ) != NULL ||
01311 strstr(new_fname,".mpeg") != NULL ||
01312 strstr(new_fname,".MPEG") != NULL ){
01313
01314 nz = mri_imcount_mpeg( new_fname ) ;
01315 if( nz > 0 ){ free(new_fname); RETURN(nz); }
01316 }
01317
01318
01319
01320 mri_dicom_seterr(0) ;
01321 nz = mri_imcount_dicom( new_fname ) ;
01322 mri_dicom_seterr(1) ;
01323 if( nz > 0 ){ free(new_fname); RETURN(nz); }
01324
01325
01326
01327 free(new_fname) ; RETURN(1) ;
01328 }
01329
01330
01331
01332
01333
01334
01335
01336
01337
01338
01339
01340
01341 MRI_IMARR * mri_read_many_files( int nf , char * fn[] )
01342 {
01343 MRI_IMARR * newar , * outar ;
01344 int kf , ii ;
01345
01346 ENTRY("mri_read_many_files") ;
01347
01348 if( nf <= 0 ) RETURN( NULL );
01349 INIT_IMARR(outar) ;
01350
01351 for( kf=0 ; kf < nf ; kf++ ){
01352 newar = mri_read_file( fn[kf] ) ;
01353
01354 if( newar == NULL ){
01355 fprintf(stderr,"cannot read images from file %s\n",fn[kf]) ;
01356 for( ii=0 ; ii < outar->num ; ii++ ) mri_free(outar->imarr[ii]) ;
01357 FREE_IMARR(outar) ;
01358 RETURN( NULL );
01359 }
01360
01361 for( ii=0 ; ii < newar->num ; ii++ )
01362 ADDTO_IMARR( outar , newar->imarr[ii] ) ;
01363
01364 FREE_IMARR(newar) ;
01365 }
01366 RETURN( outar );
01367 }
01368
01369
01370
01371
01372
01373
01374
01375
01376 MRI_IMARR * mri_read_ppm3( char * fname )
01377 {
01378 int ch , nch , nx,ny,maxval , length , npix,ii ;
01379 char buf[512] ;
01380 MRI_IMAGE *rim , *gim , *bim ;
01381 MRI_IMARR * outar ;
01382 FILE * imfile ;
01383 byte * rby , * gby , * bby , * rgby ;
01384
01385 ENTRY("mri_read_ppm3") ;
01386
01387
01388
01389 imfile = fopen( fname , "r" ) ;
01390 if( imfile == NULL ){
01391 fprintf(stderr,"couldn't open file %s in mri_read_ppm3\n",fname); RETURN(NULL) ;
01392 }
01393
01394
01395
01396 ch = getc( imfile ) ; if( ch != 'P' ) { fclose(imfile) ; RETURN(NULL); }
01397 ch = getc( imfile ) ; if( ch != '6' ) { fclose(imfile) ; RETURN(NULL); }
01398
01399
01400
01401 ch = getc(imfile) ;
01402
01403 NUMSCAN(nx) ; if( nx <= 0 ) { fclose(imfile) ; RETURN(NULL); }
01404 NUMSCAN(ny) ; if( ny <= 0 ) { fclose(imfile) ; RETURN(NULL); }
01405 NUMSCAN(maxval) ; if( maxval <= 0 ||
01406 maxval > 255 ) { fclose(imfile) ; RETURN(NULL); }
01407
01408
01409
01410 rim = mri_new( nx , ny , MRI_byte ) ; rby = mri_data_pointer( rim ) ;
01411 gim = mri_new( nx , ny , MRI_byte ) ; gby = mri_data_pointer( gim ) ;
01412 bim = mri_new( nx , ny , MRI_byte ) ; bby = mri_data_pointer( bim ) ;
01413
01414 sprintf(buf,"%s#R",fname) ; mri_add_name( buf , rim ) ;
01415 sprintf(buf,"%s#G",fname) ; mri_add_name( buf , gim ) ;
01416 sprintf(buf,"%s#B",fname) ; mri_add_name( buf , bim ) ;
01417
01418 rgby = (byte *) malloc( sizeof(byte) * 3*nx*ny ) ;
01419 if( rgby == NULL ){
01420 fprintf(stderr,"couldn't malloc workspace in mri_read_ppm3!\n") ; EXIT(1) ;
01421 }
01422
01423
01424
01425 length = fread( rgby , sizeof(byte) , 3*nx*ny , imfile ) ;
01426 fclose( imfile ) ;
01427
01428 if( length != 3*nx*ny ){
01429 free(rgby) ; mri_free(rim) ; mri_free(gim) ; mri_free(bim) ;
01430 fprintf(stderr,"couldn't read data from file %s in mri_read_ppm3\n",fname) ;
01431 RETURN(NULL);
01432 }
01433
01434
01435
01436 npix = nx*ny ;
01437 for( ii=0 ; ii < npix ; ii++ ){
01438 rby[ii] = rgby[3*ii ] ;
01439 gby[ii] = rgby[3*ii+1] ;
01440 bby[ii] = rgby[3*ii+2] ;
01441 }
01442 free( rgby ) ;
01443
01444
01445
01446 INIT_IMARR(outar) ;
01447 ADDTO_IMARR( outar , rim ) ;
01448 ADDTO_IMARR( outar , gim ) ;
01449 ADDTO_IMARR( outar , bim ) ;
01450 RETURN(outar);
01451 }
01452
01453
01454
01455
01456
01457
01458
01459
01460
01461
01462
01463
01464
01465 MRI_IMAGE *mri_read_nsize( char * fname )
01466 {
01467 MRI_IMARR *imar ;
01468 MRI_IMAGE *imout ;
01469
01470 imar = mri_read_file( fname ) ;
01471 if( imar == NULL ) return NULL ;
01472 if( imar->num != 1 ){ DESTROY_IMARR(imar) ; return NULL ; }
01473
01474 imout = mri_nsize( IMAGE_IN_IMARR(imar,0) ) ;
01475 mri_add_name( IMAGE_IN_IMARR(imar,0)->name , imout ) ;
01476
01477 DESTROY_IMARR(imar) ;
01478 return imout ;
01479 }
01480
01481
01482
01483 MRI_IMARR *mri_read_many_nsize( int nf , char * fn[] )
01484 {
01485 MRI_IMARR * newar , * outar ;
01486 MRI_IMAGE * im ;
01487 int ii ;
01488
01489 newar = mri_read_many_files( nf , fn ) ;
01490 if( newar == NULL ) return NULL ;
01491
01492 INIT_IMARR(outar) ;
01493 for( ii=0 ; ii < newar->num ; ii++ ){
01494 im = mri_nsize( IMAGE_IN_IMARR(newar,ii) ) ;
01495 mri_add_name( IMAGE_IN_IMARR(newar,ii)->name , im ) ;
01496 ADDTO_IMARR(outar,im) ;
01497 mri_free( IMAGE_IN_IMARR(newar,ii) ) ;
01498 }
01499 FREE_IMARR(newar) ;
01500 return outar ;
01501 }
01502
01503
01504
01505
01506
01507
01508
01509
01510
01511
01512
01513
01514 void init_MCW_sizes(void)
01515 {
01516 int num , count ;
01517 char ename[32] ;
01518 char * str ;
01519
01520 if( MCW_imsize_good >= 0 ) return ;
01521
01522 MCW_imsize_good = 0 ;
01523
01524 for( num=0 ; num < MAX_MCW_IMSIZE ; num++ ){
01525
01526 imsize[num].size = -1 ;
01527
01528
01529
01530 sprintf( ename , "AFNI_IMSIZE_%d" , num+1 ) ;
01531 str = my_getenv( ename ) ;
01532
01533 if( str == NULL ){
01534 sprintf( ename , "MCW_IMSIZE_%d" , num+1 ) ;
01535 str = my_getenv( ename ) ;
01536 if( str == NULL ) continue ;
01537 }
01538
01539 imsize[num].prefix = (char *) malloc( sizeof(char) * strlen(str) ) ;
01540 if( imsize[num].prefix == NULL ){
01541 fprintf(stderr,"\n*** Can't malloc in init_MCW_sizes! ***\a\n");
01542 EXIT(1) ;
01543 }
01544
01545 if( str[0] != '%' ){
01546
01547 imsize[num].head = -1 ;
01548 count = sscanf( str , "%d=%s" , &(imsize[num].size) , imsize[num].prefix ) ;
01549 if( count != 2 || imsize[num].size < 2 || strlen(imsize[num].prefix) < 2 ){
01550 free( imsize[num].prefix ) ;
01551 fprintf(stderr,"bad environment %s = %s\n" ,
01552 ename , str ) ;
01553 }
01554
01555 } else {
01556
01557 count = sscanf( str+1 , "%d+%d=%s" ,
01558 &(imsize[num].size) , &(imsize[num].head) , imsize[num].prefix ) ;
01559
01560 if( count != 3 || imsize[num].size < 2 ||
01561 imsize[num].head < 0 || strlen(imsize[num].prefix) < 2 ){
01562
01563 free( imsize[num].prefix ) ;
01564 fprintf(stderr,"bad environment %s = %s\n" ,
01565 ename , str ) ;
01566 }
01567 }
01568
01569 MCW_imsize_good ++ ;
01570 }
01571
01572 return ;
01573 }
01574
01575
01576
01577
01578 char * my_strdup( char * str )
01579 {
01580 char * new_str ;
01581 if( str == NULL ) return NULL ;
01582 new_str = (char *) malloc( sizeof(char) * (strlen(str)+1) ) ;
01583 if( new_str != NULL ) strcpy( new_str , str ) ;
01584 return new_str ;
01585 }
01586
01587
01588
01589
01590
01591
01592
01593
01594
01595
01596
01597 char * imsized_fname( char * fname )
01598 {
01599 int num , lll ;
01600 long len ;
01601 char * new_name ;
01602
01603 init_MCW_sizes() ;
01604 if( MCW_imsize_good == 0 ){
01605 new_name = my_strdup(fname) ;
01606 return new_name ;
01607 }
01608
01609 len = mri_filesize( fname ) ;
01610 if( len <= 0 ){
01611 new_name = my_strdup(fname) ;
01612 return new_name ;
01613 }
01614
01615 for( num=0 ; num < MAX_MCW_IMSIZE ; num++ ){
01616
01617 if( imsize[num].size <= 0 ) continue ;
01618
01619 if( imsize[num].head < 0 && len == imsize[num].size ){
01620
01621 lll = strlen(fname) + strlen(imsize[num].prefix) + 4 ;
01622 new_name = (char *) malloc( sizeof(char) * lll ) ;
01623 if( new_name == NULL ){
01624 fprintf(stderr,"\n*** Can't malloc in imsized_fname! ***\a\n");
01625 EXIT(1) ;
01626 }
01627 sprintf( new_name , "%s%s" , imsize[num].prefix , fname ) ;
01628 return new_name ;
01629
01630 } else if( (len-imsize[num].head) % imsize[num].size == 0 ){
01631 int count = (len-imsize[num].head) / imsize[num].size ;
01632
01633 if( count < 1 ) continue ;
01634
01635 lll = strlen(fname) + strlen(imsize[num].prefix) + 32 ;
01636 new_name = (char *) malloc( sizeof(char) * lll ) ;
01637 if( new_name == NULL ){
01638 fprintf(stderr,"\n*** Can't malloc in imsized_fname! ***\a\n");
01639 EXIT(1) ;
01640 }
01641 sprintf( new_name , "%s%d:%s" , imsize[num].prefix , count , fname ) ;
01642 return new_name ;
01643 }
01644
01645 }
01646
01647 new_name = my_strdup(fname) ;
01648 return new_name ;
01649 }
01650
01651
01652
01653
01654
01655
01656
01657
01658
01659 long mri_filesize( char * pathname )
01660 {
01661 static struct stat buf ;
01662 int ii ;
01663
01664 if( pathname == NULL ) return -1 ;
01665 ii = stat( pathname , &buf ) ; if( ii != 0 ) return -1 ;
01666 return buf.st_size ;
01667 }
01668
01669
01670
01671
01672
01673
01674
01675
01676
01677
01678
01679
01680 void mri_read_ppm_header( char *fname , int *nx, int *ny )
01681 {
01682 FILE *imfile ;
01683 int ch , nch , nxx,nyy ;
01684 char buf[256] ;
01685
01686 ENTRY("mri_read_ppm_header") ;
01687
01688 if( fname == NULL || nx == NULL || ny == NULL ) EXRETURN ;
01689
01690 *nx = *ny = 0 ;
01691
01692
01693
01694 imfile = fopen( fname , "r" ) ; if( imfile == NULL ) EXRETURN ;
01695
01696
01697
01698 ch = getc( imfile ) ; if( ch != 'P' ) { fclose(imfile) ; EXRETURN ; }
01699 ch = getc( imfile ) ; if( ch != '6' ) { fclose(imfile) ; EXRETURN ; }
01700
01701
01702
01703 ch = getc(imfile) ;
01704
01705 NUMSCAN(nxx) ; if( nxx <= 0 ){ fclose(imfile) ; EXRETURN ; }
01706 NUMSCAN(nyy) ; if( nyy <= 0 ){ fclose(imfile) ; EXRETURN ; }
01707
01708
01709
01710 fclose(imfile) ; *nx = nxx ; *ny = nyy ; EXRETURN ;
01711 }
01712
01713
01714
01715
01716
01717
01718
01719
01720
01721
01722 MRI_IMAGE * mri_read_ppm( char * fname )
01723 {
01724 int ch , nch , nx,ny,maxval , length ;
01725 MRI_IMAGE * rgbim ;
01726 FILE * imfile ;
01727 byte * rgby ;
01728 char buf[256] ;
01729
01730 ENTRY("mri_read_ppm") ;
01731
01732
01733
01734 imfile = fopen( fname , "r" ) ;
01735 if( imfile == NULL ) RETURN(NULL);
01736
01737
01738
01739 ch = getc( imfile ) ; if( ch != 'P' ) { fclose(imfile) ; RETURN(NULL); }
01740 ch = getc( imfile ) ; if( ch != '6' ) { fclose(imfile) ; RETURN(NULL); }
01741
01742
01743
01744 ch = getc(imfile) ;
01745
01746 NUMSCAN(nx) ; if( nx <= 0 ) { fclose(imfile); RETURN(NULL); }
01747 NUMSCAN(ny) ; if( ny <= 0 ) { fclose(imfile); RETURN(NULL); }
01748 NUMSCAN(maxval); if( maxval <= 0 ||
01749 maxval > 255 ){ fclose(imfile); RETURN(NULL); }
01750
01751
01752
01753 rgbim = mri_new( nx , ny , MRI_rgb ) ; mri_add_name( fname , rgbim ) ;
01754 rgby = MRI_RGB_PTR(rgbim) ;
01755
01756
01757
01758 length = fread( rgby , sizeof(byte) , 3*nx*ny , imfile ) ;
01759 fclose( imfile ) ;
01760
01761 if( length != 3*nx*ny ){ mri_free(rgbim) ; RETURN(NULL) ; }
01762
01763
01764
01765 if( maxval < 255 ){
01766 int ii ; float fac = 255.4/maxval ;
01767 for( ii=0 ; ii < 3*nx*ny ; ii++ ) rgby[ii] = (byte)( rgby[ii]*fac ) ;
01768 }
01769
01770 RETURN(rgbim) ;
01771 }
01772
01773
01774
01775
01776 #define LBUF 2524288
01777
01778
01779 #define FRB(b) do{ if( (b)!=NULL ){free((b)); (b)=NULL;} }while(0)
01780
01781 #undef USE_LASTBUF
01782
01783
01784
01785
01786
01787
01788
01789
01790
01791
01792 static char * my_fgets( char *buf , int size , FILE *fts )
01793 {
01794 char *ptr ;
01795 int nbuf , ll,ii , cflag ;
01796 static char *qbuf=NULL ;
01797
01798 #ifdef USE_LASTBUF
01799 static char *lastbuf = NULL ;
01800 static int nlastbuf = 0 ;
01801
01802 if( buf == NULL && lastbuf != NULL ){
01803 free((void *)lastbuf); lastbuf = NULL; nlastbuf = 0 ;
01804 }
01805 #endif
01806
01807 if( buf == NULL && qbuf != NULL ){ free((void *)qbuf); qbuf = NULL; }
01808
01809 if( buf == NULL || size < 1 || fts == NULL ) return NULL ;
01810
01811 if( qbuf == NULL ) qbuf = AFMALL(char, LBUF) ;
01812
01813 nbuf = 0 ;
01814 cflag = 0 ;
01815
01816 while(1){
01817
01818 ptr = fgets( qbuf , LBUF , fts ) ;
01819
01820 if( ptr == NULL ) break ;
01821
01822
01823
01824 for( ; *ptr != '\0' && isspace(*ptr) ; ptr++ ) ;
01825
01826
01827
01828 if( *ptr == '\0' ){ if(cflag) break; else continue; }
01829
01830 #ifdef USE_LASTBUF
01831
01832
01833 if( *ptr == '"' && *(ptr+1) == '"' && nlastbuf > 0 && nbuf == 0 ){
01834 ll = strlen(lastbuf) ; if( ll >= size ) ll = size-1 ;
01835 memcpy(buf,lastbuf,ll-1) ; buf[ll] = '\0' ;
01836 return buf ;
01837 }
01838 #endif
01839
01840
01841
01842 if( *ptr == '#' || (*ptr == '/' && *(ptr+1) == '/') ) continue ;
01843
01844
01845
01846 ll = strlen(ptr) ;
01847 for( ii=ll-1 ; isspace(ptr[ii]) && ii > 0 ; ii-- )
01848 ptr[ii] = '\0' ;
01849
01850 ll = strlen(ptr) ;
01851 if( ll == 0 ) continue ;
01852
01853 cflag = (ptr[ll-1] == '\\') ;
01854 if( cflag ) ptr[ll-1] = ' ' ;
01855
01856
01857
01858 if( nbuf+ll+1 > size ){
01859 ll = size - (nbuf+1) ;
01860 if( ll <= 0 ) break ;
01861 }
01862
01863 memcpy(buf+nbuf,ptr,ll+1) ; nbuf += ll ;
01864 if( !cflag ) break ;
01865
01866 }
01867
01868 #ifdef LASTBUF
01869
01870
01871 ll = strlen(buf) ;
01872 if( ll+1 > nlastbuf ){
01873 nlastbuf = ll+2 ; lastbuf = (char *)realloc((void *)lastbuf,nlastbuf) ;
01874 }
01875 memcpy(lastbuf,buf,ll+1) ;
01876 #endif
01877
01878
01879
01880 if( nbuf > 0 ) return buf ;
01881 return NULL ;
01882 }
01883
01884
01885 static float lbfill = 0.0 ;
01886
01887
01888
01889
01890 static floatvec * decode_linebuf( char *buf )
01891 {
01892 floatvec *fv=NULL ;
01893 int blen, bpos, ncol, ii, count ;
01894 char sep, vbuf[64] , *cpt ;
01895 float val ;
01896
01897 if( buf == NULL || *buf == '\0' ) return fv ;
01898
01899 blen = strlen(buf) ;
01900 ncol = 0 ;
01901
01902
01903
01904 for( ii=0 ; ii < blen ; ii++ ) if( buf[ii] == ',' ) buf[ii] = ' ' ;
01905
01906 fv = (floatvec *)malloc(sizeof(floatvec)) ;
01907 fv->nar = 0 ;
01908 fv->ar = (float *)NULL ;
01909
01910 for( bpos=0 ; bpos < blen ; ){
01911
01912
01913 for( ; bpos < blen && (isspace(buf[bpos])||buf[bpos]==',') ; bpos++ ) ;
01914 if( bpos == blen ) break ;
01915
01916 sscanf( buf+bpos , "%63s" , vbuf ) ;
01917
01918 val = 0.0 ; count = 1 ;
01919 if( vbuf[0] == '*' ){
01920 val = lbfill ;
01921 } else if( (cpt=strchr(vbuf,'@')) != NULL ){
01922 sscanf( vbuf , "%d%c%f" , &count , &sep , &val ) ;
01923 if( count < 1 ) count = 1 ;
01924 if( *(cpt+1) == '*' ) val = lbfill ;
01925 } else {
01926 sscanf( vbuf , "%f" , &val ) ;
01927 }
01928
01929 fv->ar = (float *)realloc( (void *)fv->ar , sizeof(float)*(fv->nar+count) ) ;
01930 for( ii=0 ; ii < count ; ii++ ) fv->ar[ii+fv->nar] = val ;
01931 fv->nar += count ;
01932 bpos += strlen(vbuf) ;
01933 }
01934
01935 if( fv->nar == 0 ){ KILL_floatvec(fv); fv = NULL; }
01936 return fv ;
01937 }
01938
01939
01940
01941
01942 #define INC_TSARSIZE 128
01943
01944
01945
01946
01947
01948
01949
01950
01951
01952
01953
01954
01955
01956
01957
01958
01959
01960
01961
01962
01963
01964
01965 MRI_IMAGE * mri_read_ascii( char * fname )
01966 {
01967 MRI_IMAGE * outim ;
01968 int ii,jj,val , used_tsar , alloc_tsar ;
01969 float * tsar ;
01970 float ftemp ;
01971 FILE * fts ;
01972 char * ptr ;
01973 int ncol , bpos , blen , nrow ;
01974 static char *buf=NULL ;
01975
01976 floatvec *fvec ;
01977 int incts ;
01978
01979 ENTRY("mri_read_ascii") ;
01980
01981 if( fname == NULL || fname[0] == '\0' ) RETURN(NULL) ;
01982
01983 if( strncmp(fname,"1D:",3) == 0 ){
01984 MRI_IMAGE *qim = mri_1D_fromstring( fname+3 ) ;
01985 if( qim != NULL ){
01986 outim = mri_transpose(qim); mri_free(qim); RETURN(outim);
01987 }
01988 }
01989
01990 fts = fopen( fname , "r" ); if( fts == NULL ) RETURN(NULL);
01991
01992 if( buf == NULL ) buf = AFMALL(char, LBUF) ;
01993
01994
01995
01996
01997 (void) my_fgets( NULL , 0 , NULL ) ;
01998 ptr = my_fgets( buf , LBUF , fts ) ;
01999 if( ptr==NULL || *ptr=='\0' ){ FRB(buf); fclose(fts); RETURN(NULL); }
02000
02001 lbfill = 0.0f ;
02002
02003 fvec = decode_linebuf( buf ) ;
02004 if( fvec == NULL || fvec->nar == 0 ){
02005 if( fvec != NULL ) KILL_floatvec(fvec) ;
02006 FRB(buf); fclose(fts); RETURN(NULL);
02007 }
02008 ncol = fvec->nar ; KILL_floatvec(fvec) ;
02009
02010
02011
02012 rewind( fts ) ;
02013
02014 incts = MAX(INC_TSARSIZE,ncol) ;
02015 used_tsar = 0 ;
02016 alloc_tsar = incts ;
02017 tsar = (float *) malloc( sizeof(float) * alloc_tsar ) ;
02018 if( tsar == NULL ){
02019 fprintf(stderr,"\n*** malloc error in mri_read_ascii ***\n"); EXIT(1);
02020 }
02021
02022
02023
02024 nrow = 0 ;
02025 while( 1 ){
02026 ptr = my_fgets( buf , LBUF , fts ) ;
02027 if( ptr==NULL || *ptr=='\0' ) break ;
02028
02029 fvec = decode_linebuf( buf ) ;
02030 if( fvec == NULL ) break ;
02031 if( fvec->nar == 0 ){ KILL_floatvec(fvec); break; }
02032
02033 if( used_tsar + ncol >= alloc_tsar ){
02034 alloc_tsar += incts ;
02035 tsar = (float *)realloc( (void *)tsar,sizeof(float)*alloc_tsar );
02036 if( tsar == NULL ){
02037 fprintf(stderr,"\n*** realloc error in mri_read_ascii ***\n"); EXIT(1);
02038 }
02039 }
02040 for( ii=0 ; ii < fvec->nar && ii < ncol ; ii++ )
02041 tsar[used_tsar+ii] = fvec->ar[ii] ;
02042 for( ; ii < ncol ; ii++ )
02043 tsar[used_tsar+ii] = 0.0 ;
02044 used_tsar += ncol ;
02045 KILL_floatvec(fvec) ;
02046
02047 nrow++ ;
02048 }
02049 fclose( fts ) ;
02050 (void) my_fgets( NULL , 0 , NULL ) ;
02051
02052 if( used_tsar <= 1 ){ FRB(buf); free(tsar); RETURN(NULL); }
02053
02054 tsar = (float *) realloc( tsar , sizeof(float) * used_tsar ) ;
02055 if( tsar == NULL ){
02056 fprintf(stderr,"\n*** final realloc error in mri_read_ascii ***\n"); EXIT(1);
02057 }
02058
02059 outim = mri_new_vol_empty( ncol , nrow , 1 , MRI_float ) ;
02060 mri_fix_data_pointer( tsar , outim ) ;
02061 mri_add_name( fname , outim ) ;
02062
02063 FRB(buf) ; RETURN(outim) ;
02064 }
02065
02066
02067
02068
02069
02070
02071
02072
02073
02074
02075
02076
02077
02078
02079 MRI_IMAGE * mri_read_1D( char * fname )
02080 {
02081 MRI_IMAGE *inim , *outim , *flim ;
02082 char dname[512] , *cpt , *dpt ;
02083 int ii,jj,nx,ny,nts , *ivlist , *ivl , *sslist ;
02084 float *far , *oar ;
02085
02086 ENTRY("mri_read_1D") ;
02087
02088 if( fname == NULL || fname[0] == '\0' || strlen(fname) > 511 ) RETURN(NULL) ;
02089
02090 if( strncmp(fname,"1D:",3) == 0 ){
02091 return mri_1D_fromstring( fname+3 ) ;
02092 }
02093
02094
02095
02096 cpt = strstr(fname,"[") ;
02097 dpt = strstr(fname,"{") ;
02098
02099 if( cpt == fname || dpt == fname ){
02100 fprintf(stderr,"*** Illegal filename in mri_read_1D: %s\n",fname) ;
02101 RETURN(NULL) ;
02102 } else {
02103 strcpy( dname , fname ) ;
02104 if( cpt != NULL ){ ii = cpt-fname; dname[ii] = '\0'; }
02105 if( dpt != NULL ){ ii = dpt-fname; dname[ii] = '\0'; }
02106 }
02107
02108
02109
02110 inim = mri_read_ascii(dname) ;
02111 if( inim == NULL ) RETURN(NULL) ;
02112 flim = mri_transpose(inim) ; mri_free(inim) ;
02113
02114
02115
02116 nx = flim->nx ; ny = flim->ny ;
02117
02118 ivlist = MCW_get_intlist( ny , cpt ) ;
02119 sslist = MCW_get_intlist( nx , dpt ) ;
02120
02121
02122
02123 if( ivlist != NULL && ivlist[0] > 0 ){
02124 nts = ivlist[0] ;
02125 ivl = ivlist + 1 ;
02126
02127 for( ii=0 ; ii < nts ; ii++ ){
02128 if( ivl[ii] < 0 || ivl[ii] >= ny ){
02129 fprintf(stderr,"*** Out-of-range subvector [list] in mri_read_1D: %s\n",fname) ;
02130 mri_free(flim) ; free(ivlist) ; RETURN(NULL) ;
02131 }
02132 }
02133
02134 outim = mri_new( nx , nts , MRI_float ) ;
02135 far = MRI_FLOAT_PTR( flim ) ;
02136 oar = MRI_FLOAT_PTR( outim ) ;
02137
02138 for( ii=0 ; ii < nts ; ii++ )
02139 memcpy( oar + ii*nx , far + ivl[ii]*nx , sizeof(float)*nx ) ;
02140
02141 mri_free(flim); free(ivlist); flim = outim; ny = nts;
02142 }
02143
02144
02145
02146 if( sslist != NULL && sslist[0] > 0 ){
02147 nts = sslist[0] ;
02148 ivl = sslist + 1 ;
02149
02150 for( ii=0 ; ii < nts ; ii++ ){
02151 if( ivl[ii] < 0 || ivl[ii] >= nx ){
02152 fprintf(stderr,"*** Out-of-range subsampling {list} in mri_read_1D: %s\n",fname) ;
02153 mri_free(flim) ; free(sslist) ; RETURN(NULL) ;
02154 }
02155 }
02156
02157 outim = mri_new( nts , ny , MRI_float ) ;
02158 far = MRI_FLOAT_PTR( flim ) ;
02159 oar = MRI_FLOAT_PTR( outim ) ;
02160
02161 for( ii=0 ; ii < nts ; ii++ )
02162 for( jj=0 ; jj < ny ; jj++ )
02163 oar[ii+jj*nts] = far[ivl[ii]+jj*nx] ;
02164
02165 mri_free(flim); free(sslist); flim = outim;
02166 }
02167
02168 RETURN(flim) ;
02169 }
02170
02171
02172
02173 MRI_IMAGE * mri_read_ascii_ragged( char *fname , float filler )
02174 {
02175 MRI_IMAGE *outim ;
02176 int ii,jj , ncol,nrow ;
02177 float *tsar ;
02178 FILE *fts ;
02179 char *ptr ;
02180 static char *buf=NULL ;
02181 floatvec *fvec ;
02182
02183 ENTRY("mri_read_ascii_ragged") ;
02184
02185 if( fname == NULL || *fname == '\0' ){ FRB(buf); RETURN(NULL); }
02186
02187 fts = fopen( fname , "r" ); if( fts == NULL ){ FRB(buf); RETURN(NULL); }
02188
02189 if( buf == NULL ) buf = AFMALL(char, LBUF) ;
02190
02191
02192
02193
02194 lbfill = filler ;
02195
02196 (void) my_fgets( NULL , 0 , NULL ) ;
02197 ncol = nrow = 0 ;
02198 while(1){
02199 ptr = my_fgets( buf , LBUF , fts ) ;
02200 if( ptr==NULL || *ptr=='\0' ) break ;
02201 fvec = decode_linebuf( buf ) ;
02202 if( fvec != NULL && fvec->nar > 0 ){ nrow++; ncol = MAX(ncol,fvec->nar); }
02203 if( fvec != NULL ) KILL_floatvec(fvec) ; else break ;
02204 }
02205 if( nrow == 0 || ncol == 0 ){ fclose(fts); FRB(buf); lbfill=0.0f; RETURN(NULL); }
02206
02207
02208
02209 rewind( fts ) ;
02210
02211 outim = mri_new( ncol , nrow , MRI_float ) ;
02212 tsar = MRI_FLOAT_PTR(outim) ;
02213
02214
02215
02216 nrow = 0 ;
02217 while( 1 ){
02218 ptr = my_fgets( buf , LBUF , fts ) ;
02219 if( ptr==NULL || *ptr=='\0' ) break ;
02220
02221 fvec = decode_linebuf( buf ) ;
02222 if( fvec == NULL ) break ;
02223 if( fvec->nar == 0 ){ KILL_floatvec(fvec); break; }
02224
02225 for( ii=0 ; ii < fvec->nar && ii < ncol ; ii++ )
02226 tsar[nrow*ncol+ii] = fvec->ar[ii] ;
02227 for( ; ii < ncol ; ii++ )
02228 tsar[nrow*ncol+ii] = filler ;
02229 KILL_floatvec(fvec) ;
02230 nrow++ ;
02231 }
02232 fclose( fts ) ;
02233 (void) my_fgets( NULL , 0 , NULL ) ;
02234
02235 mri_add_name( fname , outim ) ;
02236 FRB(buf) ; lbfill = 0.0f ; RETURN(outim) ;
02237 }
02238
02239
02240
02241
02242
02243 static void read_ascii_floats( char * fname, int * nff , float ** ff )
02244 {
02245 int ii,jj,val , used_tsar , alloc_tsar ;
02246 float *tsar ;
02247 float ftemp ;
02248 FILE *fts ;
02249 char *buf ;
02250 char *ptr ;
02251 int bpos , blen , nrow ;
02252
02253
02254
02255 if( nff == NULL || ff == NULL ) return ;
02256 if( fname == NULL || fname[0] == '\0' ){ *nff=0 ; *ff=NULL ; return ; }
02257
02258 fts = fopen( fname , "r" ) ;
02259 if( fts == NULL ){ *nff=0 ; *ff=NULL ; return ; }
02260
02261
02262
02263 used_tsar = 0 ;
02264 alloc_tsar = INC_TSARSIZE ;
02265 tsar = (float *) malloc( sizeof(float) * alloc_tsar ) ;
02266 if( tsar == NULL ){
02267 fprintf(stderr,"\n*** malloc fails: read_ascii_floats ***\n"); EXIT(1);
02268 }
02269
02270
02271
02272 nrow = 0 ;
02273 buf = (char *)malloc(LBUF) ;
02274 while( 1 ){
02275 ptr = fgets( buf , LBUF , fts ) ;
02276 if( ptr == NULL ) break ;
02277 blen = strlen(buf) ;
02278 if( blen <= 0 ) break ;
02279
02280 for( ii=0 ; ii < blen && isspace(buf[ii]) ; ii++ ) ;
02281
02282 if( ii == blen ) continue ;
02283 if( buf[ii] == '#' ) continue ;
02284 if( buf[ii] == '!' ) continue ;
02285
02286
02287
02288 for( jj=ii ; jj < blen ; jj++ ) if( buf[jj] == ',' ) buf[jj] = ' ' ;
02289
02290 for( bpos=ii ; bpos < blen ; ){
02291 val = sscanf( buf+bpos , "%f%n" , &ftemp , &jj ) ;
02292 if( val < 1 ) break ;
02293 bpos += jj ;
02294
02295 if( used_tsar == alloc_tsar ){
02296 alloc_tsar += INC_TSARSIZE ;
02297 tsar = (float *)realloc( tsar,sizeof(float)*alloc_tsar );
02298 if( tsar == NULL ){
02299 fprintf(stderr,"\n*** realloc fails: read_ascii_floats ***\n"); EXIT(1);
02300 }
02301 }
02302
02303 tsar[used_tsar++] = ftemp ;
02304 }
02305
02306 nrow++ ;
02307 }
02308 fclose( fts ) ;
02309 free( buf ) ;
02310
02311 if( used_tsar <= 1 ){ free(tsar); *nff=0; *ff=NULL; return; }
02312
02313 tsar = (float *) realloc( tsar , sizeof(float) * used_tsar ) ;
02314 if( tsar == NULL ){
02315 fprintf(stderr,"\n*** final realloc fails: read_ascii_floats ***\n"); EXIT(1);
02316 }
02317
02318 *nff = used_tsar; *ff = tsar; return;
02319 }
02320
02321
02322
02323
02324
02325
02326
02327 MRI_IMARR * mri_read_3A( char * tname )
02328 {
02329 int nx , ny , nz , ii , nxyz,nxy , nff ;
02330 int ngood , length , kim , datum_type ;
02331 char fname[256]="\0" , buf[512] ;
02332 MRI_IMARR * newar ;
02333 MRI_IMAGE * newim , * flim ;
02334 float * ff ;
02335
02336 ENTRY("mri_read_3A") ;
02337
02338
02339
02340 if( tname == NULL || strlen(tname) < 10 ) RETURN(NULL) ;
02341
02342 switch( tname[2] ){
02343
02344 default: ngood = 0 ; break ;
02345
02346 case 's':
02347 ngood = sscanf( tname, "3As:%d:%d:%d:%s", &nx, &ny, &nz, fname ) ;
02348 datum_type = MRI_short ;
02349 break ;
02350
02351 case 'b':
02352 ngood = sscanf( tname, "3Ab:%d:%d:%d:%s", &nx, &ny, &nz, fname ) ;
02353 datum_type = MRI_byte ;
02354 break ;
02355
02356 case 'f':
02357 ngood = sscanf( tname, "3Af:%d:%d:%d:%s", &nx, &ny, &nz, fname ) ;
02358 datum_type = MRI_float ;
02359 break ;
02360 }
02361
02362 if( ngood < 4 || nx <= 0 || ny <= 0 || nz <= 0 || strlen(fname) <= 0 ) RETURN(NULL) ;
02363
02364
02365
02366 read_ascii_floats( fname , &nff , &ff ) ;
02367
02368 if( nff <= 0 || ff == NULL ) RETURN(NULL) ;
02369
02370 nxy = nx*ny ; nxyz = nxy*nz ;
02371
02372 if( nff < nxyz ){
02373 fprintf(stderr,
02374 "\n** WARNING: %s is too short - padding with %d zeros\n",
02375 tname,nxyz-nff) ;
02376 ff = (float *) realloc( ff , sizeof(float) * nxyz ) ;
02377 for( ii=nff ; ii < nxyz ; ii++ ) ff[ii] = 0.0 ;
02378 nff = nxyz ;
02379 } else if( nff > nxyz ){
02380 fprintf(stderr,
02381 "\n** WARNING: %s is too long - truncating off last %d values\n",
02382 tname,nff-nxyz) ;
02383 }
02384
02385
02386
02387 INIT_IMARR(newar) ;
02388
02389 for( kim=0 ; kim < nz ; kim++ ){
02390 flim = mri_new( nx,ny , MRI_float ) ;
02391 memcpy( MRI_FLOAT_PTR(flim) , ff+nxy*kim , sizeof(float)*nxy ) ;
02392 switch( datum_type ){
02393 case MRI_float: newim = flim ; break ;
02394 case MRI_short: newim = mri_to_short(1.0,flim) ; mri_free(flim) ; break ;
02395 case MRI_byte: newim = mri_to_byte_scl(1.0,0.0,flim) ; mri_free(flim) ; break ;
02396 }
02397
02398 if( nz == 1 ) mri_add_name( fname , newim ) ;
02399 else {
02400 sprintf( buf , "%s#%d" , fname,kim ) ;
02401 mri_add_name( buf , newim ) ;
02402 }
02403
02404 ADDTO_IMARR(newar,newim) ;
02405 }
02406
02407 free(ff) ; RETURN(newar) ;
02408 }
02409
02410
02411
02412
02413
02414
02415
02416 #include "mayo_analyze.h"
02417
02418
02419
02420
02421 static void swap_analyze_hdr( struct dsr *pntr )
02422 {
02423 ENTRY("swap_analyze_hdr") ;
02424 swap_4(&pntr->hk.sizeof_hdr) ;
02425 swap_4(&pntr->hk.extents) ;
02426 swap_2(&pntr->hk.session_error) ;
02427 swap_2(&pntr->dime.dim[0]) ;
02428 swap_2(&pntr->dime.dim[1]) ;
02429 swap_2(&pntr->dime.dim[2]) ;
02430 swap_2(&pntr->dime.dim[3]) ;
02431 swap_2(&pntr->dime.dim[4]) ;
02432 swap_2(&pntr->dime.dim[5]) ;
02433 swap_2(&pntr->dime.dim[6]) ;
02434 swap_2(&pntr->dime.dim[7]) ;
02435 #if 0
02436 swap_2(&pntr->dime.unused1) ;
02437 #endif
02438 swap_2(&pntr->dime.datatype) ;
02439 swap_2(&pntr->dime.bitpix) ;
02440 swap_4(&pntr->dime.pixdim[0]) ;
02441 swap_4(&pntr->dime.pixdim[1]) ;
02442 swap_4(&pntr->dime.pixdim[2]) ;
02443 swap_4(&pntr->dime.pixdim[3]) ;
02444 swap_4(&pntr->dime.pixdim[4]) ;
02445 swap_4(&pntr->dime.pixdim[5]) ;
02446 swap_4(&pntr->dime.pixdim[6]) ;
02447 swap_4(&pntr->dime.pixdim[7]) ;
02448 swap_4(&pntr->dime.vox_offset) ;
02449 swap_4(&pntr->dime.funused1) ;
02450 swap_4(&pntr->dime.funused2) ;
02451 swap_4(&pntr->dime.cal_max) ;
02452 swap_4(&pntr->dime.cal_min) ;
02453 swap_4(&pntr->dime.compressed) ;
02454 swap_4(&pntr->dime.verified) ;
02455 swap_2(&pntr->dime.dim_un0) ;
02456 swap_4(&pntr->dime.glmax) ;
02457 swap_4(&pntr->dime.glmin) ;
02458 EXRETURN ;
02459 }
02460
02461
02462
02463
02464
02465
02466
02467 static int mri_imcount_analyze75( char * hname )
02468 {
02469 FILE * fp ;
02470 struct dsr hdr ;
02471 int doswap , nz ;
02472
02473 ENTRY("mri_imcount_analyze75") ;
02474
02475 fp = fopen( hname , "rb" ) ;
02476 if( fp == NULL ) RETURN(0) ;
02477 hdr.dime.dim[0] = 0 ;
02478 fread( &hdr , 1 , sizeof(struct dsr) , fp ) ;
02479 fclose(fp) ;
02480 if( hdr.dime.dim[0] == 0 ) RETURN(0) ;
02481 doswap = (hdr.dime.dim[0] < 0 || hdr.dime.dim[0] > 15) ;
02482 if( doswap ) swap_analyze_hdr( &hdr ) ;
02483
02484 switch( hdr.dime.dim[0] ){
02485 case 2: nz = 1 ; break ;
02486 case 3: nz = hdr.dime.dim[3] ; break ;
02487
02488 default:
02489 case 4: nz = hdr.dime.dim[3] * hdr.dime.dim[4] ; break ;
02490 }
02491 if( nz < 1 ) nz = 1 ;
02492
02493 RETURN(nz) ;
02494 }
02495
02496
02497
02498
02499
02500
02501
02502 MRI_IMARR * mri_read_analyze75( char * hname )
02503 {
02504 FILE * fp ;
02505 char iname[1024] , buf[1024] ;
02506 int ii , jj , doswap ;
02507 struct dsr hdr ;
02508 int ngood , length , kim , koff , datum_type , datum_len , swap ;
02509 int nx,ny,nz , hglobal=0 , himage=0 ;
02510 float dx,dy,dz ;
02511 MRI_IMARR * newar ;
02512 MRI_IMAGE * newim ;
02513 void * imar ;
02514 float fac=0.0 ;
02515 int floatize ;
02516 int spmorg=0 ;
02517
02518 ENTRY("mri_read_analyze75") ;
02519
02520
02521
02522 if( hname == NULL ) RETURN(NULL) ;
02523 jj = strlen(hname) ;
02524 if( jj < 5 ) RETURN(NULL) ;
02525 if( strcmp(hname+jj-3,"hdr") != 0 ) RETURN(NULL) ;
02526 strcpy(iname,hname) ; strcpy(iname+jj-3,"img") ;
02527
02528
02529
02530 fp = fopen( hname , "rb" ) ;
02531 if( fp == NULL ) RETURN(NULL) ;
02532 hdr.dime.dim[0] = 0 ;
02533 fread( &hdr , 1 , sizeof(struct dsr) , fp ) ;
02534 fclose(fp) ;
02535 if( hdr.dime.dim[0] == 0 ) RETURN(NULL) ;
02536
02537
02538
02539 doswap = (hdr.dime.dim[0] < 0 || hdr.dime.dim[0] > 15) ;
02540 if( doswap ) swap_analyze_hdr( &hdr ) ;
02541
02542
02543
02544 { short xyzuv[5] , xx,yy,zz ;
02545 memcpy( xyzuv , hdr.hist.originator , 10 ) ;
02546 if( xyzuv[3] == 0 && xyzuv[4] == 0 ){
02547 xx = xyzuv[0] ; yy = xyzuv[1] ; zz = xyzuv[2] ;
02548 if( doswap ){ swap_2(&xx); swap_2(&yy); swap_2(&zz); }
02549 if( xx > 0 && xx < hdr.dime.dim[1] &&
02550 yy > 0 && yy < hdr.dime.dim[2] &&
02551 zz > 0 && zz < hdr.dime.dim[3] ) spmorg = 1 ;
02552 }
02553 }
02554 if( spmorg ) strcpy( MRILIB_orients , "LRPAIS" ) ;
02555
02556
02557
02558 if( !AFNI_noenv("AFNI_ANALYZE_SCALE") ){
02559 fac = hdr.dime.funused1 ;
02560 (void) thd_floatscan( 1 , &fac ) ;
02561 if( fac < 0.0 || fac == 1.0 ) fac = 0.0 ;
02562 }
02563
02564 floatize = (fac != 0.0 || AFNI_yesenv("AFNI_ANALYZE_FLOATIZE")) ;
02565
02566
02567
02568 switch( hdr.dime.datatype ){
02569 default:
02570 fprintf(stderr,"*** %s: Unknown ANALYZE datatype=%d (%s)\n",
02571 hname,hdr.dime.datatype,ANDT_string(hdr.dime.datatype) ) ;
02572 RETURN(NULL) ;
02573
02574 case ANDT_UNSIGNED_CHAR: datum_type = MRI_byte ; break;
02575 case ANDT_SIGNED_SHORT: datum_type = MRI_short ; break;
02576 case ANDT_SIGNED_INT: datum_type = MRI_int ; break;
02577 case ANDT_FLOAT: datum_type = MRI_float ; floatize = 0; break;
02578 case ANDT_COMPLEX: datum_type = MRI_complex; floatize = 0; break;
02579 case ANDT_RGB: datum_type = MRI_rgb ; floatize = 0; break;
02580 case ANDT_DOUBLE: datum_type = MRI_double ; floatize = 0; break;
02581 }
02582
02583 datum_len = mri_datum_size(datum_type) ;
02584
02585
02586
02587 nx = hdr.dime.dim[1] ;
02588 ny = hdr.dime.dim[2] ;
02589 if( nx < 2 || ny < 2 ) RETURN(NULL) ;
02590
02591 switch( hdr.dime.dim[0] ){
02592 case 2: nz = 1 ; break ;
02593 case 3: nz = hdr.dime.dim[3] ; break ;
02594
02595 default:
02596 case 4: nz = hdr.dime.dim[3] * hdr.dime.dim[4] ; break ;
02597 }
02598 if( nz < 1 ) nz = 1 ;
02599
02600 dx = hdr.dime.pixdim[1] ;
02601 dy = hdr.dime.pixdim[2] ;
02602 dz = hdr.dime.pixdim[3] ;
02603
02604
02605
02606
02607
02608
02609 length = THD_filesize(iname) ;
02610 if( length <= 0 ){
02611 fprintf(stderr,"*** Can't find ANALYZE file %s\n",iname) ;
02612 RETURN(NULL) ;
02613 }
02614
02615 fp = fopen( iname , "rb" ) ;
02616 if( fp == NULL ){
02617 fprintf(stderr,"*** Can't open ANALYZE file %s\n",iname) ;
02618 RETURN(NULL) ;
02619 }
02620
02621 ngood = datum_len*nx*ny*nz ;
02622 if( length < ngood ){
02623 fprintf( stderr,
02624 "*** ANALYZE file %s is %d bytes long but must be at least %d bytes long\n"
02625 "*** for nx=%d ny=%d nz=%d and voxel=%d bytes\n",
02626 iname,length,ngood,nx,ny,nz,datum_len ) ;
02627 fclose(fp) ; RETURN(NULL) ;
02628 }
02629
02630
02631
02632 INIT_IMARR(newar) ;
02633
02634 for( kim=0 ; kim < nz ; kim++ ){
02635 koff = hglobal + (kim+1)*himage + datum_len*nx*ny*kim ;
02636
02637 fseek( fp , koff , SEEK_SET ) ;
02638
02639 newim = mri_new( nx , ny , datum_type ) ;
02640 imar = mri_data_pointer( newim ) ;
02641 length = fread( imar , datum_len , nx * ny , fp ) ;
02642
02643 if( doswap ){
02644 switch( datum_len ){
02645 default: break ;
02646 case 2: swap_twobytes ( nx*ny , imar ) ; break ;
02647 case 4: swap_fourbytes( nx*ny , imar ) ; break ;
02648 case 8: swap_fourbytes( 2*nx*ny , imar ) ; break ;
02649 }
02650 newim->was_swapped = 1 ;
02651 }
02652
02653
02654
02655 if( floatize ){
02656 MRI_IMAGE *qim = mri_to_float(newim) ;
02657 mri_free(newim) ; newim = qim ;
02658 }
02659
02660 if( nz == 1 ) mri_add_name( iname , newim ) ;
02661 else {
02662 sprintf( buf , "%s#%d" , iname,kim ) ;
02663 mri_add_name( buf , newim ) ;
02664 }
02665
02666 newim->dx = dx ; newim->dy = dy ; newim->dz = dz ; newim->dw = 1.0 ;
02667 ADDTO_IMARR(newar,newim) ;
02668
02669
02670
02671 if( fac != 0.0 ) mri_scale_inplace( fac , newim ) ;
02672 }
02673
02674 fclose(fp) ; RETURN(newar) ;
02675 }
02676
02677
02678
02679
02680
02681
02682
02683 MRI_IMARR * mri_read3D_analyze75( char * hname )
02684 {
02685 FILE * fp ;
02686 char iname[1024] , buf[1024] ;
02687 int ii , jj , doswap ;
02688 struct dsr hdr ;
02689 int ngood , length , kim , koff , datum_type , datum_len , swap ;
02690 int nx,ny,nz , hglobal=0 , himage=0 ;
02691 float dx,dy,dz ;
02692 MRI_IMARR * newar ;
02693 MRI_IMAGE * newim ;
02694 void * imar ;
02695 float fac=0.0 ;
02696 int floatize ;
02697 int spmorg=0 ;
02698
02699 int nt , nxyz ;
02700 float dt ;
02701
02702 ENTRY("mri_read3D_analyze75") ;
02703
02704
02705
02706 if( hname == NULL ) RETURN(NULL) ;
02707 jj = strlen(hname) ;
02708 if( jj < 5 ) RETURN(NULL) ;
02709 if( strcmp(hname+jj-3,"hdr") != 0 ) RETURN(NULL) ;
02710 strcpy(iname,hname) ; strcpy(iname+jj-3,"img") ;
02711
02712
02713
02714 fp = fopen( hname , "rb" ) ;
02715 if( fp == NULL ) RETURN(NULL) ;
02716 hdr.dime.dim[0] = 0 ;
02717 fread( &hdr , 1 , sizeof(struct dsr) , fp ) ;
02718 fclose(fp) ;
02719 if( hdr.dime.dim[0] == 0 ) RETURN(NULL) ;
02720
02721
02722
02723 doswap = (hdr.dime.dim[0] < 0 || hdr.dime.dim[0] > 15) ;
02724 if( doswap ) swap_analyze_hdr( &hdr ) ;
02725
02726
02727
02728 { short xyzuv[5] , xx,yy,zz ;
02729 memcpy( xyzuv , hdr.hist.originator , 10 ) ;
02730 if( xyzuv[3] == 0 && xyzuv[4] == 0 ){
02731 xx = xyzuv[0] ; yy = xyzuv[1] ; zz = xyzuv[2] ;
02732 if( doswap ){ swap_2(&xx); swap_2(&yy); swap_2(&zz); }
02733 if( xx > 0 && xx < hdr.dime.dim[1] &&
02734 yy > 0 && yy < hdr.dime.dim[2] &&
02735 zz > 0 && zz < hdr.dime.dim[3] ) spmorg = 1 ;
02736 }
02737 }
02738 if( spmorg ) strcpy( MRILIB_orients , "LRPAIS" ) ;
02739
02740
02741
02742 if( !AFNI_noenv("AFNI_ANALYZE_SCALE") ){
02743 fac = hdr.dime.funused1 ;
02744 (void) thd_floatscan( 1 , &fac ) ;
02745 if( fac < 0.0 || fac == 1.0 ) fac = 0.0 ;
02746 }
02747
02748 floatize = (fac != 0.0 || AFNI_yesenv("AFNI_ANALYZE_FLOATIZE")) ;
02749
02750
02751
02752 switch( hdr.dime.datatype ){
02753 default:
02754 fprintf(stderr,"*** %s: Unknown ANALYZE datatype=%d (%s)\n",
02755 hname,hdr.dime.datatype,ANDT_string(hdr.dime.datatype) ) ;
02756 RETURN(NULL) ;
02757
02758 case ANDT_UNSIGNED_CHAR: datum_type = MRI_byte ; break;
02759 case ANDT_SIGNED_SHORT: datum_type = MRI_short ; break;
02760 case ANDT_SIGNED_INT: datum_type = MRI_int ; break;
02761 case ANDT_FLOAT: datum_type = MRI_float ; floatize = 0; break;
02762 case ANDT_COMPLEX: datum_type = MRI_complex; floatize = 0; break;
02763 case ANDT_RGB: datum_type = MRI_rgb ; floatize = 0; break;
02764 }
02765
02766 datum_len = mri_datum_size(datum_type) ;
02767
02768
02769
02770 nx = hdr.dime.dim[1] ;
02771 ny = hdr.dime.dim[2] ;
02772 if( nx < 2 || ny < 2 ) RETURN(NULL) ;
02773
02774 switch( hdr.dime.dim[0] ){
02775 case 2: nz = 1 ; nt = 1 ; ; break ;
02776 case 3: nz = hdr.dime.dim[3] ; nt = 1 ; ; break ;
02777
02778 default:
02779 case 4: nz = hdr.dime.dim[3] ; nt = hdr.dime.dim[4] ; break ;
02780 }
02781 if( nz < 1 ) nz = 1 ;
02782 if( nt < 1 ) nt = 1 ;
02783
02784 dx = hdr.dime.pixdim[1] ;
02785 dy = hdr.dime.pixdim[2] ;
02786 dz = hdr.dime.pixdim[3] ;
02787 dt = hdr.dime.pixdim[4] ; if( dt <= 0.0 ) dt = 1.0 ;
02788
02789
02790
02791 length = THD_filesize(iname) ;
02792 if( length <= 0 ){
02793 fprintf(stderr,"*** Can't find ANALYZE file %s\n",iname) ;
02794 RETURN(NULL) ;
02795 }
02796
02797 fp = fopen( iname , "rb" ) ;
02798 if( fp == NULL ){
02799 fprintf(stderr,"*** Can't open ANALYZE file %s\n",iname) ;
02800 RETURN(NULL) ;
02801 }
02802
02803 ngood = datum_len*nx*ny*nz*nt ;
02804 if( length < ngood ){
02805 fprintf( stderr,
02806 "*** ANALYZE file %s is %d bytes long but must be at least %d bytes long\n"
02807 "*** for nx=%d ny=%d nz=%d nt=%d and voxel=%d bytes\n",
02808 iname,length,ngood,nx,ny,nz,nt,datum_len ) ;
02809 fclose(fp) ; RETURN(NULL) ;
02810 }
02811
02812
02813
02814 INIT_IMARR(newar) ;
02815
02816 for( kim=0 ; kim < nt ; kim++ ){
02817 koff = hglobal + (kim+1)*himage + datum_len*nxyz*kim ;
02818 fseek( fp , koff , SEEK_SET ) ;
02819
02820 newim = mri_new_vol( nx,ny,nz , datum_type ) ;
02821 imar = mri_data_pointer( newim ) ;
02822 length = fread( imar , datum_len , nxyz , fp ) ;
02823
02824 if( doswap ){
02825 switch( datum_len ){
02826 default: break ;
02827 case 2: swap_twobytes ( nxyz , imar ) ; break ;
02828 case 4: swap_fourbytes( nxyz , imar ) ; break ;
02829 case 8: swap_fourbytes( 2*nxyz , imar ) ; break ;
02830 }
02831 newim->was_swapped = 1 ;
02832 }
02833
02834
02835
02836 if( floatize ){
02837 MRI_IMAGE *qim = mri_to_float(newim) ;
02838 mri_free(newim) ; newim = qim ;
02839 }
02840
02841 if( nt == 1 ) mri_add_name( iname , newim ) ;
02842 else {
02843 sprintf( buf , "%s#%d" , iname,kim ) ;
02844 mri_add_name( buf , newim ) ;
02845 }
02846
02847 newim->dx = dx ; newim->dy = dy ; newim->dz = dz ; newim->dt = dt ; newim->dw = 1.0 ;
02848 ADDTO_IMARR(newar,newim) ;
02849
02850
02851
02852 if( fac != 0.0 ) mri_scale_inplace( fac , newim ) ;
02853 }
02854
02855 fclose(fp) ; RETURN(newar) ;
02856 }
02857
02858
02859
02860
02861
02862 #include "siemens_vision.h"
02863
02864
02865
02866
02867
02868
02869
02870
02871
02872 static int mri_imcount_siemens( char * hname )
02873 {
02874 struct Siemens_vision_header head ;
02875 FILE * fp ;
02876 int i,j,xx,yy , matrix , swap , imagesize,nxx,blank , slices ;
02877 struct stat file_stat ;
02878 short *imar ;
02879
02880
02881
02882 if( hname == NULL ) return 0 ;
02883
02884 i = stat( hname , &file_stat ) ;
02885 if( i < 0 ) return 0 ;
02886
02887
02888
02889 fp = fopen( hname , "rb" ) ;
02890 if( fp == NULL ) return 0 ;
02891 fread( &head , sizeof(struct Siemens_vision_header) , 1 , fp ) ;
02892
02893
02894
02895 swap = ( head.SiemensStudyDateMM < 0 || head.SiemensStudyDateMM > 13 ) ;
02896 if( swap ){
02897 swap_4( &(head.SiemensStudyDateMM) ) ;
02898 if( head.SiemensStudyDateMM < 0 || head.SiemensStudyDateMM > 13 ){
02899 swap = 0 ;
02900 }
02901 }
02902
02903
02904
02905 if( swap ) swap_4( &(head.DisplayMatrixSize) ) ;
02906 imagesize = head.DisplayMatrixSize ;
02907
02908
02909
02910 #undef MATRIX_MAX
02911 #define MATRIX_MAX 9
02912
02913 i = 2*imagesize*imagesize ;
02914 for( matrix=1 ; matrix < MATRIX_MAX ; matrix++ )
02915 if( file_stat.st_size == i*matrix*matrix + SIEMENS_HEADERSIZE ) break ;
02916
02917 if( matrix == MATRIX_MAX ){
02918 fclose(fp) ; return 0 ;
02919 }
02920 #undef MATRIX_MAX
02921
02922
02923
02924 imar = (short *) calloc(sizeof(short),matrix*matrix*imagesize*imagesize) ;
02925 fseek( fp , SIEMENS_HEADERSIZE , SEEK_SET ) ;
02926 fread( imar , sizeof(short) , matrix*matrix*imagesize*imagesize , fp ) ;
02927 fclose(fp) ;
02928
02929
02930
02931 slices = 0 ; nxx = matrix*imagesize ;
02932
02933 for( yy=0 ; yy < matrix ; yy++ ){
02934 for( xx=0 ; xx < matrix ; xx++ ){
02935 blank = 1 ;
02936 for( j=0 ; j < imagesize ; j++ ){
02937 for( i=0 ; i < imagesize ; i++ ){
02938 if( imar[i+xx*imagesize+(j+yy*imagesize)*nxx] ) blank = 0 ;
02939 }
02940 }
02941 if( !blank ) slices = 1 + xx + yy*matrix ;
02942 }
02943 }
02944
02945 free(imar) ; return slices ;
02946 }
02947
02948
02949
02950
02951
02952
02953
02954
02955 MRI_IMARR * mri_read_siemens( char * hname )
02956 {
02957 struct Siemens_vision_header head ;
02958 FILE * fp ;
02959 int i,j,xx,yy , matrix , swap , imagesize,nxx,blank , slices,nz ;
02960 struct stat file_stat ;
02961 short *imar ;
02962 MRI_IMARR * newar ;
02963 MRI_IMAGE * newim ;
02964 short * nar ;
02965 char buf[256] ;
02966 float dx,dy,dz ;
02967 char *eee ; int ileave=0 ;
02968
02969 ENTRY("mri_read_siemens") ;
02970
02971
02972
02973 if( hname == NULL ) RETURN(NULL) ;
02974
02975 i = stat( hname , &file_stat ) ;
02976 if( i < 0 ) RETURN(NULL) ;
02977
02978
02979
02980 fp = fopen( hname , "rb" ) ;
02981 if( fp == NULL ) RETURN(NULL) ;
02982 fread( &head , sizeof(struct Siemens_vision_header) , 1 , fp ) ;
02983
02984
02985
02986 swap = ( head.SiemensStudyDateMM < 0 || head.SiemensStudyDateMM > 13 ) ;
02987 if( swap ){
02988 swap_4( &(head.SiemensStudyDateMM) ) ;
02989 if( head.SiemensStudyDateMM < 0 || head.SiemensStudyDateMM > 13 ){
02990 swap = 0 ;
02991 }
02992 }
02993
02994
02995
02996 if( swap ) swap_4( &(head.DisplayMatrixSize) ) ;
02997 imagesize = head.DisplayMatrixSize ;
02998
02999
03000
03001 #undef MATRIX_MAX
03002 #define MATRIX_MAX 16
03003
03004 i = 2*imagesize*imagesize ;
03005 for( matrix=1 ; matrix < MATRIX_MAX ; matrix++ )
03006 if( file_stat.st_size == i*matrix*matrix + SIEMENS_HEADERSIZE ) break ;
03007
03008 if( matrix == MATRIX_MAX ){
03009 fclose(fp) ; RETURN(NULL) ;
03010 }
03011 #undef MATRIX_MAX
03012
03013
03014
03015 imar = (short *) calloc(sizeof(short),matrix*matrix*imagesize*imagesize) ;
03016 fseek( fp , SIEMENS_HEADERSIZE , SEEK_SET ) ;
03017 fread( imar , sizeof(short) , matrix*matrix*imagesize*imagesize , fp ) ;
03018 fclose(fp) ;
03019
03020 if( swap ) swap_twobytes( matrix*matrix*imagesize*imagesize , imar ) ;
03021
03022
03023
03024 slices = 0 ; nxx = matrix*imagesize ;
03025
03026 for( yy=0 ; yy < matrix ; yy++ ){
03027 for( xx=0 ; xx < matrix ; xx++ ){
03028 blank = 1 ;
03029 for( j=0 ; j < imagesize ; j++ ){
03030 for( i=0 ; i < imagesize ; i++ ){
03031 if( imar[i+xx*imagesize+(j+yy*imagesize)*nxx] ) blank = 0 ;
03032 }
03033 }
03034 if( !blank ) slices = 1 + xx + yy*matrix ;
03035 }
03036 }
03037
03038 if( slices == 0 ){ free(imar) ; RETURN(NULL) ; }
03039
03040
03041
03042 if( swap ){
03043 swap_8(&(head.FOVRow));
03044 swap_8(&(head.FOVColumn));
03045 swap_8(&(head.SliceThickness));
03046 }
03047 dx = head.FOVRow / imagesize ;
03048 dy = head.FOVColumn / imagesize ;
03049 dz = head.SliceThickness ;
03050
03051
03052
03053 MRILIB_orients[0] = head.OrientationSet1Left[0] ;
03054 MRILIB_orients[1] = head.OrientationSet2Right[0];
03055 MRILIB_orients[2] = head.OrientationSet1Top[0] ;
03056 MRILIB_orients[3] = head.OrientationSet2Down[0] ;
03057 MRILIB_orients[4] = head.OrientationSet1Back[0] ;
03058 MRILIB_orients[5] = head.OrientationSet2Front[0];
03059 for (i=0; i<6; i++) {
03060 if (MRILIB_orients[i]=='H') MRILIB_orients[i]='S';
03061 if (MRILIB_orients[i]=='F') MRILIB_orients[i]='I';
03062 }
03063 MRILIB_orients[6] = '\0' ;
03064 MRILIB_zoff = fabs(strtod(head.TextSlicePosition,NULL)) ; use_MRILIB_zoff = 1 ;
03065
03066
03067
03068 INIT_IMARR(newar) ;
03069
03070 for( yy=0 ; yy < matrix ; yy++ ){
03071 for( xx=0 ; xx < matrix ; xx++ ){
03072
03073 newim = mri_new( imagesize , imagesize , MRI_short ) ;
03074 nar = MRI_SHORT_PTR( newim ) ;
03075
03076 if( swap ) newim->was_swapped = 1 ;
03077
03078 for( j=0 ; j < imagesize ; j++ )
03079 memcpy( nar+j*imagesize ,
03080 imar+xx*imagesize+(j+yy*imagesize)*nxx , 2*imagesize ) ;
03081
03082 sprintf( buf , "%s#%d:%d" , hname,xx,yy ) ;
03083 mri_add_name( buf , newim ) ;
03084
03085 newim->dx = dx ; newim->dy = dy ; newim->dz = dz ; newim->dw = 1.0 ;
03086 ADDTO_IMARR(newar,newim) ;
03087 if( IMARR_COUNT(newar) == slices ) goto Done ;
03088 }
03089 }
03090
03091 Done:
03092
03093
03094
03095 eee = getenv("AFNI_SIEMENS_INTERLEAVE") ;
03096 ileave = ( (eee != NULL) && (*eee=='Y' || *eee=='y') ) ;
03097 if( ileave && slices > 2 ){
03098 int mid = (slices-1)/2 ;
03099 MRI_IMARR *qar ;
03100 INIT_IMARR(qar) ;
03101 for( i=0 ; i < slices ; i++ ){
03102 if( i%2 == 0 ) j = i/2 ;
03103 else j = mid + (i+1)/2 ;
03104 ADDTO_IMARR(qar,IMARR_SUBIM(newar,j)) ;
03105 }
03106 FREE_IMARR(newar) ; newar = qar ;
03107 }
03108
03109 free(imar) ; RETURN(newar) ;
03110 }
03111
03112
03113
03114
03115
03116
03117 short check_dicom_magic_num( char * fname )
03118 {
03119 FILE * fp;
03120 char test_string[5] ;
03121
03122 fp = fopen( fname, "rb" ) ;
03123 if(fp == NULL ) return 0 ;
03124 fseek( fp, 128 , SEEK_SET ) ;
03125 fread( test_string , 1 , 4 , fp ) ; test_string[4] = '\0' ;
03126 fclose( fp ) ;
03127 if( strcmp(test_string,"DICM") == 0 ) {
03128 return 1 ;
03129 } else {
03130 return 0 ;
03131 }
03132 }
03133
03134
03135
03136
03137
03138 #ifdef USE_MRI_DELAY
03139
03140
03141
03142 void mri_purge_delay( MRI_IMAGE * im )
03143 {
03144 void * ar ;
03145
03146
03147
03148
03149 if( im->fname == NULL ||
03150 (im->fondisk & INPUT_DELAY) != 0 ) return ;
03151
03152
03153
03154 ar = mri_data_pointer( im ) ;
03155 if( ar != NULL ){ free(ar) ; mri_clear_data_pointer(im) ; }
03156
03157
03158
03159 im->fondisk |= INPUT_DELAY ;
03160 return ;
03161 }
03162
03163
03164
03165 void mri_input_delay( MRI_IMAGE * im )
03166 {
03167 FILE * imfile=NULL ;
03168 void * imar ;
03169
03170
03171
03172
03173 if( im->fname == NULL ||
03174 (im->fondisk & INPUT_DELAY) == 0 ) return ;
03175
03176
03177
03178 if( strcmp(im->fname,"ALLZERO") != 0 ){
03179 imfile = fopen( im->fname , "r" ) ;
03180 if( imfile == NULL ){
03181 fprintf( stderr , "couldn't open delayed image file %s\n" , im->fname ) ;
03182 return ;
03183 }
03184 }
03185
03186
03187
03188 imar = (void *) malloc( im->nvox * im->pixel_size ) ;
03189 if( imar == NULL ){
03190 fprintf( stderr ,
03191 "malloc fails for delayed image from file %s\n" , im->fname ) ;
03192 if( imfile != NULL ) fclose( imfile ) ;
03193 return ;
03194 }
03195 mri_fix_data_pointer( imar , im ) ;
03196
03197
03198
03199 if( imfile != NULL ){
03200 fseek( imfile , im->foffset , SEEK_SET ) ;
03201 fread( imar , im->pixel_size , im->nvox , imfile ) ;
03202 fclose( imfile ) ;
03203 } else {
03204 memset( imar , 0 , im->nvox * im->pixel_size ) ;
03205 }
03206
03207
03208
03209 if( (im->fondisk & BSWAP_DELAY) ){
03210 mri_swapbytes( im ) ;
03211 im->was_swapped = 1 ;
03212 }
03213
03214
03215
03216 im->fondisk ^= INPUT_DELAY ;
03217
03218 #if 0
03219 fprintf(stderr,"delayed input from file %s at offset %d\n",im->fname,im->foffset);
03220 #endif
03221 return ;
03222 }
03223
03224
03225
03226
03227
03228 MRI_IMARR * mri_read_file_delay( char * fname )
03229 {
03230 MRI_IMARR * newar=NULL ;
03231 MRI_IMAGE * newim ;
03232 char * new_fname ;
03233
03234 new_fname = imsized_fname( fname ) ;
03235 if( new_fname == NULL ) return NULL ;
03236
03237 if( strlen(new_fname) > 9 && new_fname[0] == '3' && new_fname[1] == 'D' &&
03238 (new_fname[2] == ':' || new_fname[3] == ':') ){
03239
03240
03241 newar = mri_read_3D_delay( new_fname ) ;
03242
03243 } else if( strlen(new_fname) > 9 &&
03244 new_fname[0] == '3' && new_fname[1] == 'A' && new_fname[3] == ':' ){
03245
03246 newar = mri_read_3A( new_fname ) ;
03247
03248 } else if( strstr(new_fname,".hdr") != NULL ||
03249 strstr(new_fname,".HDR") != NULL ){
03250
03251 newar = mri_read_analyze75( new_fname ) ;
03252
03253 } else if( strstr(new_fname,".ima") != NULL ||
03254 strstr(new_fname,".IMA") != NULL ){
03255
03256 newar = mri_read_siemens( new_fname ) ;
03257
03258 } else if( strstr(new_fname,".mpg" ) != NULL ||
03259 strstr(new_fname,".MPG" ) != NULL ||
03260 strstr(new_fname,".mpeg") != NULL ||
03261 strstr(new_fname,".MPEG") != NULL ){
03262
03263 newar = mri_read_mpeg( new_fname ) ;
03264 }
03265
03266
03267
03268
03269 if ((newar == NULL) && !AFNI_yesenv("AFNI_TRY_DICOM_LAST")) {
03270 newar = mri_read_dicom( new_fname ) ;
03271 }
03272
03273
03274
03275 if( newar == NULL ){
03276 newim = mri_read( new_fname ) ;
03277 if( newim == NULL ){ free(new_fname) ; return NULL ; }
03278 INIT_IMARR(newar) ;
03279 ADDTO_IMARR(newar,newim) ;
03280 }
03281
03282 if ( (newar == NULL) && AFNI_yesenv("AFNI_TRY_DICOM_LAST")) {
03283 newar = mri_read_dicom( new_fname ) ;
03284 }
03285
03286 free(new_fname) ;
03287 return newar ;
03288 }
03289
03290
03291
03292 MRI_IMARR * mri_read_3D_delay( char * tname )
03293 {
03294 int hglobal , himage , nx , ny , nz ;
03295 char fname[256] , buf[512] ;
03296 int ngood , length , kim , datum_type , datum_len , swap ;
03297 MRI_IMARR * newar ;
03298 MRI_IMAGE * newim ;
03299 FILE * imfile ;
03300
03301
03302
03303 if( tname == NULL || strlen(tname) < 10 ) return NULL ;
03304
03305 switch( tname[2] ){
03306
03307 default:
03308 case ':':
03309 ngood = sscanf( tname , "3D:%d:%d:%d:%d:%d:%s" ,
03310 &hglobal , &himage , &nx , &ny , &nz , fname ) ;
03311
03312 swap = 0 ;
03313 datum_type = MRI_short ;
03314 datum_len = sizeof(short) ;
03315 break ;
03316
03317 case 's':
03318 ngood = sscanf( tname , "3Ds:%d:%d:%d:%d:%d:%s" ,
03319 &hglobal , &himage , &nx , &ny , &nz , fname ) ;
03320
03321 swap = 1 ;
03322 datum_type = MRI_short ;
03323 datum_len = sizeof(short) ;
03324 break ;
03325
03326 case 'b':
03327 ngood = sscanf( tname , "3Db:%d:%d:%d:%d:%d:%s" ,
03328 &hglobal , &himage , &nx , &ny , &nz , fname ) ;
03329
03330 swap = 0 ;
03331 datum_type = MRI_byte ;
03332 datum_len = sizeof(byte) ;
03333 break ;
03334
03335 case 'f':
03336 ngood = sscanf( tname , "3Df:%d:%d:%d:%d:%d:%s" ,
03337 &hglobal , &himage , &nx , &ny , &nz , fname ) ;
03338
03339 swap = 0 ;
03340 datum_type = MRI_float ;
03341 datum_len = sizeof(float) ;
03342 break ;
03343
03344 case 'd':
03345 ngood = sscanf( tname , "3Dd:%d:%d:%d:%d:%d:%s" ,
03346 &hglobal , &himage , &nx , &ny , &nz , fname ) ;
03347
03348 swap = 0 ;
03349 datum_type = MRI_float ;
03350 datum_len = sizeof(double) ;
03351 break ;
03352
03353 case 'i':
03354 ngood = sscanf( tname , "3Di:%d:%d:%d:%d:%d:%s" ,
03355 &hglobal , &himage , &nx , &ny , &nz , fname ) ;
03356
03357 swap = 0 ;
03358 datum_type = MRI_int ;
03359 datum_len = sizeof(int) ;
03360 break ;
03361
03362 case 'c':
03363 ngood = sscanf( tname , "3Dc:%d:%d:%d:%d:%d:%s" ,
03364 &hglobal , &himage , &nx , &ny , &nz , fname ) ;
03365
03366 swap = 0 ;
03367 datum_type = MRI_complex ;
03368 datum_len = sizeof(complex) ;
03369 break ;
03370
03371 case 'r':
03372 ngood = sscanf( tname , "3Dr:%d:%d:%d:%d:%d:%s" ,
03373 &hglobal , &himage , &nx , &ny , &nz , fname ) ;
03374
03375 swap = 0 ;
03376 datum_type = MRI_rgb ;
03377 datum_len = 3*sizeof(byte) ;
03378 break ;
03379 }
03380
03381 if( ngood < 6 || himage < 0 ||
03382 nx <= 0 || ny <= 0 || nz <= 0 ||
03383 strlen(fname) <= 0 ) return NULL ;
03384
03385
03386
03387 if( strcmp(fname,"ALLZERO") != 0 ){
03388 imfile = fopen( fname , "r" ) ;
03389 if( imfile == NULL ){
03390 fprintf( stderr , "couldn't open delayed image file %s\n" , fname ) ;
03391 return NULL ;
03392 }
03393 } else {
03394 imfile = NULL ;
03395 }
03396
03397 if( imfile != NULL ){
03398 fseek( imfile , 0L , SEEK_END ) ;
03399 length = ftell( imfile ) ;
03400
03401
03402
03403
03404 #if 0
03405 if( hglobal < 0 ){
03406 hglobal = length - nz*(datum_len*nx*ny+himage) ;
03407 if( hglobal < 0 ) hglobal = 0 ;
03408 }
03409 #else
03410 if( hglobal == -1 || hglobal+himage < 0 ){
03411 hglobal = length - nz*(datum_len*nx*ny+himage) ;
03412 if( hglobal < 0 ) hglobal = 0 ;
03413 }
03414 #endif
03415
03416 ngood = hglobal + nz*(datum_len*nx*ny+himage) ;
03417 if( length < ngood ){
03418 fprintf( stderr,
03419 "file %s is %d bytes long but must be at least %d bytes long\n"
03420 "for hglobal=%d himage=%d nx=%d ny=%d nz=%d and voxel=%d bytes\n",
03421 fname,length,ngood,hglobal,himage,nx,ny,nz,datum_len ) ;
03422 fclose( imfile ) ;
03423 return NULL ;
03424 }
03425 fclose( imfile ) ;
03426 }
03427
03428
03429
03430 INIT_IMARR(newar) ;
03431
03432 for( kim=0 ; kim < nz ; kim++ ){
03433 newim = mri_new_vol_empty( nx,ny,1 , datum_type ) ;
03434 mri_add_fname_delay( fname , newim ) ;
03435 newim->fondisk = (swap) ? (INPUT_DELAY | BSWAP_DELAY)
03436 : (INPUT_DELAY) ;
03437 newim->foffset = hglobal + (kim+1)*himage + datum_len*nx*ny*kim ;
03438
03439 if( nz == 1 ) mri_add_name( fname , newim ) ;
03440 else {
03441 sprintf( buf , "%s#%d" , fname,kim ) ;
03442 mri_add_name( buf , newim ) ;
03443 }
03444
03445 ADDTO_IMARR(newar,newim) ;
03446 }
03447
03448 return newar ;
03449 }
03450 #endif