/* ----------------------------- MNI Header ----------------------------------- @NAME : dim_conversion.c @DESCRIPTION: File of functions to do dimension conversion for image conversion variables (icv). These variables allow conversion of netcdf variables (the MINC image variable, in particular) to a form more convenient for a program. @METHOD : Routines included in this file : public : miicv_attach private : MI_icv_get_dim MI_get_dim_flip MI_get_dim_scale MI_get_dim_bufsize_step MI_icv_get_dim_conversion MI_icv_dimconvert MI_icv_dimconv_init @CREATED : September 9, 1992. (Peter Neelin) @MODIFIED : * $Log: dim_conversion.c,v $ * Revision 6.4 2004/10/15 13:44:52 bert * Minor changes for Windows compatibility * * Revision 6.3 2003/11/14 16:52:24 stever * More last-minute fixes. * * Revision 6.2 2003/09/18 16:16:15 bert * Use standard labs() and fabs() instead of our private macros * * Revision 6.1 1999/10/19 14:45:07 neelin * Fixed Log subsitutions for CVS * * Revision 6.0 1997/09/12 13:24:54 neelin * Release of minc version 0.6 * * Revision 5.0 1997/08/21 13:25:53 neelin * Release of minc version 0.5 * * Revision 4.0 1997/05/07 20:07:52 neelin * Release of minc version 0.4 * * Revision 3.1 1997/04/10 18:14:50 neelin * Fixed handling of invalid data when icv scale is zero. * * Revision 3.0 1995/05/15 19:33:12 neelin * Release of minc version 0.3 * * Revision 2.3 1995/02/08 19:14:44 neelin * More changes for irix 5 lint. * * Revision 2.2 1995/02/08 19:01:06 neelin * Moved private function declarations from minc_routines.h to appropriate file. * * Revision 2.1 1994/11/02 09:42:37 neelin * Fixed conversion of vector to scalar (old code simply returned the first * component of the vector - now averaging is done properly). * * Revision 2.0 94/09/28 10:37:52 neelin * Release of minc version 0.2 * * Revision 1.11 94/09/28 10:37:35 neelin * Pre-release * * Revision 1.10 94/07/05 15:31:07 neelin * Assume that MIstep is positive if it is not found (for flipping dimensions). * * Revision 1.9 94/03/16 10:30:00 neelin * Changed default MI_ICV_STEP: if not found use 1.0 (and 0.0 for start). * * Revision 1.8 93/08/11 12:06:03 neelin * Added RCS logging in source. * @COPYRIGHT : Copyright 1993 Peter Neelin, McConnell Brain Imaging Centre, Montreal Neurological Institute, McGill University. Permission to use, copy, modify, and distribute this software and its documentation for any purpose and without fee is hereby granted, provided that the above copyright notice appear in all copies. The author and McGill University make no representations about the suitability of this software for any purpose. It is provided "as is" without express or implied warranty. ---------------------------------------------------------------------------- */ #ifndef lint static char rcsid[] = "$Header: /software/source/minc/cvsroot/minc/libsrc/dim_conversion.c,v 6.4 2004/10/15 13:44:52 bert Exp $ MINC (MNI)"; #endif #include "minc_private.h" #include #include /* Private functions */ PRIVATE int MI_icv_get_dim(mi_icv_type *icvp, int cdfid, int varid); PRIVATE int MI_get_dim_flip(mi_icv_type *icvp, int cdfid, int dimvid[], int subsc[]); PRIVATE int MI_get_dim_scale(mi_icv_type *icvp, int cdfid, int dimvid[]); PRIVATE int MI_get_dim_bufsize_step(mi_icv_type *icvp, int subsc[]); PRIVATE int MI_icv_get_dim_conversion(mi_icv_type *icvp, int subsc[]); PRIVATE int MI_icv_dimconvert(int operation, mi_icv_type *icvp, long start[], long count[], void *values, long bufstart[], long bufcount[], void *buffer); PRIVATE int MI_icv_dimconv_init(int operation, mi_icv_type *icvp, mi_icv_dimconv_type *dcp, long start[], long count[], void *values, long bufstart[], long bufcount[], void *buffer); /* ----------------------------- MNI Header ----------------------------------- @NAME : miicv_attach @INPUT : icvid - icv id cdfid - cdf file id varid - cdf variable id @OUTPUT : (none) @RETURNS : MI_ERROR if an error occurs @DESCRIPTION: Attaches an open cdf file and variable to an image conversion variable for subsequent access through miicvget and miicvput. File must be in data mode. @METHOD : @GLOBALS : @CALLS : NetCDF routines @CREATED : September 9, 1992 (Peter Neelin) @MODIFIED : ---------------------------------------------------------------------------- */ MNCAPI int miicv_attach(int icvid, int cdfid, int varid) { mi_icv_type *icvp; /* Pointer to icv structure */ long size_diff, user_dim_size; int idim; MI_SAVE_ROUTINE_NAME("miicv_attach"); /* Check icv id */ if ((icvp=MI_icv_chkid(icvid)) == NULL) MI_RETURN_ERROR(MI_ERROR); /* Call routine to set variables for everything except dimension conversion */ {MI_CHK_ERR(miicv_ndattach(icvid, cdfid, varid))} /* Check to see if we need to worry about dimension conversion */ if (!icvp->user_do_dimconv) { MI_RETURN(MI_NOERROR); } /* Reset cdfid and varid in icv structure in case something goes wrong in dimension conversion calculations */ icvp->cdfid = MI_ERROR; icvp->varid = MI_ERROR; /* Get dimensions info */ {MI_CHK_ERR(MI_icv_get_dim(icvp, cdfid, varid))} /* Set the do_dimconvert field of icv structure We do dimension conversion if any dimension needs flipping, scaling or offset, or if we have to convert from vector to scalar. */ icvp->do_dimconvert = (icvp->user_do_scalar && icvp->var_is_vector); for (idim=0; idimuser_num_imgdims; idim++) { if (icvp->derv_dim_flip[idim] || (icvp->derv_dim_scale[idim] != 1) || (icvp->derv_dim_off[idim] != 0)) icvp->do_dimconvert = TRUE; } icvp->dimconvert_func = MI_icv_dimconvert; /* Check if we have to zero user's buffer on GETs */ icvp->derv_do_zero = FALSE; for (idim=0; idimuser_num_imgdims; idim++) { user_dim_size = ((icvp->user_dim_size[idim]<=0) ? icvp->var_dim_size[idim] : icvp->user_dim_size[idim]); if (icvp->derv_dim_grow[idim]) size_diff = user_dim_size - icvp->var_dim_size[idim] * icvp->derv_dim_scale[idim]; else size_diff = user_dim_size - 1 - (icvp->var_dim_size[idim] - 1) / icvp->derv_dim_scale[idim]; if ((icvp->derv_dim_off[idim]!=0) || (size_diff!=0)) icvp->derv_do_zero = TRUE; } /* Set the cdfid and varid fields */ icvp->cdfid = cdfid; icvp->varid = varid; MI_RETURN(MI_NOERROR); } /* ----------------------------- MNI Header ----------------------------------- @NAME : MI_icv_get_dim @INPUT : icvp - pointer to icv structure cdfid - cdf file id varid - variable id @OUTPUT : (none) @RETURNS : MI_ERROR if an error occurs @DESCRIPTION: Gets dimension info for the icv @METHOD : @GLOBALS : @CALLS : NetCDF routines @CREATED : August 10, 1992 (Peter Neelin) @MODIFIED : ---------------------------------------------------------------------------- */ PRIVATE int MI_icv_get_dim(mi_icv_type *icvp, int cdfid, int varid) /* ARGSUSED */ { int oldncopts; /* For saving value of ncopts */ char dimname[MAX_NC_NAME]; /* Dimensions name */ int idim; /* Looping counter for fastest image dims */ int subsc[MI_MAX_IMGDIMS]; /* Subscripts for fastest image dims */ int dimvid[MI_MAX_IMGDIMS]; /* Variable ids for dimensions */ MI_SAVE_ROUTINE_NAME("MI_icv_get_dim"); /* Check that the variable has at least icvp->user_num_imgdims dimensions */ if (icvp->var_ndims < icvp->user_num_imgdims) { MI_LOG_PKG_ERROR2(MI_ERR_TOOFEWDIMS, "Variable has too few dimensions to be an image"); MI_RETURN_ERROR(MI_ERROR); } /* Check the first dimensions of the variable */ MI_CHK_ERR(ncdiminq(cdfid, icvp->var_dim[icvp->var_ndims-1], dimname, &(icvp->var_vector_size))) icvp->var_is_vector = STRINGS_EQUAL(dimname, MIvector_dimension); /* Check that the variable has at least icvp->user_num_imgdims+1 dimensions if it is a vector field */ if (icvp->var_is_vector && (icvp->var_ndims < icvp->user_num_imgdims+1)) { MI_LOG_PKG_ERROR2(MI_ERR_TOOFEWDIMS, "Variable has too few dimensions to be an image"); MI_RETURN_ERROR(MI_ERROR); } /* Check for dimension flipping and get dimension sizes */ /* Get subscripts for first icvp->user_num_imgdims dimensions */ subsc[0] = (icvp->var_is_vector) ? icvp->var_ndims-2 : icvp->var_ndims-1; for (idim=1; idim < icvp->user_num_imgdims; idim++) subsc[idim]=subsc[idim-1]-1; /* Get dimension variable ids */ for (idim=0; idim < icvp->user_num_imgdims; idim++) { {MI_CHK_ERR(ncdiminq(cdfid, icvp->var_dim[subsc[idim]], dimname, &(icvp->var_dim_size[idim])))}; oldncopts = ncopts; ncopts = 0; dimvid[idim] = ncvarid(cdfid, dimname); ncopts = oldncopts; } /* Check for flipping */ {MI_CHK_ERR(MI_get_dim_flip(icvp, cdfid, dimvid, subsc))} /* Check for scaling of dimension */ {MI_CHK_ERR(MI_get_dim_scale(icvp, cdfid, dimvid))} /* Check for variable buffer size increments */ {MI_CHK_ERR(MI_get_dim_bufsize_step(icvp, subsc))} /* Get information for dimension conversion */ {MI_CHK_ERR(MI_icv_get_dim_conversion(icvp, subsc))} MI_RETURN(MI_NOERROR); } /* ----------------------------- MNI Header ----------------------------------- @NAME : MI_icv_get_dim_flip @INPUT : icvp - icv pointer cdfid - cdf file id dimvid - variable id subsc - array of dimension subscripts for fastest varying dimensions @OUTPUT : (none) @RETURNS : MI_ERROR if an error occurs @DESCRIPTION: Checks for flipping of icv. @METHOD : @GLOBALS : @CALLS : NetCDF routines @CREATED : September 1, 1992 (Peter Neelin) @MODIFIED : ---------------------------------------------------------------------------- */ PRIVATE int MI_get_dim_flip(mi_icv_type *icvp, int cdfid, int dimvid[], int subsc[]) { int oldncopts; /* For saving value of ncopts */ char dimname[MAX_NC_NAME]; /* Dimensions name */ int dim_dir; /* Desired direction for current dimension */ double dimstep; /* Dimension step size (and direction) */ int idim; MI_SAVE_ROUTINE_NAME("MI_get_dim_flip"); /* Loop through fast dimensions */ for (idim=0; idim < icvp->user_num_imgdims; idim++) { /* Get name of the dimension */ {MI_CHK_ERR(ncdiminq(cdfid, icvp->var_dim[subsc[idim]], dimname, NULL))} /* Should we look for dimension flipping? */ icvp->derv_dim_flip[idim]=FALSE; if (STRINGS_EQUAL(dimname, MIxspace) || STRINGS_EQUAL(dimname, MIxfrequency)) dim_dir = icvp->user_xdim_dir; else if (STRINGS_EQUAL(dimname, MIyspace) || STRINGS_EQUAL(dimname, MIyfrequency)) dim_dir = icvp->user_ydim_dir; else if (STRINGS_EQUAL(dimname, MIzspace) || STRINGS_EQUAL(dimname, MIzfrequency)) dim_dir = icvp->user_zdim_dir; else dim_dir = MI_ICV_ANYDIR; /* Look for variable corresponding to dimension */ if (dim_dir != MI_ICV_ANYDIR) { /* Look for flipping? */ /* Get the value of the MIstep attribute to determine whether flipping is needed. Assume that direction is positive if no step is provided. */ dimstep = 1.0; if (dimvid[idim] != MI_ERROR) { /* if dimension exists */ oldncopts = ncopts; ncopts = 0; (void) miattget1(cdfid, dimvid[idim], MIstep, NC_DOUBLE, &dimstep); ncopts = oldncopts; } /* if dimension exists */ if (dim_dir == MI_ICV_POSITIVE) icvp->derv_dim_flip[idim] = (dimstep<0.0); else if (dim_dir == MI_ICV_NEGATIVE) icvp->derv_dim_flip[idim] = (dimstep>0.0); } /* if look for flipping */ } /* for each dimension */ MI_RETURN(MI_NOERROR); } /* ----------------------------- MNI Header ----------------------------------- @NAME : MI_icv_get_dim_scale @INPUT : icvp - icv pointer cdfid - cdf file id dimvid - dimension variable id @OUTPUT : (none) @RETURNS : MI_ERROR if an error occurs @DESCRIPTION: Checks for scaling of images @METHOD : @GLOBALS : @CALLS : NetCDF routines @CREATED : September 1, 1992 (Peter Neelin) @MODIFIED : ---------------------------------------------------------------------------- */ PRIVATE int MI_get_dim_scale(mi_icv_type *icvp, int cdfid, int dimvid[]) { int oldncopts; /* For saving value of ncopts */ int min_grow, dim_grow; int min_scale, dim_scale; double dimstep, dimstart; int idim; long user_dim_size; MI_SAVE_ROUTINE_NAME("MI_get_dim_scale"); /* Loop through dimensions, calculating scale and looking for smallest one. For each dimension, check to see if we need to shrink or grow the image. (This is get-oriented: grow is TRUE if the variable dimension has to be expanded to fit the user's dimensions). */ for (idim=0; idim < icvp->user_num_imgdims; idim++) { /* Check to see if we user wants resize */ if (icvp->user_dim_size[idim] <= 0) { icvp->derv_dim_grow[idim]=TRUE; icvp->derv_dim_scale[idim]=1; } else { /* Check for growing or shrinking */ icvp->derv_dim_grow[idim] = (icvp->var_dim_size[idim] <= icvp->user_dim_size[idim]); /* If growing */ if (icvp->derv_dim_grow[idim]) { /* Get scale so that whole image fits in user array */ icvp->derv_dim_scale[idim] = icvp->user_dim_size[idim] / icvp->var_dim_size[idim]; } /* Otherwise, shrinking. Things are complicated by the fact that the external variable must fit completely in the user's array */ else { icvp->derv_dim_scale[idim] = 1 + (icvp->var_dim_size[idim] - 1) / icvp->user_dim_size[idim]; } } /* if user wants resizing */ /* Check for smallest scale */ if (idim==0) { min_grow = icvp->derv_dim_grow[idim]; min_scale = icvp->derv_dim_scale[idim]; } else { dim_grow = icvp->derv_dim_grow[idim]; dim_scale = icvp->derv_dim_scale[idim]; /* Check for one of three conditions : (1) smallest so far is growing, but this dim is shrinking (2) both are growing and this dim has smaller scale (3) both are shrinking and this dim has larger scale */ if ((min_grow && !dim_grow) || ((min_grow && dim_grow) && (min_scale > dim_scale)) || ((!min_grow && !dim_grow) && (min_scale < dim_scale))) { min_grow = dim_grow; min_scale = dim_scale; } } /* if not first dim */ } /* for each dimension, get scale */ /* Loop through dimensions, resetting scale if needed, setting offset and pixel step and start */ for (idim=0; idim < icvp->user_num_imgdims; idim++) { /* Check for aspect ratio */ if (icvp->user_keep_aspect) { icvp->derv_dim_grow[idim] = min_grow; icvp->derv_dim_scale[idim] = min_scale; } /* Get user's buffer size */ user_dim_size = ((icvp->user_dim_size[idim]<=0) ? icvp->var_dim_size[idim] : icvp->user_dim_size[idim]); /* Set offset of variable into user's image */ /* If growing */ if (icvp->derv_dim_grow[idim]) { /* Calculate remainder and split it in half */ icvp->derv_dim_off[idim] = ( user_dim_size - icvp->var_dim_size[idim] * icvp->derv_dim_scale[idim] ) / 2; } /* Otherwise, shrinking. Things are complicated by the fact that the external variable must fit completely in the user's array */ else { /* Calculate remainder and split it in half */ icvp->derv_dim_off[idim] = ( user_dim_size - 1 - (icvp->var_dim_size[idim] - 1) / icvp->derv_dim_scale[idim] ) / 2 ; } /* Get pixel step and start for variable and calculate for user. Look for them in the dimension variable (if MIstep is not there, then use defaults step = 1.0, start = 0.0 */ oldncopts = ncopts; ncopts = 0; dimstep = 1.0; (void) miattget1(cdfid, dimvid[idim], MIstep, NC_DOUBLE, &dimstep); /* Flip dimstep if needed */ if (icvp->derv_dim_flip[idim]) dimstep *= (-1); /* Get step size for user's image */ icvp->derv_dim_step[idim] = icvp->derv_dim_grow[idim] ? dimstep / icvp->derv_dim_scale[idim] : dimstep * icvp->derv_dim_scale[idim]; /* Get start position for user's image - if no MIstart for dimension, then assume 0.0 */ dimstart = 0.0; (void) miattget1(cdfid, dimvid[idim], MIstart, NC_DOUBLE, &dimstart); /* Flip dimstart if needed */ if (icvp->derv_dim_flip[idim]) dimstart -= dimstep * (icvp->var_dim_size[idim]-1); /* Calculate start position */ icvp->derv_dim_start[idim] = dimstart + (icvp->derv_dim_step[idim] - dimstep) / 2.0 - icvp->derv_dim_off[idim] * icvp->derv_dim_step[idim]; ncopts = oldncopts; } /* for each dimension */ MI_RETURN(MI_NOERROR); } /* ----------------------------- MNI Header ----------------------------------- @NAME : MI_icv_get_dim_bufsize_step @INPUT : icvp - icv pointer subsc - array of dimension subscripts for fastest varying dimensions @OUTPUT : (none) @RETURNS : MI_ERROR if an error occurs @DESCRIPTION: Sets the variables giving variable buffer size @METHOD : @GLOBALS : @CALLS : NetCDF routines @CREATED : September 3, 1992 (Peter Neelin) @MODIFIED : ---------------------------------------------------------------------------- */ PRIVATE int MI_get_dim_bufsize_step(mi_icv_type *icvp, int subsc[]) { int idim; MI_SAVE_ROUTINE_NAME("MI_get_dim_bufsize_step"); /* Set default buffer size step */ for (idim=0; idim < MAX_VAR_DIMS; idim++) icvp->derv_bufsize_step[idim]=1; /* Check for converting vector to scalar */ icvp->derv_do_bufsize_step = (icvp->var_is_vector && icvp->user_do_scalar); if (icvp->derv_do_bufsize_step) icvp->derv_bufsize_step[icvp->var_ndims-1] = icvp->var_vector_size; /* Check each dimension to see if we need to worry about the variable buffer size. This only occurs if we are shrinking the dimension from the variable buffer to the user buffer. */ for (idim=0; idim < icvp->user_num_imgdims; idim++) { if (!icvp->derv_dim_grow[idim]) icvp->derv_bufsize_step[subsc[idim]]=icvp->derv_dim_scale[idim]; if (icvp->derv_bufsize_step[subsc[idim]] != 1) icvp->derv_do_bufsize_step = TRUE; } MI_RETURN(MI_NOERROR); } /* ----------------------------- MNI Header ----------------------------------- @NAME : MI_icv_get_dim_conversion @INPUT : icvp - icv pointer subsc - array of dimension subscripts for fastest varying dimensions @OUTPUT : (none) @RETURNS : MI_ERROR if an error occurs @DESCRIPTION: Sets the variables for dimensions converions @METHOD : @GLOBALS : @CALLS : NetCDF routines @CREATED : September 8, 1992 (Peter Neelin) @MODIFIED : ---------------------------------------------------------------------------- */ PRIVATE int MI_icv_get_dim_conversion(mi_icv_type *icvp, int subsc[]) /* ARGSUSED */ { int idim; MI_SAVE_ROUTINE_NAME("MI_icv_get_dim_conversion"); /* Find out whether we need to compress variable or user buffer */ icvp->derv_var_compress = (icvp->var_is_vector && icvp->user_do_scalar); icvp->derv_usr_compress = FALSE; for (idim=0; idimuser_num_imgdims; idim++) { if (icvp->derv_dim_scale[idim]!=1) { if (icvp->derv_dim_grow[idim]) icvp->derv_usr_compress = TRUE; else icvp->derv_var_compress = TRUE; } } /* Get the fastest varying dimension in user's buffer */ icvp->derv_dimconv_fastdim = icvp->var_ndims-1; if (icvp->var_is_vector && icvp->user_do_scalar) icvp->derv_dimconv_fastdim--; /* Find out how many pixels to compress/expand for variable and user buffers and allocate arrays */ if (icvp->var_is_vector && icvp->user_do_scalar) icvp->derv_var_pix_num=icvp->var_vector_size; else icvp->derv_var_pix_num=1; icvp->derv_usr_pix_num=1; for (idim=0; idimuser_num_imgdims; idim++) { if (icvp->derv_dim_grow[idim]) icvp->derv_usr_pix_num *= icvp->derv_dim_scale[idim]; else icvp->derv_var_pix_num *= MIN(icvp->var_dim_size[idim], icvp->derv_dim_scale[idim]); } icvp->derv_var_pix_off = MALLOC(icvp->derv_var_pix_num, long); icvp->derv_usr_pix_off = MALLOC(icvp->derv_usr_pix_num, long); if ((icvp->derv_var_pix_off == NULL) || (icvp->derv_usr_pix_off == NULL)) { MI_LOG_SYS_ERROR1("MI_icv_get_dim_conversion"); MI_RETURN_ERROR(MI_ERROR); } MI_RETURN(MI_NOERROR); } /* ----------------------------- MNI Header ----------------------------------- @NAME : MI_icv_dimconvert @INPUT : operation - MI_PRIV_GET or MI_PRIV_PUT icvp - icv structure pointer start - start passed by user count - count passed by user values - pointer to user's data area (for put) bufstart - start of variable buffer bufcount - count of variable buffer buffer - pointer to variable buffer (for get) @OUTPUT : values - pointer to user's data area (for get) buffer - pointer to variable buffer (for put) @RETURNS : MI_ERROR if an error occurs @DESCRIPTION: Converts values and dimensions from an input buffer to the user's buffer. Called by MI_var_action. @METHOD : @GLOBALS : @CALLS : NetCDF routines @CREATED : August 27, 1992 (Peter Neelin) @MODIFIED : ---------------------------------------------------------------------------- */ PRIVATE int MI_icv_dimconvert(int operation, mi_icv_type *icvp, long start[], long count[], void *values, long bufstart[], long bufcount[], void *buffer) { mi_icv_dimconv_type dim_conv_struct; mi_icv_dimconv_type *dcp; double sum0, sum1; /* Counters for averaging values */ double dvalue; /* Pixel value */ long counter[MAX_VAR_DIMS]; /* Dimension loop counter */ void *ptr, *iptr, *optr; /* Pointers for looping through fastest dim */ void *ivecptr[MAX_VAR_DIMS]; /* Pointers to start of each dimension */ void *ovecptr[MAX_VAR_DIMS]; long *end; /* Pointer to array of dimension ends */ int fastdim; /* Dimension that varies fastest */ long ipix; /* Buffer subscript */ int idim; /* Dimension subscript */ int notmodified; /* First dimension not reset */ int out_of_range; /* Flag indicating one pixel of sum out of range */ double dmin, dmax, epsilon; /* Range limits */ MI_SAVE_ROUTINE_NAME("MI_icv_dimconvert"); /* Initialize variables */ dcp = &dim_conv_struct; {MI_CHK_ERR(MI_icv_dimconv_init(operation, icvp, dcp, start, count, values, bufstart, bufcount, buffer))} /* Initialize local variables */ iptr = dcp->istart; optr = dcp->ostart; end = dcp->end; fastdim = icvp->derv_dimconv_fastdim; dmax = icvp->fill_valid_max; dmin = icvp->fill_valid_min; epsilon = (dmax - dmin) * FILLVALUE_EPSILON; epsilon = fabs(epsilon); dmax += epsilon; dmin -= epsilon; /* Initialize counters */ for (idim=0; idim<=fastdim; idim++) { counter[idim] = 0; ivecptr[idim] = iptr; ovecptr[idim] = optr; } /* Loop through data */ while (counter[0] < end[0]) { /* Compress data by averaging if needed */ if (!dcp->do_compress) { {MI_TO_DOUBLE(dvalue, dcp->intype, dcp->insign, iptr)} out_of_range = (icvp->do_fillvalue && ((dvalue < dmin) || (dvalue > dmax))); } else { sum1 = 0.0; sum0 = 0.0; out_of_range=FALSE; for (ipix=0; ipixin_pix_num; ipix++) { ptr=(void *) ((char *)iptr + dcp->in_pix_off[ipix]); /* Check if we are outside the buffer. If we are looking before the buffer, then we need to add in the previous result to do averaging properly. If we are looking after the buffer, then break. */ if (ptrin_pix_first) { /* Get the output value and re-scale it */ {MI_TO_DOUBLE(dvalue, dcp->outtype, dcp->outsign, optr)} if (icvp->do_scale) { dvalue = ((icvp->scale==0.0) ? 0.0 : (dvalue - icvp->offset) / icvp->scale); } } else if (ptr>dcp->in_pix_last) { continue; } else { {MI_TO_DOUBLE(dvalue, dcp->intype, dcp->insign, ptr)} } /* Add in the value, checking for range if needed */ if (icvp->do_fillvalue && ((dvalue < dmin) || (dvalue > dmax))) { out_of_range = TRUE; } else { sum1 += dvalue; sum0++; } } /* Foreach pixel to compress */ /* Average values */ if (sum0!=0.0) dvalue = sum1/sum0; else dvalue = 0.0; } /* If compress */ /* Check for out of range values and scale result */ if (out_of_range) { dvalue = icvp->user_fillvalue; } else if (icvp->do_scale) { dvalue = icvp->scale * dvalue + icvp->offset; } /* Expand data if needed */ if (!dcp->do_expand) { {MI_FROM_DOUBLE(dvalue, dcp->outtype, dcp->outsign, optr)} } else { for (ipix=0; ipixout_pix_num; ipix++) { ptr=(void *) ((char *)optr + dcp->out_pix_off[ipix]); /* Check if we are outside the buffer. */ if ((ptr>=dcp->out_pix_first) && (ptr<=dcp->out_pix_last)) { {MI_FROM_DOUBLE(dvalue, dcp->outtype, dcp->outsign, ptr)} } } /* Foreach pixel to expand */ } /* if expand */ /* Increment the counter and the pointers */ if ((++counter[fastdim]) < end[fastdim]) { optr = (void *) ((char *) optr + dcp->ostep[fastdim]); iptr = (void *) ((char *) iptr + dcp->istep[fastdim]); } else { /* If we reach the end of fastdim, then reset the counter and increment the next dimension down - keep going as needed. The vectors ovecptr and ivecptr give the starting values of optr and iptr for that dimension. */ idim = fastdim; while ((idim>0) && (counter[idim] >= end[idim])) { counter[idim] = 0; idim--; counter[idim]++; ovecptr[idim] = (void *)((char *)ovecptr[idim]+dcp->ostep[idim]); ivecptr[idim] = (void *)((char *)ivecptr[idim]+dcp->istep[idim]); } notmodified = idim; /* Copy the starting index up the vector */ for (idim=notmodified+1; idim<=fastdim; idim++) { ovecptr[idim]=ovecptr[notmodified]; ivecptr[idim]=ivecptr[notmodified]; } optr = ovecptr[fastdim]; iptr = ivecptr[fastdim]; } /* if at end of row */ } /* while more pixels to process */ MI_RETURN(MI_NOERROR); } /* ----------------------------- MNI Header ----------------------------------- @NAME : MI_icv_dimconv_init @INPUT : operation - MI_PRIV_GET or MI_PRIV_PUT icvp - icv structure pointer dcp - dimconvert structure pointer start - start passed by user count - count passed by user values - pointer to user's data area (for put) bufstart - start of variable buffer bufcount - count of variable buffer buffer - pointer to variable buffer (for get) @OUTPUT : @RETURNS : MI_ERROR if an error occurs @DESCRIPTION: Sets up stuff for MI_icv_dimconvert. @METHOD : @GLOBALS : @CALLS : NetCDF routines @CREATED : September 4, 1992 (Peter Neelin) @MODIFIED : ---------------------------------------------------------------------------- */ PRIVATE int MI_icv_dimconv_init(int operation, mi_icv_type *icvp, mi_icv_dimconv_type *dcp, long start[], long count[], void *values, long bufstart[], long bufcount[], void *buffer) /* ARGSUSED */ { long buffer_len, values_len; /* Buffer lengths, offsets and indices */ long buffer_off, values_off; long buffer_index, values_index; int imgdim_high, imgdim_low; /* Range of subscripts of image dimensions */ int scale, offset; /* Dimension scale and offset */ int idim, jdim; int fastdim; /* Variables for calculating pixel offsets for compress/expand */ long var_dcount[MI_MAX_IMGDIMS+1], var_dend[MI_MAX_IMGDIMS+1]; long usr_dcount[MI_MAX_IMGDIMS+1], usr_dend[MI_MAX_IMGDIMS+1]; long pixcount; int var_fd, usr_fd, dshift; long ipix; MI_SAVE_ROUTINE_NAME("MI_icv_dimconv_init"); /* Check to see if any compression or expansion needs to be done. Work it out for a GET and then swap if a PUT. */ if (operation==MI_PRIV_GET) { dcp->do_compress = icvp->derv_var_compress; dcp->do_expand = icvp->derv_usr_compress; } else { dcp->do_expand = icvp->derv_var_compress; dcp->do_compress = icvp->derv_usr_compress; } fastdim = icvp->derv_dimconv_fastdim; /* Get the indices of high and low image dimensions */ imgdim_high=icvp->var_ndims-1; if (icvp->var_is_vector) imgdim_high--; imgdim_low = imgdim_high - icvp->user_num_imgdims + 1; /* Get the buffer sizes */ buffer_len = icvp->var_typelen; values_len = icvp->user_typelen; for (idim=0; idim < icvp->var_ndims; idim++) { buffer_len *= bufcount[idim]; if (idim<=fastdim) values_len *= icvp->derv_icv_count[idim]; } /* Calculate step size for variable and user buffers. This does not allow for growing or shrinking pixels. That correction is done below. */ if (icvp->var_is_vector && icvp->user_do_scalar) { dcp->buf_step[fastdim+1] = icvp->var_typelen; dcp->buf_step[fastdim] = dcp->buf_step[fastdim+1] * bufcount[fastdim+1]; } else { dcp->buf_step[fastdim] = icvp->var_typelen; } dcp->usr_step[fastdim] = icvp->user_typelen; for (idim=fastdim-1; idim>=0; idim--) { dcp->buf_step[idim] = dcp->buf_step[idim+1] * bufcount[idim+1]; dcp->usr_step[idim] = dcp->usr_step[idim+1] * icvp->derv_icv_count[idim+1]; } /* Set sign of user steps for flipping, if needed */ for (idim=imgdim_low; idim <= imgdim_high; idim++) { if (icvp->derv_dim_flip[imgdim_high-idim]) dcp->usr_step[idim] *= (-1); } /* Get the pointers to the start of buffers and the number of pixels in each dimension (count a pixel as one expansion/compression - one time through the loop below) */ buffer_off = 0; values_off = 0; for (idim=0; idim <= fastdim; idim++) { if ((idim < imgdim_low) || (idim > imgdim_high)) { dcp->end[idim] = bufcount[idim]; buffer_index = 0; values_index = bufstart[idim] - icvp->derv_icv_start[idim]; } else { jdim = imgdim_high - idim; scale = icvp->derv_dim_scale[jdim]; offset = icvp->derv_dim_off[jdim]; if (icvp->derv_dim_grow[jdim]) { dcp->end[idim] = bufcount[idim]; buffer_index = 0; if (!icvp->derv_dim_flip[jdim]) values_index = bufstart[idim]*scale - icvp->derv_icv_start[idim] + offset; else values_index = (icvp->var_dim_size[jdim] - bufstart[idim])*scale - 1 - icvp->derv_icv_start[idim] + offset; } else { dcp->end[idim] = (bufcount[idim] - 1 + bufstart[idim]%scale) / scale + 1; buffer_index = -(bufstart[idim] % scale); if (!icvp->derv_dim_flip[jdim]) values_index = bufstart[idim]/scale - icvp->derv_icv_start[idim] + offset; else values_index = (icvp->var_dim_size[jdim] - bufstart[idim] - 1)/scale - icvp->derv_icv_start[idim] + offset; } } buffer_off += buffer_index * labs(dcp->buf_step[idim]); values_off += values_index * labs(dcp->usr_step[idim]); } /* Calculate arrays of offsets for compress/expand. */ if (dcp->do_compress || dcp->do_expand) { /* Initialize counters */ var_fd = icvp->user_num_imgdims-1; usr_fd = icvp->user_num_imgdims-1; if (icvp->var_is_vector && icvp->user_do_scalar) { var_fd++; var_dcount[var_fd]=0; var_dend[var_fd]=icvp->var_vector_size; } for (jdim=0; jdimuser_num_imgdims; jdim++) { idim=icvp->user_num_imgdims - jdim - 1; var_dcount[idim] = 0; usr_dcount[idim] = 0; var_dend[idim] = (icvp->derv_dim_grow[jdim] ? 1 : MIN(icvp->var_dim_size[jdim], icvp->derv_dim_scale[jdim])); usr_dend[idim] = (icvp->derv_dim_grow[jdim] ? icvp->derv_dim_scale[jdim] : 1); } /* Loop through variable buffer pixels */ pixcount=0; dshift = imgdim_low; for (ipix=0; ipixderv_var_pix_num; ipix++) { icvp->derv_var_pix_off[ipix] = pixcount; pixcount += dcp->buf_step[var_fd+dshift]; if ((++var_dcount[var_fd]) >= var_dend[var_fd]) { idim=var_fd; while ((idim>0) && (var_dcount[idim]>=var_dend[idim])) { var_dcount[idim]=0; idim--; var_dcount[idim]++; } for (idim=0, pixcount=0; idim<=var_fd; idim++) { pixcount += var_dcount[idim] * dcp->buf_step[idim+dshift]; } } } /* Loop through user buffer pixels */ pixcount=0; dshift = imgdim_low; for (ipix=0; ipixderv_usr_pix_num; ipix++) { icvp->derv_usr_pix_off[ipix] = pixcount; pixcount += dcp->usr_step[usr_fd+dshift]; if ((++usr_dcount[usr_fd]) >= usr_dend[usr_fd]) { idim=usr_fd; while ((idim>0) && (usr_dcount[idim]>=usr_dend[idim])) { usr_dcount[idim]=0; idim--; usr_dcount[idim]++; } for (idim=0, pixcount=0; idim<=var_fd; idim++) { pixcount += usr_dcount[idim] * dcp->usr_step[idim+dshift]; } } } /* Correct buffer steps for compress/expand */ for (idim=imgdim_low; idim <= imgdim_high; idim++) { jdim = imgdim_high-idim; if (icvp->derv_dim_grow[jdim]) dcp->usr_step[idim] *= icvp->derv_dim_scale[jdim]; else dcp->buf_step[idim] *= icvp->derv_dim_scale[jdim]; } } /* if compress/expand */ /* Set input and output variables */ if (operation==MI_PRIV_GET) { /* For a GET */ dcp->in_pix_num = icvp->derv_var_pix_num; dcp->in_pix_off = icvp->derv_var_pix_off; dcp->in_pix_first = buffer; dcp->in_pix_last = (void *) ((char *)buffer + buffer_len - 1); dcp->out_pix_num = icvp->derv_usr_pix_num; dcp->out_pix_off = icvp->derv_usr_pix_off; dcp->out_pix_first = values; dcp->out_pix_last = (void *) ((char *)values + values_len - 1); dcp->intype = icvp->var_type; dcp->insign = icvp->var_sign; dcp->outtype = icvp->user_type; dcp->outsign = icvp->user_sign; dcp->istep = dcp->buf_step; dcp->ostep = dcp->usr_step; dcp->istart = (void *) ((char *) buffer + buffer_off); dcp->ostart = (void *) ((char *) values + values_off); } /* if GET */ else { /* For a PUT */ dcp->out_pix_num = icvp->derv_var_pix_num; dcp->out_pix_off = icvp->derv_var_pix_off; dcp->out_pix_first = buffer; dcp->out_pix_last = (void *) ((char *)buffer + buffer_len - 1); dcp->in_pix_num = icvp->derv_usr_pix_num; dcp->in_pix_off = icvp->derv_usr_pix_off; dcp->in_pix_first = values; dcp->in_pix_last = (void *) ((char *)values + values_len - 1); dcp->outtype = icvp->var_type; dcp->outsign = icvp->var_sign; dcp->intype = icvp->user_type; dcp->insign = icvp->user_sign; dcp->ostep = dcp->buf_step; dcp->istep = dcp->usr_step; dcp->ostart = (void *) ((char *) buffer + buffer_off); dcp->istart = (void *) ((char *) values + values_off); } /* if PUT */ MI_RETURN(MI_NOERROR); }