render_z.c

Go to the documentation of this file.

/* Converts point precipitation data in struct pdata to hrap grid.  The
   20 closest stations are precalculated for each grid point - this
   speeds computations for data sets with many precipitation points.
   If there are no good precipitation points for a grid point, then
   a recalculation is made using all precipitation points. 1/R**2
   interpolation is used.  If requested, the final interpolation
   is scaled using seasonal isohyets. The grid is saved as a disk
   file and used later for a raster or vector (HRAP) plot. */
 
#include "prototypes_new.h"
 
void render_z(int pcpn_day,int pcpn_time,int pcpn_time_step)
 
{
 
extern int isom;
extern Widget top_level;
extern int pcp_in_use[500];
extern Cursor watch_cursor;
extern int method;
extern struct zdata zdata[10];
extern struct station zstation[500];
int i,j,h,hh,time_pos,htotal;
double distance,dist1,dist2,dist,value;
double temp;
extern struct hrap_grid *hrap_grid;
extern struct pcp *pcp;
extern Window window;
extern Display *display;
extern int max_zstations;
int totals[5];
 
for(i=0;i<5;i++)
        totals[i]=0;
 
XDefineCursor(display,window,watch_cursor); 
XmUpdateDisplay(top_level);
 
time_pos=pcpn_time;
 
for(i=0;i<hrap_grid->maxi;i++) {
 
for(j=0;j<hrap_grid->maxj;j++) {
 
           if(hrap_grid->owner[i][j]==-1) {
           
                                 pcp->value[i][j]=-9999;
                                 continue;
                                 
                                 }
 
           value=0.0;
           distance=0.0;
           htotal=0;
 
           for(h=0;h<5;h++) {
 
               hh=hrap_grid->gage[i][j].zindex[h];
 
               if(zdata[pcpn_day].stn[hh].zlevel2[time_pos].data < 0 )
                              continue;
 
               if(zdata[pcpn_day].stn[hh].zlevel2[time_pos].qual != 5 &&
                  zdata[pcpn_day].stn[hh].zlevel2[time_pos].qual != 8 &&
                  zdata[pcpn_day].stn[hh].zlevel2[time_pos].qual != 2 )
                          continue;
 
               dist1=(double)(hrap_grid->coord[i][j].x-
                      zstation[hh].x);    
               dist2=(double)(hrap_grid->coord[i][j].y-
                      zstation[hh].y);  
 
               dist=pow(dist1,2)+pow(dist2,2);
 
               if(dist < .00001)
                        dist=.00001;
               
               dist=1/dist;
 
               temp=zdata[pcpn_day].stn[hh].zlevel2[time_pos].data*dist;
 
               value=value+temp;
               
               distance=distance+dist;                    
 
               htotal++;
               
         }
                
 
           if(htotal < 4) {
          
           value=0.0;
           distance=0.0;
           htotal=0;
 
           for(h=0;h<max_zstations;h++) {
 
               if(zdata[pcpn_day].stn[h].zlevel2[time_pos].data < 0 )
                                                 continue;
 
               if(zdata[pcpn_day].stn[h].zlevel2[time_pos].qual != 5 &&
                  zdata[pcpn_day].stn[h].zlevel2[time_pos].qual != 8 &&
                  zdata[pcpn_day].stn[h].zlevel2[time_pos].qual != 2)
                              continue;
 
               dist1=(double)(hrap_grid->coord[i][j].x-
                      zstation[h].x);    
               dist2=(double)(hrap_grid->coord[i][j].y-
                      zstation[h].y);  
 
               dist=pow(dist1,2)+pow(dist2,2);
 
               if(dist < .00001)
                        dist=.00001;
               
               dist=1/dist;
 
               temp=zdata[pcpn_day].stn[h].zlevel2[time_pos].data*dist;
 
               value=value+temp;
 
               distance=distance+dist;                    
 
               htotal++;
 
         }
 
           }
 
           if(htotal == 0) 
               pcp->value[i][j]=-9999;
  
           else
               pcp->value[i][j]=(int)(value/distance*100);
  
 
     }
 
}
 
time_pos=100+pcpn_day*4+3-pcpn_time;
 
pcp_in_use[time_pos]=1;
 
write_file("pcp",time_pos,pcp);
 
}