SHAPE file for draw polyf

Graziano Giuliani giuliani at LAMMA.RETE.TOSCANA.IT
Wed Dec 14 10:21:53 EST 2005


Dear all,
I have made a little hack to grads code (version 1.9b4 codebase, with some 
patches applied of whom I did speak about in previous message of 4 November, 
just some cosmetic to build against recent OpenDAP and some clean-up for 
Debian install with gcc 4.0) to directly use a shapefile for the draw polyf 
function.
The Shapefile format is a working and interchange format promulagated by ESRI 
for simple vector data with attributes. The file format actually consists of 
three files.
XXX.shp - holds the actual vertices.
XXX.shx - hold index data pointing to the structures in the .shp file.
XXX.dbf - holds the attributes in xBase (dBase) format.
There exist a simple C library to access data available from 
http://shapelib.maptools.org/ maintained by Frank Warmerdam which is LGPL (or 
MIT style, You can choose).
As a lot of "shapefiles" are available on the web, (freegis for example, and 
also gshhs) I wanted to bypass the ASCII lpolyxxxx files and the basemap gs 
script.
Sooo...

43a44,47
> #if USESHAPE == 1
> #include <shapefil.h>
> #endif
>
1081a1086,1166
> #if USESHAPE == 1
>     if (i==3) {
>       SHPHandle hSHP;
>       int nShapeType, nshp, nobj, npart, pstart, nvert, ivert;
>       double adfMinBound[4], adfMaxBound[4];
>       float shplon, shplat;
>       SHPObject *psShape;
>       char shpname[1024];
>       char message[2048];
>
>       sscanf(nxtwrd(cmd), "%s", shpname);
>       sprintf(message, "Using shapefile %s\n", gxgnam(shpname));
>       gaprnt (0,message);
>       hSHP = SHPOpen(gxgnam(shpname), "rb");
>       if (hSHP == NULL)
>       {
>         gaprnt (0,"SHAPEFILE not found.\n");
>         gaprnt (0,"DRAW error: Syntax is DRAW POLYF shapefile object 
part\n");
>         return (1);
>       }
>       SHPGetInfo(hSHP, &nshp, &nShapeType, adfMinBound, adfMaxBound);
>       cmd = nxtwrd(cmd);
>       cmd = nxtwrd(cmd);
>       if (intprs(cmd, &nobj) == NULL)
>       {
>         gaprnt (0,"DRAW error: Invalid object number in shapefile\n");
>         return (1);
>       }
>       cmd = nxtwrd(cmd);
>       if (intprs(cmd, &npart) == NULL)
>       {
>         gaprnt (0,"DRAW error: Invalid part of object number\n");
>         return (1);
>       }
>       if (nobj < 0 || nobj >= nshp)
>       {
>         gaprnt (0,"DRAW error: Invalid object number in shapefile\n");
>         return (1);
>       }
>       psShape = SHPReadObject(hSHP, nobj);
>       if (npart < 0 || npart >= psShape->nParts)
>       {
>         gaprnt (0,"DRAW error: Invalid part of object number\n");
>         return (1);
>       }
>       pstart = psShape->panPartStart[npart];
>       if (npart == psShape->nParts-1)
>       {
>         nvert = psShape->nVertices-pstart;
>       }
>       else
>       {
>         nvert = psShape->panPartStart[npart+1]-pstart;
>       }
>       xy = (float *)malloc(2*sizeof(float)*(nvert+2));
>       if (xy==NULL) {
>         gaprnt (0,"DRAW error: Memory allocation error\n");
>         return (1);
>       }
>       i = 0;
>       for (ivert = 0; ivert < nvert; ivert ++)
>       {
>         shplon = psShape->padfX[ivert+pstart];
>         shplat = psShape->padfY[ivert+pstart];
>         if ((shplon >= pcm->dmin[0] && shplon <= pcm->dmax[0]) &&
>             (shplat >= pcm->dmin[1] && shplat <= pcm->dmax[1]))
>         {
>           gxconv (shplon,shplat,&x,&y,2);
>           xy[i]   = x;
>           xy[i+1] = y;
>           i += 2;
>         }
>       }
>       if (i == 0)
>       {
>         gaprnt (0,"DRAW warning: Outside window.\n");
>         return (0);
>       }
>     }
>     else if (i<6) {
> #else
1082a1168
> #endif
1086,1096c1172,1185
<     xy = (float *)malloc(sizeof(float)*(i+2));
<     if (xy==NULL) {
<       gaprnt (0,"DRAW error: Memory allocation error\n");
<       return (1);
<     }
<     i = 0;
<     while ( (cmd = nxtwrd(cmd)) != NULL) {
<       if ( valprs(cmd,xy+i) == NULL ) {
<         gaprnt (0,"DRAW error: Invalid polyf coordinate\n");
<         free (xy);
<         return(1);
---
>     else {
>       xy = (float *) malloc(sizeof(float)*(i+2));
>       if (xy==NULL) {
>         gaprnt (0,"DRAW error: Memory allocation error\n");
>         return (1);
>       }
>       i = 0;
>       while ( (cmd = nxtwrd(cmd)) != NULL) {
>         if ( valprs(cmd,xy+i) == NULL ) {
>           gaprnt (0,"DRAW error: Invalid polyf coordinate\n");
>           free (xy);
>           return(1);
>         }
>         i++;
1098d1186
<       i++;

this way a new syntax is added to draw polyf:

draw polyf shapefile object partition

as in

draw polyf shapes/gshhs/gshhs_land 128 0

using shapefiles from 
http://www.ngdc.noaa.gov/mgg/shorelines/data/gshhs/version1.3/shapefiles

Just in case I attach another simple program to be linked with same shapelib 
to convert any shapefile or directly the gshhs format or any simple ascii 
format into a map data for the GrADS program.

Cheers,
                 Graziano.


-- 
                             \ | /
                             (@ @)
 -------------------------o00-(_)-00o -----------------------------

 LaMMA - Laboratorio per la Meteorologia e la Modellistica Ambientale
         Laboratory for Meteorology and Environmental Modelling
 Via Madonna del Piano, 50019 Sesto Fiorentino (FI)

     tel: + 39 055 4483049
     fax: + 39 055 444083
     web: www.lamma.rete.toscana.it
  e-mail: giuliani at lamma.rete.toscana.it
-------------- next part --------------
A non-text attachment was scrubbed...
Name: create_map.c
Type: text/x-csrc
Size: 14238 bytes
Desc: not available
Url : http://gradsusr.org/pipermail/gradsusr/attachments/20051214/124fa573/attachment.bin 
-------------- next part --------------
A non-text attachment was scrubbed...
Name: Makefile
Type: text/x-makefile
Size: 463 bytes
Desc: not available
Url : http://gradsusr.org/pipermail/gradsusr/attachments/20051214/124fa573/attachment-0001.bin 
-------------- next part --------------
A non-text attachment was scrubbed...
Name: not available
Type: application/pgp-signature
Size: 189 bytes
Desc: not available
Url : http://gradsusr.org/pipermail/gradsusr/attachments/20051214/124fa573/attachment-0002.bin 


More information about the gradsusr mailing list