Rich,
The same web page offers a C program that generates all of the
lat/lon values as one step in translating to an alternate file
format. I found it useful to change this program to output the
coordinates as a free standing Netcdf file, which can then be easily
combined or applied with with the precip data.
Please see the attached program. If needed, it should be fairly
easy to add the grid_mapping elements that you mentioned.
--Dave A.
NOAA/PSD/CIRES
On 12/17/2010 1:45 PM, Rich Signell wrote:
NetCDF folks,
We are looking for a netcdf file that contains the lon/lat values and
perhaps the grid_mapping parameters for the precip data available
from:
http://water.weather.gov/precip/download.php
You can download gridded precip data as NetCDF, but the files do not
contain the lon/lat values, and no such netcdf "grid" file seems
available on the site. There is a Shapefile with all the lon/lat
points that we could convert to NetCDF, but surely someone must have
already done this, and perhaps even added the x,y coordinates and the
polar_stereographic parameters in a grid_mapping to make it all
CF-compliant?
Thanks,
Rich
/* Program hrap-coordinates-to-netcdf.c
This is a modified version of program nctoasc.c
written by Ken Pavelle
National Weather Service
Arkansas-Red Basin River Forecast Center
2010-sep-30 nctoasc.c:
Original version, by Ken Pavelle, NWS.
Downloaded from NWS/AHPS website:
http://water.weather.gov/precip/download.php
This version converts a Netcdf precip file WITHOUT
coordinates, into a text file WITH lat/lon coordinates.
2010-dec-20 hrap-coordinates-to-netcdf.c:
Modified to extract coordinates only, by Dave Allured,
NOAA/PSD/CIRES.
Read a Netcdf precip file WITHOUT coordinates,
only to get the HRAP grid location parameters.
Output a Netcdf file containing ONLY the lat and lon
coordinates.
Promote computations to double precision.
Input -- Netcdf file - Containing HRAP grid location parameters.
Output -- Netcdf file - Containing ONLY lat and lon coordinate
arrays matching the specified input file.
USAGE: "hrap-coordinates-to-netcdf infile outfile"
NOTE: This version computes coordinates for ANY HRAP grid,
not just CONUS and Puerto Rico.
There is one constraint. The HRAP grid coordinate system
is defined on a polar stereographic projection with origin
at 60°N / 105°W. The origin is currently hard coded, and
all HRAP grid parameters for this program must be relative
to this origin. */
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <malloc.h>
#include <math.h>
#include <netcdf.h>
typedef struct _HRAP {
double x;
double y;
} HRAP;
/*--------------------------------------------------------------------------*/
HRAP HrapToLatLong(HRAP hrap)
{
double raddeg = 57.29577951;
double earthrad = 6371.2;
double stdlon = 105.;
double mesh_len = 4.7625;
double tlat, rr, gi, ang, x, y;
HRAP ll;
tlat=60./raddeg;
x = hrap.x - 401.;
y = hrap.y - 1601.;
rr = x*x + y*y;
gi = ((earthrad * (1+sin(tlat)))/mesh_len);
gi=gi*gi;
ll.y = asin((gi-rr)/(gi+rr))*raddeg;
ang = atan2(y,x)*raddeg;
if (ang < 0) ang=ang+360.;
ll.x = 270+stdlon-ang;
if (ll.x < 0) ll.x=ll.x+360;
if (ll.x > 360) ll.x=ll.x-360;
return ll;
}
/*--------------------------------------------------------------------------*/
void check (int nc_status, const char *func_name)
{
if (nc_status != NC_NOERR)
{
printf ("*** Netcdf error, %s returned error no. %d:\n", func_name,
nc_status);
printf ("*** %s\n", nc_strerror (nc_status));
exit (-1);
}
}
/*--------------------------------------------------------------------------*/
void *malloc_check (size_t size, const char *varname)
{
void *ptr = malloc (size);
if (ptr == NULL)
{
printf ("*** Malloc error, can't allocate %lu bytes for %s.\n", size,
varname);
exit(1);
}
return (ptr);
}
/*--------------------------------------------------------------------------*/
int main (int argc, char *argv [])
{
int i, j;
char infile[150], outfile[150];
int hrap_xor_id, hrap_yor_id, hrapx_id, hrapy_id; /*netCDF stuff*/
int ncid, xid, yid, lon2d_id, lat2d_id;
int ndims, dimids[2];
double hrap_xor, hrap_yor;
size_t hrapx, hrapy, vsize, index;
HRAP hrap, latlon;
float *lat2d, *lon2d;
static char title[] = "HRAP 2-D lat and lon coordinates";
static char history[] = "Generated by hrap-coordinates-to-netcdf.c";
static char lat_lname[] = "latitude";
static char lon_lname[] = "longitude";
static char lat_units[] = "degrees_north";
static char lon_units[] = "degrees_east";
/* Get command line args */
if (argc == 3) {
sprintf(infile, "%s", argv[1]);
sprintf(outfile, "%s", argv[2]);
} else {
printf ("USAGE: hrap-coordinates-to-netcdf infile outfile\n");
return 1;
}
/* Open Netcdf input file. */
printf ("\n");
printf ("Read input file: %s\n", infile);
check (nc_open (infile, NC_NOWRITE, &ncid), "nc_open");
/* Read dimensions and HRAP grid parameters. */
/* Use Netcdf automatic type conversion, as needed. */
check (nc_inq_dimid (ncid,"hrapx",&hrapx_id), "nc_inq_dimid");
check (nc_inq_dimlen (ncid,hrapx_id,&hrapx), "nc_inq_dimlen");
check (nc_inq_dimid (ncid,"hrapy",&hrapy_id), "nc_inq_dimid");
check (nc_inq_dimlen (ncid,hrapy_id,&hrapy), "nc_inq_dimlen");
check (nc_inq_varid (ncid,"hrap_xor",&hrap_xor_id), "nc_inq_varid");
check (nc_get_var1_double (ncid,hrap_xor_id,0,&hrap_xor), "nc_get_var1_float");
check (nc_inq_varid (ncid,"hrap_yor",&hrap_yor_id), "nc_inq_varid");
check (nc_get_var1_double (ncid,hrap_yor_id,0,&hrap_yor), "nc_get_var1_float");
/* Close input file. */
check (nc_close(ncid), "nc_close");
printf ("nx = hrapx = %lu\n", (long) hrapx);
printf ("ny = hrapy = %lu\n", (long) hrapy);
printf ("Grid X offset = hrap_xor = %lu\n", (long) hrap_xor);
printf ("Grid Y offset = hrap_yor = %lu\n", (long) hrap_yor);
/* Compute 2-D coordinates for HRAP grid. */
printf ("\n");
printf ("Compute 2-D coordinates for HRAP grid.\n");
vsize = hrapx * hrapy * sizeof(*lat2d); /* allocate arrays */
lat2d = malloc_check (vsize, "lat2d"); /* simple 1-D indexing */
lon2d = malloc_check (vsize, "lon2d");
for (i = 0; i< hrapy; i++) {
for (j = 0; j< hrapx; j++) {
hrap.x = j + hrap_xor + 0.5;
hrap.y = i + hrap_yor + 0.5;
latlon = HrapToLatLong (hrap);
index = j + hrapx * i; /* index for 2-D contiguous */
lon2d[index] = (float) -latlon.x; /* convert lons WEST to lons EAST */
lat2d[index] = (float) latlon.y;
}
}
/* Write 2-D coordinates to output file. */
printf ("Write 2-D coordinates to output file: %s\n", outfile);
/* Note: No overwrite protect in this version. */
check (nc_create (outfile, NC_CLOBBER, &ncid), "nc_create"); /* create new
file */
check (nc_def_dim (ncid, "hrapy", hrapy, &yid), "nc_def_dim"); /* define
dimensions */
check (nc_def_dim (ncid, "hrapx", hrapx, &xid), "nc_def_dim");
/* Define 2-D lat and lon arrays. */
ndims = 2;
dimids[0] = yid;
dimids[1] = xid;
check (nc_def_var (ncid, "lat2d", NC_FLOAT, ndims, dimids, &lat2d_id),
"nc_def_var");
check (nc_def_var (ncid, "lon2d", NC_FLOAT, ndims, dimids, &lon2d_id),
"nc_def_var");
/* Add lat and lon attributes. */
check (nc_put_att_text (ncid, lat2d_id, "long_name", strlen(lat_lname),
lat_lname),
"nc_put_att_text");
check (nc_put_att_text (ncid, lon2d_id, "long_name", strlen(lon_lname),
lon_lname),
"nc_put_att_text");
check (nc_put_att_text (ncid, lat2d_id, "units", strlen(lat_units), lat_units),
"nc_put_att_text");
check (nc_put_att_text (ncid, lon2d_id, "units", strlen(lon_units), lon_units),
"nc_put_att_text");
/* Add global attributes. */
check (nc_put_att_text (ncid, NC_GLOBAL, "title", strlen(title), title),
"nc_put_att_text");
check (nc_put_att_text (ncid, NC_GLOBAL, "history", strlen(history), history),
"nc_put_att_text");
check (nc_put_att_text (ncid, NC_GLOBAL, "grid_specification_from",
strlen(infile), infile),
"nc_put_att_text");
check (nc_enddef (ncid), "nc_enddef"); /* leave define mode, enter
data mode */
/* Write lat and lon arrays to output file. */
check (nc_put_var_float (ncid, lat2d_id, lat2d), "nc_put_var_float");
check (nc_put_var_float (ncid, lon2d_id, lon2d), "nc_put_var_float");
check (nc_close(ncid), "nc_close"); /* flush buffers and close output file
*/
printf("Done.\n\n");
free (lon2d); /* clean up dynamic memory */
free (lat2d);
return 0;
}
# 2009-oct-06 Original Makefile for nctoasc.c, from NWS/AHPS website.
# 2010-dec-20 Configured for local compile. (D. Allured)
# Add main program hrap-coordinates-to-netcdf.c.
Netcdf_root = /Users/dallured/lib/gfortran/netcdf3
CFLAGS = -I$(Netcdf_root)/include -I/usr/local/include -Wall -Werror -g
LFLAGS = -L$(Netcdf_root)/lib
LIBS = -lnetcdf -lz -lm -lc
CC = /usr/local/bin/gcc
hrap-coordinates-to-netcdf: hrap-coordinates-to-netcdf.o makefile
$(CC) -o $@ $(LFLAGS) $@.o $(LIBS);
nctoasc: nctoasc.o makefile
$(CC) -o $@.exe $(LFLAGS) $@.o $(LIBS);
clean:
rm *.o
