home *** CD-ROM | disk | FTP | other *** search
/ Photo CD Demo 1 / Demo.bin / hdf / unix / hdf3_2r2 / src / dfufp2i.c < prev    next >
Encoding:
C/C++ Source or Header  |  1992-10-29  |  26.7 KB  |  800 lines

  1. /***************************************************************************
  2. *
  3. *
  4. *                         NCSA HDF version 3.2r2
  5. *                            October 30, 1992
  6. *
  7. * NCSA HDF Version 3.2 source code and documentation are in the public
  8. * domain.  Specifically, we give to the public domain all rights for future
  9. * licensing of the source code, all resale rights, and all publishing rights.
  10. *
  11. * We ask, but do not require, that the following message be included in all
  12. * derived works:
  13. *
  14. * Portions developed at the National Center for Supercomputing Applications at
  15. * the University of Illinois at Urbana-Champaign, in collaboration with the
  16. * Information Technology Institute of Singapore.
  17. *
  18. * THE UNIVERSITY OF ILLINOIS GIVES NO WARRANTY, EXPRESSED OR IMPLIED, FOR THE
  19. * SOFTWARE AND/OR DOCUMENTATION PROVIDED, INCLUDING, WITHOUT LIMITATION,
  20. * WARRANTY OF MERCHANTABILITY AND WARRANTY OF FITNESS FOR A PARTICULAR PURPOSE
  21. *
  22. ****************************************************************************
  23. */
  24.  
  25. #ifdef RCSID
  26. static char RcsId[] = "@(#)$Revision: 1.2 $";
  27. #endif
  28. /*
  29. $Header: /hdf/hdf/v3.2r2/src/RCS/dfufp2i.c,v 1.2 1992/10/23 00:14:11 koziol beta koziol $
  30.  
  31. $Log: dfufp2i.c,v $
  32.  * Revision 1.2  1992/10/23  00:14:11  koziol
  33.  * Changed all DFIstr*() and DFImem*() calls to HDstr*() and HDmem*() calls
  34.  * #ifdef'd out the macros Jason defined for Hopen, Hclose, etc. for Vsets
  35.  * Replaced Vset VFREESPACE and VGETSPACE calls with actual calls to HDfreespace
  36.  * and HDgetspace
  37.  * Added a MS-Windows lower lower for file I/O (which may not be completely working
  38.  *
  39.  * Revision 1.1  1992/08/25  21:40:44  koziol
  40.  * Initial revision
  41.  *
  42. */
  43. /*-----------------------------------------------------------------------------c
  44.  *    dfufp2i.c
  45.  *
  46.  * Purpose:  Utility function to convert floating point data to 8-bit 
  47.  *           raster image set (RIS8) format, storing the results in 
  48.  *           an hdf file.
  49.  *                                   -----------
  50.  *                                  |           | ----------> RIS8
  51.  *       floating point data   ---> | fp_to_hdf |   and/or
  52.  *         (in an array)            |           | ----------> SDS
  53.  *                                   -----------
  54.  *
  55.  * Invokes:  libdf.a
  56.  *
  57.  * Includes: stdio.h, ctype.h, "df.h"
  58.  *
  59.  * Public function:
  60.  *      DFUfptoimage: sets up structs with input params, calls process()
  61.  * 
  62.  * Private functions:
  63.  *      process: main driver routine: transforms the data to an image and 
  64.  *               stores it in the file
  65.  *      generate_scale: generates a scale, if none provided
  66.  *      convert_interp: creates an interpolated image
  67.  *      pixrep_scaled: creates an expanded image using scales provided
  68.  *      compute_offsets: called by pixrep_scaled 
  69.  *      pixrep_simple: creates an expanded image assuming equal gaps in scales
  70.  *
  71.  * Fortran stub function:
  72.  *      duif2i - intermediate, called by fortran functions in DFUfptoimage.f
  73.  *
  74.  * Remarks:
  75.  *      This routine is very similar to the utility fp_to_hdf, which
  76.  *      takes its input from one or more files, rather than from internal
  77.  *      memory. 
  78.  *      Another difference is that this routine allows compression (run
  79.  *      length encoding), whereas fp_to_hdf does not at present (8/31/89).
  80.  *      Since this routine is meant to mimic many of the features of
  81.  *      NCSA DataScope, much of the code has been taken directly from
  82.  *      the DataScope source.
  83.  *
  84.  *  National Center for Supercomputing Applications
  85.  *  University of Illinois, Urbana-Champaign
  86.  *
  87.  *  by Mike Folk (mfolk@ncsa.uiuc.edu)
  88.  *  Beta version: 9/1/89
  89.  *  Released:     6/5/90
  90.  *
  91.  *  This program is in the public domain
  92.  *
  93.  *--------------------------------------------------------------------------*/
  94.  
  95. #include <ctype.h>
  96. #include "hdf.h"
  97. #include "dfufp2i.h"
  98. #include "dfsd.h"
  99.  
  100. /**********************************************************************
  101. *
  102. *  Header information
  103. *
  104. ***********************************************************************/
  105.  
  106. /*
  107. *  global definitions
  108. */
  109.  
  110. #ifndef FALSE
  111. #define FALSE 0
  112. #endif
  113.  
  114. #ifndef TRUE
  115. #define TRUE    1
  116. #endif /* TRUE */
  117. #define EXPAND  1   /* -e: expand image with pixel replication */
  118. #define INTERP  2   /* -i: expand image with bilinear interpolation */
  119.  
  120.  
  121. #ifndef DFUFP2I_FNAMES
  122. #   define DFUFP2I_FNAMES
  123. #ifdef DF_CAPFNAMES
  124. #   define nduif2i       FNAME(DUIF2I)
  125. #else
  126. #   define nduif2i       FNAME(duif2i)
  127. #endif /* DF_CAPFNAMES */
  128. #endif /* DFUFP2I_FNAMES */
  129.  
  130. /*-----------------------------------------------------------------------------
  131.  * Name:    duif2i
  132.  * Purpose: Intermediate Fortran callable version of DFUfptoimage
  133.  *          (See DFUfptoimage for details)
  134.  *
  135.  * Invokes: DFUfptoimage
  136.  *---------------------------------------------------------------------------*/
  137.  
  138.     FRETVAL(int)
  139. #ifdef PROTOTYPE
  140. nduif2i(int32 *hdim, int32 *vdim, float32 *max, float32 *min, float32 hscale[],
  141.     float32 vscale[], float32 data[], uint8 *palette, _fcd outfile,
  142.     int *ct_method, int32 *hres, int32 *vres, int *compress, int *lenfn)
  143. #else
  144. nduif2i(hdim,vdim,max,min,hscale,vscale,data,palette,
  145.                                    outfile,ct_method,hres,vres,compress,lenfn)
  146.  
  147.     int32 *hdim, *vdim; /* horizontal and vertical dimensions of input data */
  148.     float32 *max, *min,/* maximum and minimum values in the data */
  149.             hscale[]   /* horizontal and vertical scales */,
  150.             vscale[],
  151.             data[];    /* input data */
  152.     uint8 *palette;    /* palette to be stored with the image */
  153.     _fcd outfile;      /* name of file to store image in */
  154.     int  *ct_method;   /* color transform method: 1=EXPAND; 2=INTERP */
  155.     int32 *hres, *vres; /* resolutions desired for output image      */
  156.     int *compress,      /* flag: 0 = don't compress; 1=do compression */
  157.          *lenfn;        /* length of outfile string */
  158. #endif /* PROTOTYPE */
  159. {
  160.     char *fn;
  161.     int  ret;
  162.     
  163.     fn = HDf2cstring(outfile, *lenfn);
  164.     ret = DFUfptoimage( *hdim, *vdim, *max, *min, hscale, vscale, data,
  165.             _fcdtocp(palette), fn, *ct_method, *hres, *vres, *compress );
  166.     HDfreespace(fn);
  167.     return(ret);
  168. }
  169.  
  170.  
  171. /*-----------------------------------------------------------------------------s
  172.  * DFUfptoimage()
  173.  *  
  174.  * Purpose:sets up structs with input params, calls process()
  175.  * Inputs:
  176.  *     hdim, vdim: horizontal and vertical dimensions of input data
  177.  *     max, min:   maximum and minimum values in the data
  178.  *     hscale,vscale: optional horizontal and vertical scales
  179.  *     data:       input data
  180.  *     palette:    optional palette to be stored with the image
  181.  *     outfile:n   name of hdf file to store image in 
  182.  *     ct_method:  color transform method: 1=EXPAND; 2=INTERP
  183.  *     hres, vres: resolutions desired for output image
  184.  *     compress:   compression flag: 0=don't; 1=do
  185.  * Returns: 0 on success, -1 on failure with error set
  186.  * Users:       HDF HLL (high-level library) users, utilities, other routines
  187.  * Invokes: process
  188.  * Remarks: none
  189.  *----------------------------------------------------------------------------*/
  190.  
  191. #ifdef PROTOTYPE
  192. int DFUfptoimage(int32 hdim, int32 vdim, float32 max, float32 min,
  193.         float32 *hscale, float32 *vscale, float32 *data, uint8 *palette,
  194.         char *outfile, int ct_method, int32 hres, int32 vres, int compress)
  195. #else
  196. int DFUfptoimage(hdim,vdim,max,min,hscale,vscale,data,palette,
  197.                                    outfile,ct_method,hres,vres,compress)
  198.  
  199. int32   hdim, vdim,      /* horizontal and vertical dimensions of input data */
  200.         hres, vres;      /* resolutions desired for output image      */
  201. float32 max, min,        /* maximum and minimum values in the data */
  202.         *hscale,*vscale, /* horizontal and vertical scales */
  203.         *data;           /* input data */
  204. uint8   *palette;        /* palette to be stored with the image */
  205. char    *outfile;        /* name of file to store image in */
  206. int     ct_method,       /* color transform method: 1=EXPAND; 2=INTERP */
  207.         compress;        /* flag: 0 = don't compress; 1=do compression */
  208. #endif /* PROTOTYPE */
  209. {
  210.     struct Input in;
  211.     struct Output out;
  212.  
  213.     in.hdim = hdim;
  214.     in.vdim = vdim;
  215.     in.max = max;
  216.     in.min = min;
  217.     in.is_hscale = (hscale == NULL) ? FALSE : TRUE;
  218.     in.is_vscale = (vscale == NULL) ? FALSE : TRUE;
  219.     in.hscale = hscale;
  220.     in.vscale = vscale;
  221.     in.data = data;
  222.     in.is_pal = (palette == NULL) ? FALSE : TRUE;
  223.     in.ct_method = ct_method;
  224.     HDstrcpy(out.outfile, outfile);  /* get outfile name */
  225.     out.palette = palette;         /* get palette address (may be NULL) */
  226.     out.hres = hres;
  227.     out.vres = vres;
  228.     out.compress = compress ? 11 : 0; /* 0=>don't; 11=>RLE compression */
  229.  
  230.     /* tloc1 = time((long *) 0); */  /* these 4 lines for debugging */
  231.     process(&in,&out);
  232.     /* tloc2 = time((long *) 0); */
  233.     /* printf("Time:    %ld\n",tloc2-tloc1); */
  234.     return 0;
  235. } /* end of DFUfptoimage */
  236.  
  237. /*-----------------------------------------------------------------------------s
  238.  * process
  239.  *
  240.  * Purpose:   to transform the data to an image and stores it in the file
  241.  * Inputs:
  242.  *     in:   structure with information about data to be converted to image
  243.  *     out:  structure with information about image
  244.  * Returns: 0 on success, -1 on failure with error set
  245.  * Users:    DFUfptoimage
  246.  * Invokes: from libdf.a: DFR8setpalette, Hopen, Hclose, DFR8addimage
  247.  *          local: generate_scale, pixrep_scaled, pixrep_simple,convert_interp 
  248.  * Remarks: none
  249.  *----------------------------------------------------------------------------*/
  250.  
  251. #ifdef PROTOTYPE
  252. int process(struct Input *in, struct Output *out)
  253. #else
  254. int process(in, out)
  255. struct Input *in;
  256. struct Output *out;
  257. #endif /* PROTOTYPE */
  258. {
  259.     int ret;
  260.     int32   file_id;
  261.         
  262. /*    printinput(in);*/   /* for debugging */
  263.  
  264.     if (in->is_pal) {
  265.         ret = DFR8setpalette((uint8*)out->palette);  /* output as HDF palette */
  266.         if (ret < 0) return ret;
  267.     }
  268.  
  269.     file_id = Hopen(out->outfile, DFACC_WRITE, 0);
  270.     Hclose(file_id);
  271.  
  272.     /*
  273.     *  allocate buffers for output and scales 
  274.     */
  275.     if (!in->is_hscale)
  276.         in->hscale = (float32 *) HDgetspace((uint32)(1+in->hdim)*sizeof(float32));
  277.     if (!in->is_vscale)
  278.         in->vscale = (float32 *) HDgetspace((uint32)(1+in->vdim)*sizeof(float32));
  279.     out->hres = (out->hres <= in->hdim) ? in->hdim : out->hres ;
  280.     out->vres = (out->vres <= in->vdim) ? in->vdim : out->vres ;
  281.     out->image = (uint8 *) HDgetspace((uint32)out->hres * out->vres);
  282.  
  283.     /*
  284.     *  if necessary, generate x and y scales
  285.     */
  286.     if (!in->is_hscale)
  287.         generate_scale(in->hdim, in->hscale);
  288.     if (!in->is_vscale)
  289.         generate_scale(in->vdim, in->vscale);
  290.     /*
  291.     *  output raster hdf file 
  292.     */
  293.     if (in->ct_method == EXPAND) {
  294.         if (in->is_hscale || in->is_vscale)
  295.             pixrep_scaled(in, out);
  296.         else
  297.             pixrep_simple(in, out);
  298.     } else
  299.         convert_interp(in, out);
  300.  
  301. /*    printoutput(out);*/   /* for debugging */
  302.  
  303.     ret = DFR8addimage(out->outfile, (char *)out->image,
  304.                              out->hres,out->vres,out->compress);
  305.     if (ret < 0) return ret;
  306.     /*
  307.     *  free allocated space
  308.     */
  309.     if (!in->is_hscale) HDfreespace((char *)in->hscale);
  310.     if (!in->is_vscale) HDfreespace((char *)in->vscale);
  311.     HDfreespace((char *)out->image);
  312.     return 0;
  313. } /* end of process */
  314.  
  315.  
  316. /*-----------------------------------------------------------------------------
  317.  * generate_scale
  318.  *
  319.  * Purpose:   to generate the scale 1 2 3 ... dim
  320.  * Input: 
  321.  *     dim:   length of scale 
  322.  * Output:
  323.  *     scale: array of floating point numbers from 1 to dim
  324.  * Returns: 0 on success, -1 on failure with error set
  325.  * Users:    process
  326.  * Invokes: none
  327.  * Remarks: none
  328.  *---------------------------------------------------------------------------*/
  329.  
  330. #ifdef PROTOTYPE
  331. int generate_scale(int32 dim, float32 *scale)
  332. #else
  333. int generate_scale(dim, scale)
  334. int32 dim;
  335. float32 *scale;
  336. #endif /* PROTOTYPE */
  337. {
  338.     int32 i;
  339.  
  340.     for (i=0; i <= dim; i++)
  341.         *scale++ = (float32) i;
  342.     return 0;
  343. }
  344.  
  345.  
  346. /*-----------------------------------------------------------------------------
  347.  * printinput
  348.  *
  349.  * Purpose:   debugging: prints input values to stdout
  350.  * Input: 
  351.  *     in:    struct with all input values
  352.  * Returns: 0 on success, -1 on failure with error set
  353.  * Users:    process and other local routines
  354.  * Invokes: none
  355.  * Remarks: none
  356.  *---------------------------------------------------------------------------*/
  357. /*  This function is commented out of the code!! */
  358. #ifdef DEBUG_HDF
  359. #ifdef PROTOTYPE
  360. int printinput(struct Input *in)
  361. #else
  362. int printinput(in)
  363. struct Input *in;
  364. #endif /* PROTOTYPE */
  365. {
  366.     int i,j; 
  367.  
  368.     printf("\nmax: %8.2f   min: %8.2f\n", in->max, in->min);
  369.  
  370.     printf("There %s a palette\n", in->is_pal ? "IS" : "is NOT");
  371.     printf("color tranform method = %s\n",
  372.                            (in->ct_method==EXPAND)?"expand":"interpolate");
  373.     if ( in->hscale != NULL) {
  374.         printf("\nHorizontal scale:\n");
  375.         for (i=0; i<(int)(in->hdim); i++)
  376.           printf("%8.2f",in->hscale[i]);
  377.     }
  378.     else
  379.         printf("\nNo horizontal scale\n");
  380.  
  381.     if ( in->vscale != NULL) {
  382.         printf("\nVertical scale:\n");
  383.         for (i=0; i<(int)(in->vdim); i++)
  384.           printf("%8.2f",in->vscale[i]);
  385.     }
  386.     else
  387.         printf("\nNo vertical scale.\n");
  388.  
  389.     printf("\n");
  390.     printf("Data:");
  391.     for(i=0;i<(int)(in->vdim) && i<11; i++) {
  392.         printf("\n");
  393.         for (j=0; j< (int)(in->hdim); j++)
  394.             printf("%6.1f ",in->data[i*in->hdim+j]);
  395.     }
  396.     printf("\n");
  397.     return 0;
  398. } /* end of print_input */
  399. #endif  /* DEBUG_HDF */
  400.  
  401. /*-----------------------------------------------------------------------------
  402.  * printoutput
  403.  *
  404.  * Purpose:   debugging: prints input values to stdout
  405.  * Input: 
  406.  *     out:    struct with all output values
  407.  * Returns: 0 on success, -1 on failure with error set
  408.  * Users:    process and other local routines
  409.  * Invokes: none
  410.  * Remarks: none
  411.  *---------------------------------------------------------------------------*/
  412. /*  This function is commented out of the code!! */
  413. #ifdef DEBUG_HDF
  414. #ifdef PROTOTYPE
  415. int printoutput(struct Output *out)
  416. #else
  417. int printoutput(out)
  418. struct Output *out;
  419. #endif /* PROTOTYPE */
  420. {
  421.     int i,j;
  422.     printf("\n");
  423.       for(i=0;i<(int)(out->vres) && i<20; i++) {
  424.           printf("\n");
  425.           for (j=0; j< (int)(out->hres); j++)
  426.              if (j<19)
  427.                 printf("%4d",(uint8) out->image[i*out->hres + j]);
  428.       }
  429.      printf("\n");
  430.     return 0;
  431. } /* end of printoutput */
  432. #endif  /* DEBUG_HDF */
  433.  
  434.  
  435.  
  436. /***************************************************************************
  437. *
  438. *  Next comes the routine for performing bilinear interpolation
  439. *
  440. ****************************************************************************/
  441.  
  442.  
  443. /*-----------------------------------------------------------------------------
  444.  * convert_interp
  445.  *
  446.  * Purpose:   Create an interpolated image from the data array
  447.  * Input: 
  448.  *     in:    struct with all input values
  449.  *     out:   struct with all output values
  450.  * Returns: 0 on success, -1 on failure with error set
  451.  * Users:   process
  452.  * Invokes: none
  453.  * Remarks: Uses a bilinear interpolation method to fill in the picture.
  454.  *---------------------------------------------------------------------------*/
  455.  
  456. #ifdef PROTOTYPE
  457. int convert_interp(struct Input *in, struct Output *out)
  458. #else
  459. int convert_interp(in, out)
  460. struct Input   *in;
  461. struct Output  *out;
  462. #endif /* PROTOTYPE */
  463. {
  464.     register int j,theval;
  465.     float32 *f,*dxs,*dys,*xv,*yv,*lim,delx,dely,pt,xrange,yrange,range,zy,
  466.         *z1,*z2,*z3,*z4,z;
  467.     register uint8 *p;
  468.     uint8 *xinc;
  469.     int32 i,*yoffs;
  470.  
  471.     p = (uint8 *) out->image; /* space for interpolated image */
  472.  
  473.     range = in->max - in->min;
  474.     xrange = *(in->hscale + in->hdim-1) - *in->hscale;
  475.     yrange = *(in->vscale + in->vdim-1) - *in->vscale;
  476.     delx = xrange / out->hres;              /* x axis increment in image */
  477.     dely = yrange / out->vres;              /* y axis increment in image */
  478.     
  479.     dxs = (float32 *)HDgetspace((uint32)sizeof(float32)*out->hres);
  480.                                             /* temp space for dx's */
  481.     dys = (float32 *)HDgetspace((uint32)sizeof(float32)*out->vres);
  482.                                             /* temp space for dy's */
  483.     xinc = (uint8 *)HDgetspace((uint32)out->hres);
  484.     yoffs = (int32 *)HDgetspace((uint32)(out->vres+1)*sizeof(int32));
  485.     yoffs[0] = 0;
  486.     
  487.     if (range < 0) range = -range;          /* max must be > min */
  488.     
  489.     f = dys;                                /* beginning of dys to fill in */
  490.     yv = in->vscale;                        /* beginning and end of yvals */
  491.     lim = in->vscale + in->vdim-2;
  492.     
  493.     if (yrange > 0) {
  494.         for (i=0; i<out->vres; i++) {        /* fill in dy's */
  495.             pt = dely*(float32)i + *in->vscale;/* scaled pos in new image */
  496.             
  497.             while (*(yv+1) < pt && yv < lim) {        /* move y pointer */
  498.                 yv++; yoffs[i]++;
  499.             }
  500.             *f++ = (*(yv+1) - pt)/(*(yv+1) - *yv);/* calc dy pcnt and put in */
  501.             
  502.             yoffs[i+1] = yoffs[i];
  503.         }
  504.     } 
  505.     else {                                        /* decrementing instead */
  506.         yrange = -yrange;    
  507.         
  508.         for (i=0; i<out->vres; i++) {        /* fill in dy's */
  509.             pt = dely*(float32)i + *in->vscale;
  510.  
  511.             while (*(yv+1) > pt && yv < lim) {    /* move y pointer */
  512.                 yv++; yoffs[i]++;
  513.             }
  514.             *f++ = -(*(yv+1) - pt)/(*(yv+1) - *yv);/* calc dy pcnt and put in */
  515.             
  516.             yoffs[i+1] = yoffs[i];
  517.         }
  518.     }
  519.  
  520.     f = dxs;                            /* beginning of dxs to fill in */
  521.     xv = in->hscale;                    /* beginning and end of xvals */
  522.     lim = in->hscale + in->hdim-2;
  523.     
  524.     if (xrange > 0) {
  525.         for (i=0; i<out->hres; i++) {        /* fill in dx's */
  526.             pt = delx*(float32)i + *in->hscale;
  527.             xinc[i] = 0;
  528.             
  529.             while (*(xv+1) < pt && xv < lim) {    /* move xv pointer */
  530.                 xv++;  xinc[i]++;
  531.             }
  532.             *f++ = (*(xv+1) - pt)/(*(xv+1) - *xv);/* calc dy prct and put in */
  533.         }
  534.     } 
  535.     else {                                        /* decrementing instead */
  536.         xrange = -xrange;    
  537.         
  538.         for (i=0; i<out->hres; i++) {        /* fill in dx's */
  539.             pt = delx*(float32)i + *in->hscale;
  540.             xinc[i] = 0;
  541.     
  542.             while (*(xv+1) > pt && xv < lim) {    /* move y pointer */
  543.                 xv++;  xinc[i]++;
  544.             }
  545.             *f++ = -(*(xv+1) - pt)/(*(xv+1) - *xv);/* calc dy pcnt and put in */
  546.         }
  547.     }
  548.  
  549. /*
  550. *   Do the interpolation for each point in the target image.
  551. *   We take advantage of the fact that we know the target is evenly spaced 
  552. *   along both axes.
  553. */
  554.     yv = dys;
  555.     
  556.     for (i= 0; i < out->vres; i++,yv++) {
  557.         
  558.         z1 = in->data + in->hdim*(yoffs[i]);
  559.         z2 = z1 + 1;
  560.         z3 = z1 + in->hdim;
  561.         z4 = z3 + 1;
  562.     
  563.         xv = dxs;
  564.         zy = *yv;
  565.         
  566.         for (j=0; j < (int)(out->hres); j++,xv++) {    /* for each target point */
  567.             
  568.             z1 += xinc[j];         /* xinc == 0 when we don't need to shift */
  569.             z2 += xinc[j];
  570.             z3 += xinc[j];
  571.             z4 += xinc[j];
  572.             
  573.             z = (*z1 - *z3 - *z2 + *z4)*(*xv)*zy +        /* weighted sum */
  574.                 (*z3 - *z4)*(*xv) + (*z2 - *z4)*zy + *z4;
  575.                 
  576.             theval = (int)(1.0 + 237.9*(z - in->min)/range);/* scaled value  */
  577.             if (theval >= 240 || theval < 1)
  578.                 *p++ = 0;
  579.             else
  580.                 *p++ = (uint8)theval;
  581.         }
  582.     }
  583.     HDfreespace( (char *) dxs);
  584.     HDfreespace( (char *) dys);
  585.     HDfreespace( (char *) xinc);
  586.     HDfreespace( (char *) yoffs);
  587.     return 0;
  588. } /* end of convert_interp */
  589.  
  590.  
  591. /****************************************************************************
  592. *
  593. *  Next come the routines for pixel replication
  594. *
  595. *  Two routines to create expanded image via pixel replication
  596. *
  597. *  pixrep_scaled replicates pixels according to given scales
  598. *  pixrep_simple replicates the same number of pixels for each point
  599. *
  600. ******************************************************************************/
  601.  
  602. /*-----------------------------------------------------------------------------
  603.  * pixrep_scaled
  604.  *
  605.  * Purpose:   Create an expanded image from the data array
  606.  * Input: 
  607.  *     in:    struct with all input values
  608.  *     out:   struct with all output values
  609.  * Returns: 0 on success, -1 on failure with error set
  610.  * Users:   process
  611.  * Invokes: compute_offsets
  612.  * Remarks: Uses pixel replication to fill in the picture. Replicates
  613.  *          pixels according to in->vscale and in->hscale
  614.  *---------------------------------------------------------------------------*/
  615.  
  616. #ifdef  PROTOTYPE
  617. int pixrep_scaled(struct Input *in, struct Output *out)
  618. #else
  619. int pixrep_scaled(in, out)
  620. struct Input   *in;
  621. struct Output  *out;
  622. #endif /* PROTOTYPE */
  623. {
  624.     register int32 j;
  625.     float *data, range, ratio;
  626.     register uint8 *image, *prevrow;
  627.     uint8 *pixvals;
  628.     int32 i, theval, *hoffsets, *voffsets, prevoffset;
  629.  
  630.     data = in->data;                          /* space for data */
  631.     image = (uint8 *) out->image;  /* space for image */
  632.  
  633.     range = in->max - in->min;
  634.     if (range < 0) range = -range;            /* max must be > min */
  635.     
  636.     hoffsets = (int32 *)HDgetspace((uint32)(out->hres+1)*sizeof(int32));
  637.     voffsets = (int32 *)HDgetspace((uint32)(out->vres+1)*sizeof(int32));
  638.     pixvals = (uint8 *)HDgetspace((uint32)in->hdim+1);
  639.     compute_offsets(in->hscale,in->hdim,hoffsets,out->hres);
  640.     compute_offsets(in->vscale,in->vdim,voffsets,out->vres);
  641.     
  642.     prevoffset = voffsets[0] - 1;
  643.     ratio = (float32)237.9 / range;
  644.  
  645.     for (i=0; i<out->vres; i++) {   /* for each row, store pixel vals */
  646.  
  647.         if (voffsets[i] > prevoffset) {  /* if new data row, compute pix vals */
  648.  
  649.             for (j=0; j<in->hdim; j++) {  /* compute vals for each data point*/
  650.                 theval = (int) (1.5 + ratio*(*data++ - in->min));
  651.                 if (theval >= 240 || theval < 1) theval = 0;
  652.                 pixvals[j] = (uint8) theval;
  653.             }
  654.  
  655.             for (j=0; j<out->hres; j++)       /* put row of pix vals into */
  656.                 *image++ = pixvals[hoffsets[j]]; /* next row of image */
  657.         } 
  658.  
  659.         else {                         /* else repeating a previous row */
  660.             prevrow = image - out->hres;
  661.             for (j=0; j<out->hres; j++)   /* put previous  row of pix vals*/
  662.                 *image++ = *prevrow++;       /* into next row of image */
  663.         }
  664.         prevoffset = voffsets[i];
  665.     }
  666.     HDfreespace( (char *) hoffsets);
  667.     HDfreespace( (char *) voffsets);
  668.     HDfreespace( (char *) pixvals);
  669.     return 0;
  670. } /* end of pixrep_scaled */
  671.             
  672.  
  673. /*-----------------------------------------------------------------------------
  674.  * compute_offsets
  675.  *
  676.  * Purpose:  For each pixel position on the horizontal or vertical 
  677.  *           dimension, compute the offet of the corresponding value 
  678.  *           in the scale array.
  679.  * Input: 
  680.  *     scale: the scale
  681.  *     dim:   length of scale
  682.  *     res:   resolution: length of the array 'offsets'
  683.  * Output:
  684.  *     offsets: the set of offsets that were computed
  685.  *     out:   struct with all output values
  686.  * Returns: 0 on success, -1 on failure with error set
  687.  * Users:   pixrep_scaled
  688.  * Invokes: compute_offsets
  689.  * Remarks: The array 'offsets' can be used to determine which scaled 
  690.  *          pixel value to place in the final picture.
  691.  *---------------------------------------------------------------------------*/
  692.  
  693. #ifdef PROTOTYPE
  694. int compute_offsets(float32 *scale, int32 dim, int32 *offsets, int32 res)
  695. #else
  696. int compute_offsets(scale, dim, offsets, res)
  697. float32 *scale;
  698. int32 dim,         /* dimension of scale */ 
  699.       res;         /* resolution: length of array 'offsets' */
  700. int32 *offsets;
  701. #endif /* PROTOTYPE */
  702. {
  703.     int32 i,j;
  704.     float32 *midpt, pt, delta;
  705.  
  706.     midpt = (float32 *) HDgetspace((uint32)sizeof(float32)*dim);
  707.  
  708.     for (i=0; i < dim-1; i++) {                  /* compute all midpoints */
  709.         midpt[i] = (scale[i] + scale[i+1])/2.0;
  710. /*        printf("midpt[%d]=%8.1f\tscale[%d]=%8.1f\n",i,midpt[i],i,scale[i]);*/
  711.     }
  712.     midpt[i] = scale[i] + scale[i] - midpt[i-1]; /* tack one onto end */
  713.  
  714.     delta = (*(scale+dim-1) - *scale)/(res-1); /* amt of change along scale */
  715.                                                /* per pixel position */
  716.     offsets[0] = 0;
  717.     pt = *scale;          /* base point has value of 1st scale item */
  718.  
  719.     for (i=1, j=0; i < res; i++) {  /* compute & store offsets of pix vals */
  720.         pt += delta;
  721.         offsets[i] = offsets[i-1];  /* keep offsets same until past midpt */
  722.         while (pt >= midpt[j]) {
  723.             offsets[i]++;
  724.             j++;
  725.         }
  726.     }
  727.     HDfreespace( (char *) midpt);
  728.     return 0;
  729. } /* end of compute_offsets */
  730.     
  731.  
  732. /*-----------------------------------------------------------------------------
  733.  * pixrep_simple
  734.  *
  735.  * Purpose:   Create an expanded image from the data array
  736.  * Input: 
  737.  *     in:    struct with all input values
  738.  *     out:   struct with all output values
  739.  * Returns: 0 on success, -1 on failure with error set
  740.  * Users:   process
  741.  * Invokes: compute_offsets
  742.  * Remarks: Uses pixel replication to fill in the picture. Replicates
  743.  *          the same number of pixels for each point
  744.  *---------------------------------------------------------------------------*/
  745.  
  746. #ifdef PROTOTYPE
  747. int pixrep_simple(struct Input *in, struct Output *out)
  748. #else
  749. int pixrep_simple(in, out)
  750. struct Input   *in;
  751. struct Output  *out;
  752. #endif /* PROTOTYPE */
  753. {
  754.     int32 i,j;
  755.     uint8 raster_val;
  756.     register uint8 *image, *row_buf;
  757.     float32 *in_row_ptr, *in_buf;
  758.     float32 ratio, delh, delv, hblockend, vblockend;
  759.  
  760.     ratio = (float32)237.9 / (in->max - in->min);
  761.     image = (uint8 *) out->image;
  762.     in_buf = in->data;
  763.  
  764.     delh = ( (float32) out->hres) / in->hdim;  /* horiz block size */
  765.     delv = ( (float32) out->vres) / in->vdim;  /* vert block size  */
  766.  
  767. /*
  768. * Compute expanded image
  769. * Do it a vertical block at a time
  770. * (Note the trick with the counters i and j vis-a-vis the blockends.)
  771. */
  772.     vblockend = delv;
  773.     for(i=0; i<out->vres; i++, vblockend += delv) { 
  774.         in_row_ptr = in_buf;          
  775.         row_buf = image;              /* start of next NEW row of output */
  776.  
  777.         /* compute raster values for this row */
  778.         hblockend = delh;
  779.         for(j=0; j<out->hres; j++, hblockend += delh) { 
  780.  
  781.             raster_val = (uint8)
  782.                              (1.5 + ratio * (float32) (*in_row_ptr++-in->min));
  783.             *image++ = raster_val;
  784.  
  785.             for (; j<(int32)hblockend-1; j++) /* store vals for this blk of this row*/
  786.                 *image++ = raster_val;
  787.         }
  788.  
  789.         /* repeat same row for whole vertical block*/
  790.         for(; i<(int32)vblockend-1; i++)
  791.             for (j=0; j<out->hres; j++)
  792.                 *image++ = row_buf[j];
  793.          
  794.         in_buf += in->hdim;            /* move to next row in input array */
  795.     }
  796.     return 0;
  797. } /* end of pixrep_simple() */
  798.  
  799.