#include "mrilib.h"
/****************************************************************************
Functions for dealing with data extracted/inserted from rows of a 3D
dataset - parallel to x, y, or z axes.
-- 08 Mar 2001 - RWCox
*****************************************************************************/
/*---------------------------------------------------------------------------*/
/*! Get the number of elements in a row of the dataset in a particular
direction given by dcode:
- +1 = x direction -1 = reversed x direction
- +2 = y direction -2 = reversed y direction
- +3 = z direction -3 = reversed z direction
This routine is mostly for the pitiful programmer who can't be bothered
to do this calculation himself.
-----------------------------------------------------------------------------*/
int THD_get_dset_rowcount( THD_3dim_dataset *dset, int dcode )
{
if( !ISVALID_DSET(dset) ) return 0 ; /* bad */
switch( dcode ){
case 1: case -1: return DSET_NX(dset) ; /* good */
case 2: case -2: return DSET_NY(dset) ; /* good */
case 3: case -3: return DSET_NZ(dset) ; /* good */
}
return 0 ; /* bad */
}
/*---------------------------------------------------------------------------*/
/*! Extract a row from a dataset sub-brick in the direction given by dcode.
The 3-index of a voxel from the row is in (xx,yy,zz). The return value
is a pointer to a malloc()-ed array, whose type is given by
DSET_BRICK_TYPE(dset,ival) and whose length is given by
THD_get_dset_rowcount(dset,dcode). If NULL is returned, something bad
happened.
N.B.: dcode < 0 ==> data is extracted in the reverse direction.
-----------------------------------------------------------------------------*/
void * THD_get_dset_row( THD_3dim_dataset *dset, int ival,
int dcode , int xx,int yy,int zz )
{
void *row , *brick ;
int nrow , kind , nx,ny,nz,nxy , kbot,kdel,kk,ii ;
ENTRY("THD_get_dset_row") ;
nrow = THD_get_dset_rowcount( dset , dcode ) ;
if( nrow < 1 ) RETURN(NULL) ;
nx = DSET_NX(dset) ;
ny = DSET_NY(dset) ; nxy = nx*ny ;
nz = DSET_NZ(dset) ;
/*-- We will extract brick[kbot+i*kdel] for i=0..nrow-1 --*/
switch( dcode ){
case 1: case -1:
if( yy < 0 || yy >= ny || zz < 0 || zz >= nz ) RETURN(NULL) ;
kbot = yy*nx + zz*nxy ; kdel = 1 ;
break ;
case 2: case -2:
if( xx < 0 || xx >= nx || zz < 0 || zz >= nz ) RETURN(NULL) ;
kbot = xx + zz*nxy ; kdel = nx ;
break ;
case 3: case -3:
if( xx < 0 || xx >= nx || yy < 0 || yy >= ny ) RETURN(NULL) ;
kbot = xx + yy*nx ; kdel = nxy ;
break ;
}
kind = DSET_BRICK_TYPE(dset,ival) ;
brick = DSET_ARRAY(dset,ival) ;
row = AFMALL(void, mri_datum_size((MRI_TYPE)kind) * nrow ) ;
/*-- extract row, based on kind of data in sub-brick --*/
switch( kind ){
default: free(row) ; RETURN(NULL) ; /* bad */
case MRI_short:{
short *rr = (short *)row , *bb = (short *)brick ;
if( dcode > 0 )
for( ii=kk=0 ; ii < nrow ; ii++,kk+=kdel ) rr[ii] = bb[kk+kbot] ;
else
for( ii=kk=0 ; ii < nrow ; ii++,kk+=kdel ) rr[nrow-1-ii] = bb[kk+kbot] ;
}
break ;
case MRI_byte:{
byte *rr = (byte *)row , *bb = (byte *)brick ;
if( dcode > 0 )
for( ii=kk=0 ; ii < nrow ; ii++,kk+=kdel ) rr[ii] = bb[kk+kbot] ;
else
for( ii=kk=0 ; ii < nrow ; ii++,kk+=kdel ) rr[nrow-1-ii] = bb[kk+kbot] ;
}
break ;
case MRI_float:{
float *rr = (float *)row , *bb = (float *)brick ;
if( dcode > 0 )
for( ii=kk=0 ; ii < nrow ; ii++,kk+=kdel ) rr[ii] = bb[kk+kbot] ;
else
for( ii=kk=0 ; ii < nrow ; ii++,kk+=kdel ) rr[nrow-1-ii] = bb[kk+kbot] ;
}
break ;
case MRI_complex:{
complex *rr = (complex *)row , *bb = (complex *)brick ;
if( dcode > 0 )
for( ii=kk=0 ; ii < nrow ; ii++,kk+=kdel ) rr[ii] = bb[kk+kbot] ;
else
for( ii=kk=0 ; ii < nrow ; ii++,kk+=kdel ) rr[nrow-1-ii] = bb[kk+kbot] ;
}
break ;
case MRI_rgb:{
byte *rr = (byte *)row , *bb = (byte *)brick ;
if( dcode > 0 )
for( ii=kk=0 ; ii < nrow ; ii++,kk+=kdel ){
rr[3*ii ] = bb[3*(kk+kbot) ] ;
rr[3*ii+1] = bb[3*(kk+kbot)+1] ;
rr[3*ii+2] = bb[3*(kk+kbot)+2] ;
}
else
for( ii=kk=0 ; ii < nrow ; ii++,kk+=kdel ){
rr[3*(nrow-1-ii) ] = bb[3*(kk+kbot) ] ;
rr[3*(nrow-1-ii)+1] = bb[3*(kk+kbot)+1] ;
rr[3*(nrow-1-ii)+2] = bb[3*(kk+kbot)+2] ;
}
}
break ;
}
RETURN(row) ;
}
/*---------------------------------------------------------------------------*/
/*! Get a row, but always return a float MRI_IMAGE, and scaled if need be.
-----------------------------------------------------------------------------*/
MRI_IMAGE * mri_get_dset_row( THD_3dim_dataset *dset, int ival,
int dcode , int xx,int yy,int zz )
{
void *rawrow ;
MRI_IMAGE *im ;
float *fim , fac ;
int ii , nrow , kind ;
ENTRY("MRI_get_dset_row") ;
rawrow = THD_get_dset_row( dset , ival , dcode , xx,yy,zz ) ;
if( rawrow == NULL ) RETURN(NULL) ;
nrow = THD_get_dset_rowcount( dset , dcode ) ;
kind = DSET_BRICK_TYPE(dset,ival) ;
fac = DSET_BRICK_FACTOR(dset,ival); if( fac <= 0.0f ) fac = 1.0f;
switch( kind ){
case MRI_float:{
fim = (float *)rawrow ;
if( fac != 1.0f ) for( ii=0 ; ii < nrow ; ii++ ) fim[ii] *= fac ;
im = mri_new_vol_empty( nrow,1,1 , MRI_float ) ;
mri_fix_data_pointer( fim , im ) ;
}
break ;
case MRI_short:{
short *rr = (short *)rawrow ;
im = mri_new( nrow , 1 , MRI_float ); fim = MRI_FLOAT_PTR(im);
for( ii=0 ; ii < nrow ; ii++ ) fim[ii] = rr[ii] * fac ;
}
break ;
case MRI_byte:{
byte *rr = (byte *)rawrow ;
im = mri_new( nrow , 1 , MRI_float ); fim = MRI_FLOAT_PTR(im);
for( ii=0 ; ii < nrow ; ii++ ) fim[ii] = rr[ii] * fac ;
}
break ;
case MRI_complex:{
complex *rr = (complex *)rawrow ;
im = mri_new( nrow , 1 , MRI_float ); fim = MRI_FLOAT_PTR(im);
for( ii=0 ; ii < nrow ; ii++ ) fim[ii] = CABS(rr[ii]) ;
}
break ;
case MRI_rgb:{
byte *rr = (byte *)rawrow ;
im = mri_new( nrow , 1 , MRI_float ); fim = MRI_FLOAT_PTR(im);
for( ii=0 ; ii < nrow ; ii++ ) fim[ii] = 0.299 * rr[3*ii ]
+ 0.587 * rr[3*ii+1]
+ 0.114 * rr[3*ii+2] ;
}
break ;
}
if( rawrow != (void *)fim ) free(rawrow) ;
RETURN(im) ;
}
/*---------------------------------------------------------------------------
The inverse of THD_get_dset_row: put data back into a dataset.
-----------------------------------------------------------------------------*/
void THD_put_dset_row( THD_3dim_dataset *dset, int ival,
int dcode, int xx,int yy,int zz, void *row )
{
void *brick ;
int nrow , kind , nx,ny,nz,nxy , kbot,kdel,kk,ii ;
ENTRY("THD_put_dset_row") ;
nrow = THD_get_dset_rowcount( dset , dcode ) ;
if( nrow < 1 ) EXRETURN ;
nx = DSET_NX(dset) ;
ny = DSET_NY(dset) ; nxy = nx*ny ;
nz = DSET_NZ(dset) ;
/*-- We will insert brick[kbot+i*kdel] for i=0..nrow-1 --*/
switch( dcode ){
case 1: case -1:
if( yy < 0 || yy >= ny || zz < 0 || zz >= nz ) EXRETURN ;
kbot = yy*nx + zz*nxy ; kdel = 1 ;
break ;
case 2: case -2:
if( xx < 0 || xx >= nx || zz < 0 || zz >= nz ) EXRETURN ;
kbot = xx + zz*nxy ; kdel = nx ;
break ;
case 3: case -3:
if( xx < 0 || xx >= nx || yy < 0 || yy >= ny ) EXRETURN ;
kbot = xx + yy*nx ; kdel = nxy ;
break ;
}
kind = DSET_BRICK_TYPE(dset,ival) ;
brick = DSET_ARRAY(dset,ival) ;
/*-- insert row, based on kind of data in sub-brick --*/
switch( kind ){
default: EXRETURN ; /* bad */
case MRI_short:{
short *rr = (short *)row , *bb = (short *)brick ;
if( dcode > 0 )
for( ii=kk=0 ; ii < nrow ; ii++,kk+=kdel ) bb[kk+kbot] = rr[ii] ;
else
for( ii=kk=0 ; ii < nrow ; ii++,kk+=kdel ) bb[kk+kbot] = rr[nrow-1-ii] ;
}
break ;
case MRI_byte:{
byte *rr = (byte *)row , *bb = (byte *)brick ;
if( dcode > 0 )
for( ii=kk=0 ; ii < nrow ; ii++,kk+=kdel ) bb[kk+kbot] = rr[ii] ;
else
for( ii=kk=0 ; ii < nrow ; ii++,kk+=kdel ) bb[kk+kbot] = rr[nrow-1-ii] ;
}
break ;
case MRI_complex:{
complex *rr = (complex *)row , *bb = (complex *)brick ;
if( dcode > 0 )
for( ii=kk=0 ; ii < nrow ; ii++,kk+=kdel ) bb[kk+kbot] = rr[ii] ;
else
for( ii=kk=0 ; ii < nrow ; ii++,kk+=kdel ) bb[kk+kbot] = rr[nrow-1-ii] ;
}
break ;
case MRI_float:{
float *rr = (float *)row , *bb = (float *)brick ;
if( dcode > 0 )
for( ii=kk=0 ; ii < nrow ; ii++,kk+=kdel ) bb[kk+kbot] = rr[ii] ;
else
for( ii=kk=0 ; ii < nrow ; ii++,kk+=kdel ) bb[kk+kbot] = rr[nrow-1-ii] ;
}
break ;
case MRI_rgb:{
byte *rr = (byte *)row , *bb = (byte *)brick ;
if( dcode > 0 )
for( ii=kk=0 ; ii < nrow ; ii++,kk+=kdel ){
bb[3*(kk+kbot) ] = rr[3*ii ] ;
bb[3*(kk+kbot)+1] = rr[3*ii+1] ;
bb[3*(kk+kbot)+2] = rr[3*ii+2] ;
}
else
for( ii=kk=0 ; ii < nrow ; ii++,kk+=kdel ){
bb[3*(kk+kbot) ] = rr[3*(nrow-1-ii) ] ;
bb[3*(kk+kbot)+1] = rr[3*(nrow-1-ii)+1] ;
bb[3*(kk+kbot)+2] = rr[3*(nrow-1-ii)+2] ;
}
}
break ;
}
EXRETURN ;
}
syntax highlighted by Code2HTML, v. 0.9.1