00001 #include "mrilib.h"
00002 #include "netcdf.h"
00003 
00004 /*******************************************************************/
00005 /******* 29 Oct 2001: Read a 3D MINC file as an AFNI dataset *******/
00006 /*******************************************************************/
00007 
00008 static int first_err = 1 ;
00009 static char *fname_err ;
00010 
00011 #define EPR(x,s1,s2)                                                       \
00012  do{ if( first_err ){ fprintf(stderr,"\n"); first_err=0; }                 \
00013      fprintf(stderr,"** NetCDF error [%s,%s]: %s\n",s1,s2,nc_strerror(x)); \
00014  } while(0)
00015 
00016 /*-----------------------------------------------------------------
00017   Read the attributes relevant to a MINC dimension variable
00018 -------------------------------------------------------------------*/
00019 
00020 typedef struct {
00021    int dimid , varid , good , afni_orient ;
00022    size_t len ;
00023    float start , step , xcos,ycos,zcos ;
00024    char spacetype[32] ;
00025 } mincdim ;
00026 
00027 static mincdim read_mincdim( int ncid , char *dname )
00028 {
00029    mincdim ddd ;
00030    int code ;
00031    float ccc[3] ;
00032    size_t lll ;
00033    static char *ename=NULL ;
00034 
00035    ddd.good = 0 ;  /* flag for bad result */
00036 
00037    ddd.step  = 1.0 ;                       /* defaults */
00038    ddd.xcos  = ddd.ycos = ddd.zcos = 0.0 ;
00039    ddd.spacetype[0] = '\0' ;
00040 
00041    lll = strlen(fname_err) + strlen(dname) + 4 ;
00042    ename = AFREALL(ename, char, lll) ;
00043    sprintf(ename,"%s:%s",fname_err,dname) ;
00044 
00045    /* get ID of this dimension name */
00046 
00047    code = nc_inq_dimid( ncid , dname , &ddd.dimid ) ;
00048    if( code != NC_NOERR ){ EPR(code,ename,"dimid"); return ddd; }
00049 
00050    /* get length of this dimension */
00051 
00052    code = nc_inq_dimlen( ncid , ddd.dimid , &ddd.len ) ;
00053    if( code != NC_NOERR ){ EPR(code,ename,"dimlen"); return ddd; }
00054 
00055    /* get ID of corresponding variable */
00056 
00057    code = nc_inq_varid( ncid , dname , &ddd.varid ) ;
00058    if( code != NC_NOERR ){   /* this is bad: try to fake it */
00059       if( first_err ){ fprintf(stderr,"\n"); first_err=0; }
00060       fprintf(stderr,"** MINC warning: %s variable missing\n",ename);
00061       ddd.start = -0.5*ddd.step*ddd.len ;
00062       ddd.good = 1 ; return ddd ;
00063    }
00064 
00065    /* get step attribute of this variable */
00066 
00067    code = nc_get_att_float( ncid , ddd.varid , "step" , &ddd.step ) ;
00068    if( code != NC_NOERR ){
00069       if( first_err ){ fprintf(stderr,"\n"); first_err=0; }
00070       fprintf(stderr,"** MINC warning: %s:step missing\n",ename);
00071    } else if( ddd.step == 0.0 ){
00072       if( first_err ){ fprintf(stderr,"\n"); first_err=0; }
00073       fprintf(stderr,"** MINC warning: %s:step=0\n",ename);
00074       ddd.step = 1.0 ;
00075    }
00076 
00077    /* get start attribute of this variable */
00078 
00079    code = nc_get_att_float( ncid , ddd.varid , "start" , &ddd.start ) ;
00080    if( code != NC_NOERR ){
00081       if( first_err ){ fprintf(stderr,"\n"); first_err=0; }
00082       fprintf(stderr,"** MINC warning: %s:start missing\n",ename) ;
00083       ddd.start = -0.5*ddd.step*ddd.len ;
00084    }
00085 
00086    /* get direction_cosines attribute of this variable [not used yet] */
00087 
00088    code = nc_get_att_float( ncid,ddd.varid , "direction_cosines" , ccc ) ;
00089    if( code == NC_NOERR ){
00090       ddd.xcos = ccc[0] ; ddd.ycos = ccc[1] ; ddd.zcos = ccc[2] ;
00091    }
00092 
00093    /* get spacetype attribute of this variable [Talairach or not?] */
00094 
00095    lll = 0 ;
00096    code = nc_inq_attlen( ncid , ddd.varid , "spacetype" , &lll ) ;
00097    if( code == NC_NOERR && lll > 0 && lll < 32 ){
00098       code = nc_get_att_text( ncid,ddd.varid , "spacetype" , ddd.spacetype ) ;
00099       if( code == NC_NOERR ){
00100          ddd.spacetype[lll] = '\0' ;  /* make sure ends in NUL */
00101       } else {
00102          ddd.spacetype[0] = '\0' ;    /* make sure is empty */
00103       }
00104    }
00105 
00106    ddd.good = 1 ; return ddd ;
00107 }
00108 
00109 /*-----------------------------------------------------------------
00110   Open a MINC file as an AFNI dataset
00111 -------------------------------------------------------------------*/
00112 
00113 THD_3dim_dataset * THD_open_minc( char *pathname )
00114 {
00115    THD_3dim_dataset *dset=NULL ;
00116    int ncid , code ;
00117    mincdim xspace , yspace , zspace , *xyz[3] ;
00118    int im_varid , im_rank , im_dimid[4] ;
00119    nc_type im_type ;
00120 
00121    char prefix[1024] , *ppp , tname[12] ;
00122    THD_ivec3 nxyz , orixyz ;
00123    THD_fvec3 dxyz , orgxyz ;
00124    int datum , nvals=1 , ibr , iview ;
00125    size_t len ;
00126    float dtime=1.0 ;
00127 
00128 ENTRY("THD_open_minc") ;
00129 
00130    /*-- open input file --*/
00131 
00132    if( pathname == NULL || pathname[0] == '\0' ) RETURN(NULL);
00133 
00134    first_err = 1 ;        /* put a newline before 1st error message */
00135    fname_err = pathname ; /* filename for errors in read_mincdim()  */
00136 
00137    code = nc_open( pathname , NC_NOWRITE , &ncid ) ;
00138    if( code != NC_NOERR ){ EPR(code,pathname,"open"); RETURN(NULL); }
00139 
00140    /*---------- get info about spatial axes ---------*/
00141    /*[[ MINC x and y are reversed relative to AFNI ]]*/
00142    /*[[ MINC x=L-R y=P-A;  AFNI x=R-L y=A-P; z=I-S ]]*/
00143 
00144    xspace = read_mincdim( ncid , "xspace" ) ;
00145    xspace.step = -xspace.step ; xspace.start = -xspace.start ;
00146    xspace.afni_orient = (xspace.step < 0.0) ? ORI_L2R_TYPE : ORI_R2L_TYPE ;
00147 
00148    yspace = read_mincdim( ncid , "yspace" ) ;
00149    yspace.step = -yspace.step ; yspace.start = -yspace.start ;
00150    yspace.afni_orient = (yspace.step < 0.0) ? ORI_P2A_TYPE : ORI_A2P_TYPE ;
00151 
00152    zspace = read_mincdim( ncid , "zspace" ) ;
00153    zspace.afni_orient = (zspace.step > 0.0) ? ORI_I2S_TYPE : ORI_S2I_TYPE ;
00154 
00155    /*-- if don't have all 3 named dimensions, quit --*/
00156 
00157    if( !xspace.good || !yspace.good || !zspace.good ){
00158       nc_close(ncid) ; RETURN(NULL) ;
00159    }
00160 
00161    /*-- get info about the image data --*/
00162 
00163    code = nc_inq_varid( ncid , "image" , &im_varid ) ;
00164    if( code != NC_NOERR ){ EPR(code,pathname,"image varid"); nc_close(ncid); RETURN(NULL); }
00165 
00166    if( !AFNI_yesenv("AFNI_MINC_FLOATIZE") ){  /* determine type of dataset values */
00167       code = nc_inq_vartype( ncid , im_varid , &im_type ) ;
00168       if( code != NC_NOERR ){ EPR(code,pathname,"image vartype"); nc_close(ncid); RETURN(NULL); }
00169       switch( im_type ){
00170          default:
00171             if( first_err ){ fprintf(stderr,"\n"); first_err=0; }
00172             fprintf(stderr,"** MINC error: can't handle image type code %d\n",im_type);
00173             nc_close(ncid); RETURN(NULL);
00174          break ;
00175 
00176          case NC_BYTE:  datum = MRI_byte  ; break ;
00177 
00178          case NC_SHORT: datum = MRI_short ; break ;
00179 
00180          case NC_FLOAT:
00181          case NC_DOUBLE:
00182          case NC_INT:   datum = MRI_float ; break ;
00183       }
00184    } else {                               /* always read in as floats if */
00185       datum = MRI_float ;                 /* AFNI_MINC_FLOATIZE is Yes   */
00186    }
00187 
00188    /* we allow nD mages only for n=3 or n=4 */
00189 
00190    code = nc_inq_varndims( ncid , im_varid , &im_rank ) ;
00191    if( code != NC_NOERR ){ EPR(code,pathname,"image varndims"); nc_close(ncid); RETURN(NULL); }
00192    if( im_rank < 3 || im_rank > 4 ){
00193       if( first_err ){ fprintf(stderr,"\n"); first_err=0; }
00194       fprintf(stderr,"** MINC error: image rank=%d\n",im_rank) ;
00195       nc_close(ncid); RETURN(NULL);
00196    }
00197 
00198    code = nc_inq_vardimid( ncid , im_varid , im_dimid ) ;
00199    if( code != NC_NOERR ){ EPR(code,pathname,"image vardimid"); nc_close(ncid); RETURN(NULL); }
00200 
00201    /*-- find axes order --*/
00202 
00203      /* fastest variation */
00204 
00205         if( im_dimid[im_rank-1] == xspace.dimid ) xyz[0] = &xspace ;
00206    else if( im_dimid[im_rank-1] == yspace.dimid ) xyz[0] = &yspace ;
00207    else if( im_dimid[im_rank-1] == zspace.dimid ) xyz[0] = &zspace ;
00208    else {
00209       if( first_err ){ fprintf(stderr,"\n"); first_err=0; }
00210       fprintf(stderr,"** MINC error: image dim[2] illegal\n") ;
00211       nc_close(ncid) ; RETURN(NULL) ;
00212    }
00213 
00214      /* next fastest variation */
00215 
00216         if( im_dimid[im_rank-2] == xspace.dimid ) xyz[1] = &xspace ;
00217    else if( im_dimid[im_rank-2] == yspace.dimid ) xyz[1] = &yspace ;
00218    else if( im_dimid[im_rank-2] == zspace.dimid ) xyz[1] = &zspace ;
00219    else {
00220       if( first_err ){ fprintf(stderr,"\n"); first_err=0; }
00221       fprintf(stderr,"** MINC error: image dim[1] illegal\n") ;
00222       nc_close(ncid) ; RETURN(NULL) ;
00223    }
00224 
00225      /* slowest variation */
00226 
00227         if( im_dimid[im_rank-3] == xspace.dimid ) xyz[2] = &xspace ;
00228    else if( im_dimid[im_rank-3] == yspace.dimid ) xyz[2] = &yspace ;
00229    else if( im_dimid[im_rank-3] == zspace.dimid ) xyz[2] = &zspace ;
00230    else {
00231       if( first_err ){ fprintf(stderr,"\n"); first_err=0; }
00232       fprintf(stderr,"** MINC error: image dim[0] illegal\n") ;
00233       nc_close(ncid) ; RETURN(NULL) ;
00234    }
00235 
00236    /*-- determine the length of the 4th dimension, if any --*/
00237 
00238    strcpy(tname,"MINC") ;  /* default name for 4th dimension */
00239 
00240    if( im_rank == 4 ){
00241       char fname[NC_MAX_NAME+4] ;
00242       code = nc_inq_dimlen( ncid , im_dimid[0] , &len ) ;
00243       if( code != NC_NOERR ){ EPR(code,pathname,"dimid[0] dimlen"); nc_close(ncid); RETURN(NULL); }
00244 
00245       nvals = len ;  /* number of volumes in this file */
00246 
00247       /* get the name of this dimension */
00248 
00249       code = nc_inq_dimname( ncid , im_dimid[0] , fname ) ;
00250       if( code == NC_NOERR ){
00251          MCW_strncpy(tname,fname,10) ;
00252       }
00253 
00254       /* get the stepsize of this dimension */
00255 
00256       code = nc_get_att_float( ncid , im_dimid[0] , "step" , &dtime ) ;
00257       if( code != NC_NOERR || dtime <= 0.0 ) dtime = 1.0 ;
00258    }
00259 
00260    /*-- make a dataset --*/
00261 
00262    dset = EDIT_empty_copy(NULL) ;
00263 
00264    ppp  = THD_trailname(pathname,0) ;               /* strip directory */
00265    MCW_strncpy( prefix , ppp , THD_MAX_PREFIX ) ;   /* to make prefix */
00266 
00267    nxyz.ijk[0] = xyz[0]->len ; dxyz.xyz[0] = xyz[0]->step ;  /* setup axes */
00268    nxyz.ijk[1] = xyz[1]->len ; dxyz.xyz[1] = xyz[1]->step ;
00269    nxyz.ijk[2] = xyz[2]->len ; dxyz.xyz[2] = xyz[2]->step ;
00270 
00271    orixyz.ijk[0] = xyz[0]->afni_orient ; orgxyz.xyz[0] = xyz[0]->start ;
00272    orixyz.ijk[1] = xyz[1]->afni_orient ; orgxyz.xyz[1] = xyz[1]->start ;
00273    orixyz.ijk[2] = xyz[2]->afni_orient ; orgxyz.xyz[2] = xyz[2]->start ;
00274 
00275    if( strstr(xyz[0]->spacetype,"alairach") != NULL ){  /* coord system */
00276       iview = VIEW_TALAIRACH_TYPE ;
00277    } else {
00278       iview = VIEW_ORIGINAL_TYPE ;
00279    }
00280 
00281    dset->idcode.str[0] = 'M' ;  /* overwrite 1st 3 bytes with something special */
00282    dset->idcode.str[1] = 'N' ;
00283    dset->idcode.str[2] = 'C' ;
00284 
00285    MCW_hash_idcode( pathname , dset ) ;  /* 06 May 2005 */
00286 
00287    EDIT_dset_items( dset ,
00288                       ADN_prefix      , prefix ,
00289                       ADN_datum_all   , datum ,
00290                       ADN_nxyz        , nxyz ,
00291                       ADN_xyzdel      , dxyz ,
00292                       ADN_xyzorg      , orgxyz ,
00293                       ADN_xyzorient   , orixyz ,
00294                       ADN_malloc_type , DATABLOCK_MEM_MALLOC ,
00295                       ADN_nvals       , nvals ,
00296                       ADN_type        , HEAD_ANAT_TYPE ,
00297                       ADN_view_type   , iview ,
00298                       ADN_func_type   , (nvals==1) ? ANAT_MRAN_TYPE
00299                                                    : ANAT_BUCK_TYPE ,
00300                     ADN_none ) ;
00301 
00302    if( nvals > 9 )              /* pretend it is 3D+time */
00303       EDIT_dset_items( dset ,
00304                          ADN_func_type, ANAT_EPI_TYPE ,
00305                          ADN_ntt      , nvals ,
00306                          ADN_ttorg    , 0.0 ,
00307                          ADN_ttdel    , dtime ,
00308                          ADN_ttdur    , 0.0 ,
00309                          ADN_tunits   , UNITS_SEC_TYPE ,
00310                        ADN_none ) ;
00311 
00312    /*-- flag to read data from disk using MINC/netCDF functions --*/
00313 
00314    dset->dblk->diskptr->storage_mode = STORAGE_BY_MINC ;
00315    strcpy( dset->dblk->diskptr->brick_name , pathname ) ;
00316 
00317    for( ibr=0 ; ibr < nvals ; ibr++ ){     /* make sub-brick labels */
00318       sprintf(prefix,"%s[%d]",tname,ibr) ;
00319       EDIT_BRICK_LABEL( dset , ibr , prefix ) ;
00320    }
00321 
00322    /*-- read and save the global:history attribute --*/
00323 
00324 #if 0
00325    sprintf(prefix,"THD_open_minc(%s)",pathname) ;
00326    tross_Append_History( dset , prefix ) ;
00327 #endif
00328 
00329    code = nc_inq_attlen( ncid , NC_GLOBAL , "history" , &len ) ;
00330    if( code == NC_NOERR && len > 0 ){
00331       ppp = AFMALL(char,len+4) ;
00332       code = nc_get_att_text( ncid , NC_GLOBAL , "history" , ppp ) ;
00333       if( code == NC_NOERR ){  /* should always happen */
00334          ppp[len] = '\0' ;
00335          tross_Append_History( dset , ppp ) ;
00336       }
00337       free(ppp) ;
00338    }
00339 
00340    /*-- close file and return the empty dataset */
00341 
00342    nc_close(ncid) ; RETURN(dset) ;
00343 }
00344 
00345 /*-----------------------------------------------------------------
00346   Load a MINC dataset's data into memory
00347   (called from THD_load_datablock in thd_loaddblk.c)
00348 -------------------------------------------------------------------*/
00349 
00350 void THD_load_minc( THD_datablock *dblk )
00351 {
00352    THD_diskptr *dkptr ;
00353    int im_varid , code , ncid ;
00354    int nx,ny,nz,nxy,nxyz,nxyzv , nbad,ibr,nv, nslice ;
00355    void *ptr ;
00356    nc_type im_type=NC_SHORT ;
00357    float im_valid_range[2], fac,denom , intop,inbot , outbot,outtop , sfac=1.0;
00358    float *im_min=NULL  , *im_max=NULL  ;
00359    int    im_min_varid ,  im_max_varid , do_scale ;
00360 
00361 ENTRY("THD_load_minc") ;
00362 
00363    /*-- open input [these errors should never occur] --*/
00364 
00365    if( !ISVALID_DATABLOCK(dblk)                       ||
00366        dblk->diskptr->storage_mode != STORAGE_BY_MINC ||
00367        dblk->brick == NULL                              ) EXRETURN ;
00368 
00369    dkptr = dblk->diskptr ;
00370 
00371    code = nc_open( dkptr->brick_name , NC_NOWRITE , &ncid ) ;
00372    if( code != NC_NOERR ){ EPR(code,dkptr->brick_name,"open"); EXRETURN; }
00373 
00374    code = nc_inq_varid( ncid , "image" , &im_varid ) ;
00375    if( code != NC_NOERR ){ EPR(code,dkptr->brick_name,"image varid"); nc_close(ncid); EXRETURN; }
00376 
00377    /*-- allocate space for data --*/
00378 
00379    nx = dkptr->dimsizes[0] ;
00380    ny = dkptr->dimsizes[1] ;  nxy   = nx * ny   ;
00381    nz = dkptr->dimsizes[2] ;  nxyz  = nxy * nz  ;
00382    nv = dkptr->nvals       ;  nxyzv = nxyz * nv ; nslice = nz*nv ;
00383 
00384    dblk->malloc_type = DATABLOCK_MEM_MALLOC ;
00385 
00386    /** malloc space for each brick separately **/
00387 
00388    for( nbad=ibr=0 ; ibr < nv ; ibr++ ){
00389       if( DBLK_ARRAY(dblk,ibr) == NULL ){
00390          ptr = AFMALL(void, DBLK_BRICK_BYTES(dblk,ibr) ) ;
00391          mri_fix_data_pointer( ptr ,  DBLK_BRICK(dblk,ibr) ) ;
00392          if( ptr == NULL ) nbad++ ;
00393       }
00394    }
00395 
00396    /** if couldn't get them all, take our ball and go home in a snit **/
00397 
00398    if( nbad > 0 ){
00399       fprintf(stderr,"\n** failed to malloc %d MINC bricks out of %d\n\a",nbad,nv) ;
00400       for( ibr=0 ; ibr < nv ; ibr++ ){
00401         if( DBLK_ARRAY(dblk,ibr) != NULL ){
00402            free(DBLK_ARRAY(dblk,ibr)) ;
00403            mri_fix_data_pointer( NULL , DBLK_BRICK(dblk,ibr) ) ;
00404         }
00405       }
00406       nc_close(ncid) ; EXRETURN ;
00407    }
00408 
00409    /** get scaling attributes/variables **/
00410 
00411    (void) nc_inq_vartype( ncid , im_varid , &im_type ) ;
00412 
00413    /* N.B.: we don't scale if input data is stored as floats */
00414 
00415    do_scale = (im_type==NC_BYTE || im_type==NC_SHORT || im_type==NC_INT) ;
00416 
00417    if( do_scale && AFNI_noenv("AFNI_MINC_SLICESCALE") ) do_scale = 0 ;
00418 
00419    code = nc_get_att_float( ncid,im_varid , "valid_range" , im_valid_range ) ;
00420 
00421    /** if can't get valid_range, make something up **/
00422 
00423    if( code != NC_NOERR || im_valid_range[0] >= im_valid_range[1] ){
00424 
00425       im_valid_range[0] = 0.0 ;
00426 
00427       switch( im_type ){
00428          case NC_BYTE:   im_valid_range[1] = 255.0        ; break ;
00429          case NC_SHORT:  im_valid_range[1] = 32767.0      ; break ;
00430          case NC_INT:    im_valid_range[1] = 2147483647.0 ; break ;
00431          case NC_FLOAT:
00432          case NC_DOUBLE: im_valid_range[1] = 1.0          ; break ;
00433       break ;
00434       }
00435    }
00436    inbot = im_valid_range[0] ;
00437    intop = im_valid_range[1] ;
00438    denom = intop - inbot  ;  /* always positive */
00439 
00440    /** Get range of image (per 2D slice) to which valid_range will be scaled **/
00441    /** Scaling will only be done if we get both image-min and image-max      **/
00442 
00443    if( do_scale ){
00444 
00445      code = nc_inq_varid( ncid , "image-min" , &im_min_varid ) ;
00446      if( code == NC_NOERR ){
00447        im_min = (float *) calloc(sizeof(float),nslice) ;
00448        code = nc_get_var_float( ncid , im_min_varid , im_min ) ;
00449        if( code != NC_NOERR ){ free(im_min); im_min = NULL; }
00450      }
00451 
00452      /* got im_min ==> try for im_max */
00453 
00454      if( im_min != NULL ){
00455        code = nc_inq_varid( ncid , "image-max" , &im_max_varid ) ;
00456        if( code == NC_NOERR ){
00457          im_max = (float *) calloc(sizeof(float),nslice) ;
00458          code = nc_get_var_float( ncid , im_max_varid , im_max ) ;
00459          if( code != NC_NOERR ){
00460            free(im_min); im_min = NULL; free(im_max); im_max = NULL;
00461          }
00462        }
00463      }
00464 
00465      /* 19 Mar 2003: make sure we don't scale out of range! */
00466 
00467      if( im_max != NULL && MRI_IS_INT_TYPE(DBLK_BRICK_TYPE(dblk,0)) ){
00468        float stop=0.0 , vbot,vtop ;
00469        int ii ;
00470        for( ii=0 ; ii < nslice ; ii++ ){
00471          if( im_min[ii] < im_max[ii] ){
00472            vbot = fabs(im_min[ii]) ; vtop = fabs(im_max[ii]) ;
00473            vtop = MAX(vtop,vbot) ; stop = MAX(vtop,stop) ;
00474          }
00475        }
00476        if( stop > MRI_TYPE_maxval[DBLK_BRICK_TYPE(dblk,0)] ){
00477          sfac = MRI_TYPE_maxval[DBLK_BRICK_TYPE(dblk,0)] / stop ;
00478          fprintf(stderr,"++ Scaling %s by %g to avoid overflow to %g of %s data!\n",
00479                  dkptr->brick_name, sfac, stop, MRI_TYPE_name[DBLK_BRICK_TYPE(dblk,0)] ) ;
00480        } else if( stop == 0.0 ){
00481          free(im_min); im_min = NULL; free(im_max); im_max = NULL;
00482        }
00483      }
00484 
00485    } /* end of do_scale */
00486 
00487    /** read data! **/
00488 
00489    switch( DBLK_BRICK_TYPE(dblk,0) ){  /* output datum */
00490 
00491       /* this case should never happen */
00492 
00493       default:
00494          fprintf(stderr,"\n** MINC illegal datum %d found\n",
00495                  DBLK_BRICK_TYPE(dblk,0)                     ) ;
00496          for( ibr=0 ; ibr < nv ; ibr++ ){
00497            if( DBLK_ARRAY(dblk,ibr) != NULL ){
00498               free(DBLK_ARRAY(dblk,ibr)) ;
00499               mri_fix_data_pointer( NULL , DBLK_BRICK(dblk,ibr) ) ;
00500            }
00501          }
00502       break ;
00503 
00504       /* data is stored as bytes in AFNI */
00505 
00506       case MRI_byte:{
00507          int kk,ii,qq ; byte *br ;
00508 
00509          if( nv == 1 ){
00510             code = nc_get_var_uchar( ncid , im_varid , DBLK_ARRAY(dblk,0) ) ;
00511          } else {
00512             size_t start[4] , count[4] ;
00513             start[1] = start[2] = start[3] = 0 ;
00514             count[0] = 1 ; count[1] = nz ; count[2] = ny ; count[3] = nx ;
00515             for( ibr=0 ; ibr < nv ; ibr++ ){
00516                start[0] = ibr ;
00517                code = nc_get_vara_uchar( ncid,im_varid ,
00518                                          start,count , DBLK_ARRAY(dblk,ibr) ) ;
00519             }
00520          }
00521          if( code != NC_NOERR ) EPR(code,dkptr->brick_name,"image get_var_uchar") ;
00522 
00523          if( im_max != NULL ){                 /* must scale values */
00524            for( ibr=0 ; ibr < nv ; ibr++ ){     /* loop over bricks */
00525              br = DBLK_ARRAY(dblk,ibr) ;
00526              for( kk=0 ; kk < nz ; kk++ ){      /* loop over slices */
00527                 outbot = im_min[kk+ibr*nz] ;
00528                 outtop = im_max[kk+ibr*nz] ;
00529                 if( outbot >= outtop) continue ;            /* skip */
00530                 fac = (outtop-outbot) / denom ;
00531                 qq = kk*nxy ; outbot += 0.499 ;
00532                 for( ii=0 ; ii < nxy ; ii++ )        /* scale slice */
00533                    br[ii+qq] = (byte) ((fac*(br[ii+qq]-inbot) + outbot)*sfac) ;
00534              }
00535            }
00536          }
00537       }
00538       break ;
00539 
00540       /* data is stored as shorts in AFNI */
00541 
00542       case MRI_short:{
00543          int kk,ii,qq ; short *br ;
00544 
00545          if( nv == 1 ){
00546             code = nc_get_var_short( ncid , im_varid , DBLK_ARRAY(dblk,0) ) ;
00547          } else {
00548             size_t start[4] , count[4] ;
00549             start[1] = start[2] = start[3] = 0 ;
00550             count[0] = 1 ; count[1] = nz ; count[2] = ny ; count[3] = nx ;
00551             for( ibr=0 ; ibr < nv ; ibr++ ){
00552                start[0] = ibr ;
00553                code = nc_get_vara_short( ncid,im_varid ,
00554                                          start,count , DBLK_ARRAY(dblk,ibr) ) ;
00555             }
00556          }
00557          if( code != NC_NOERR ) EPR(code,dkptr->brick_name,"image get_var_short") ;
00558 
00559          if( im_max != NULL ){                 /* must scale values */
00560            for( ibr=0 ; ibr < nv ; ibr++ ){     /* loop over bricks */
00561              br = DBLK_ARRAY(dblk,ibr) ;
00562              for( kk=0 ; kk < nz ; kk++ ){      /* loop over slices */
00563                 outbot = im_min[kk+ibr*nz] ;
00564                 outtop = im_max[kk+ibr*nz] ;
00565                 if( outbot >= outtop) continue ;            /* skip */
00566                 fac = (outtop-outbot) / denom ;
00567                 qq = kk*nxy ; outbot += 0.499 ;
00568                 for( ii=0 ; ii < nxy ; ii++ )        /* scale slice */
00569                    br[ii+qq] = (short) ((fac*(br[ii+qq]-inbot) + outbot)*sfac) ;
00570              }
00571            }
00572          }
00573       }
00574       break ;
00575 
00576       /* data is stored as floats in AFNI */
00577 
00578       case MRI_float:{
00579          int kk,ii,qq ; float *br ;
00580 
00581          if( nv == 1 ){
00582             code = nc_get_var_float( ncid , im_varid , DBLK_ARRAY(dblk,0) ) ;
00583          } else {
00584             size_t start[4] , count[4] ;
00585             start[1] = start[2] = start[3] = 0 ;
00586             count[0] = 1 ; count[1] = nz ; count[2] = ny ; count[3] = nx ;
00587             for( ibr=0 ; ibr < nv ; ibr++ ){
00588                start[0] = ibr ;
00589                code = nc_get_vara_float( ncid,im_varid ,
00590                                          start,count , DBLK_ARRAY(dblk,ibr) ) ;
00591             }
00592          }
00593          if( code != NC_NOERR ) EPR(code,dkptr->brick_name,"image get_var_float") ;
00594 
00595          if( im_type == NC_BYTE ){              /* fix sign of bytes */
00596            for( ibr=0 ; ibr < nv ; ibr++ ){     /* loop over bricks */
00597              br = DBLK_ARRAY(dblk,ibr) ;
00598              for( ii=0 ; ii < nxyz ; ii++ )
00599                 if( br[ii] < 0.0 ) br[ii] += 256.0 ;
00600            }
00601          }
00602 
00603          if( im_max != NULL ){                 /* must scale values */
00604            for( ibr=0 ; ibr < nv ; ibr++ ){     /* loop over bricks */
00605              br = DBLK_ARRAY(dblk,ibr) ;
00606              for( kk=0 ; kk < nz ; kk++ ){      /* loop over slices */
00607                 outbot = im_min[kk+ibr*nz] ;
00608                 outtop = im_max[kk+ibr*nz] ;
00609                 if( outbot >= outtop) continue ;            /* skip */
00610                 fac = (outtop-outbot) / denom ;
00611                 qq = kk*nxy ;
00612                 for( ii=0 ; ii < nxy ; ii++ )        /* scale slice */
00613                    br[ii+qq] = (fac*(br[ii+qq]-inbot) + outbot) ;
00614              }
00615            }
00616          }
00617       }
00618       break ;
00619 
00620    } /* end of switch on output datum */
00621 
00622    /*-- throw away the trash and return --*/
00623 
00624    if( im_min != NULL ) free(im_min) ;
00625    if( im_max != NULL ) free(im_max) ;
00626    nc_close(ncid) ; EXRETURN ;
00627 }