Salome HOME
refs# 1917 + protection against null shapes wires in GetMiddlePoint()
[modules/hydro.git] / src / shapelib / shputils.c
1 /******************************************************************************
2  * $Id: shputils.c,v 1.10 2007-12-13 19:59:23 fwarmerdam Exp $
3  *
4  * Project:  Shapelib
5  * Purpose:  
6  *   Altered "shpdump" and "dbfdump" to allow two files to be appended.
7  *   Other Functions:
8  *     Selecting from the DBF before the write occurs.
9  *     Change the UNITS between Feet and Meters and Shift X,Y.
10  *     Clip and Erase boundary.  The program only passes thru the
11  *     data once.
12  *
13  *   Bill Miller   North Carolina - Department of Transporation 
14  *   Feb. 1997 -- bmiller@dot.state.nc.us
15  *         There was not a lot of time to debug hidden problems;
16  *         And the code is not very well organized or documented.
17  *         The clip/erase function was not well tested.
18  *   Oct. 2000 -- bmiller@dot.state.nc.us
19  *         Fixed the problem when select is using numbers
20  *         larger than short integer.  It now reads long integer.
21  *   NOTE: DBF files created using windows NT will read as a string with
22  *         a length of 381 characters.  This is a bug in "dbfopen".
23  *
24  *
25  * Author:   Bill Miller (bmiller@dot.state.nc.us)
26  *
27  ******************************************************************************
28  * Copyright (c) 1999, Frank Warmerdam
29  *
30  * This software is available under the following "MIT Style" license,
31  * or at the option of the licensee under the LGPL (see LICENSE.LGPL).  This
32  * option is discussed in more detail in shapelib.html.
33  *
34  * --
35  * 
36  * Permission is hereby granted, free of charge, to any person obtaining a
37  * copy of this software and associated documentation files (the "Software"),
38  * to deal in the Software without restriction, including without limitation
39  * the rights to use, copy, modify, merge, publish, distribute, sublicense,
40  * and/or sell copies of the Software, and to permit persons to whom the
41  * Software is furnished to do so, subject to the following conditions:
42  *
43  * The above copyright notice and this permission notice shall be included
44  * in all copies or substantial portions of the Software.
45  *
46  * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
47  * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
48  * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
49  * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
50  * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
51  * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
52  * DEALINGS IN THE SOFTWARE.
53  ******************************************************************************
54  *
55  * $Log: shputils.c,v $
56  * Revision 1.10  2007-12-13 19:59:23  fwarmerdam
57  * reindent code, avoid some warnings.
58  *
59  * Revision 1.9  2004/01/14 14:56:00  fwarmerdam
60  * Some cleanlyness improvements.
61  *
62  * Revision 1.8  2004/01/14 14:40:22  fwarmerdam
63  * Fixed exit() call to include code.
64  *
65  * Revision 1.7  2003/02/25 17:20:22  warmerda
66  * Set psCShape to NULL after SHPDestroyObject() to avoid multi-frees of
67  * the same memory ... as submitted by Fred Fox.
68  *
69  * Revision 1.6  2001/08/28 13:57:14  warmerda
70  * fixed DBFAddField return value check
71  *
72  * Revision 1.5  2000/11/02 13:52:48  warmerda
73  * major upgrade from Bill Miller
74  *
75  * Revision 1.4  1999/11/05 14:12:05  warmerda
76  * updated license terms
77  *
78  * Revision 1.3  1998/12/03 15:47:39  warmerda
79  * Did a bunch of rewriting to make it work with the V1.2 API.
80  *
81  * Revision 1.2  1998/06/18 01:19:49  warmerda
82  * Made C++ compilable.
83  *
84  * Revision 1.1  1997/05/27 20:40:27  warmerda
85  * Initial revision
86  */
87
88 #include "shapefil.h"
89 #include "string.h"
90 #include <stdlib.h>
91
92 SHP_CVSID("$Id: shputils.c,v 1.10 2007-12-13 19:59:23 fwarmerdam Exp $")
93
94 #ifndef FALSE
95 #  define FALSE         0
96 #  define TRUE          1
97 #endif
98
99 char            infile[80], outfile[80], temp[400];
100
101 /* Variables for shape files */
102 SHPHandle       hSHP;
103 SHPHandle       hSHPappend;
104 int             nShapeType, nEntities, iPart;
105 int             nShapeTypeAppend, nEntitiesAppend;
106 SHPObject       *psCShape;
107 double          adfBoundsMin[4], adfBoundsMax[4];
108
109
110 /* Variables for DBF files */
111 DBFHandle       hDBF;
112 DBFHandle       hDBFappend;
113     
114 DBFFieldType    iType;
115 DBFFieldType    jType;
116     
117 char    iszTitle[12];
118 char    jszTitle[12];
119
120 int     *pt;
121 char    iszFormat[32], iszField[1024];
122 char    jszFormat[32], jszField[1024];
123 int     i, ti, iWidth, iDecimals, iRecord;
124 int     j, tj, jWidth, jDecimals, jRecord;
125
126
127 int clip_boundary();
128 double findunit(char *unit);
129 void openfiles(void);
130 void setext(char *pt, char *ext);
131 int strncasecmp2(char *s1, char *s2, int n);
132 void mergefields(void);
133 void findselect(void);
134 void showitems(void);
135 int selectrec();
136 void check_theme_bnd();
137 int clip_boundary();
138 void error();
139
140
141 /* -------------------------------------------------------------------- */
142 /* Variables for the DESCRIBE function */
143 /* -------------------------------------------------------------------- */
144    int       ilist = FALSE, iall = FALSE;
145 /* -------------------------------------------------------------------- */
146 /* Variables for the SELECT function */
147 /* -------------------------------------------------------------------- */
148    int       found = FALSE, newdbf = FALSE;
149    char      selectitem[40], *cpt;
150    long int  selectvalues[150], selcount=0;
151    int       iselect = FALSE, iselectitem = -1;
152    int       iunselect = FALSE;
153
154 /* -------------------------------------------------------------------- */
155 /* Variables for the CLIP and ERASE functions */
156 /* -------------------------------------------------------------------- */
157    double  cxmin, cymin, cxmax, cymax; 
158    int     iclip  = FALSE, ierase = FALSE;
159    int     itouch = FALSE, iinside = FALSE, icut = FALSE;
160    int     ibound = FALSE, ipoly = FALSE;
161    char    clipfile[80];
162
163 /* -------------------------------------------------------------------- */
164 /* Variables for the FACTOR function */
165 /* -------------------------------------------------------------------- */
166    double  infactor,outfactor,factor = 0;  /* NO FACTOR */
167    int     iunit = FALSE;
168    int     ifactor = FALSE;
169
170    
171 /* -------------------------------------------------------------------- */
172 /* Variables for the SHIFT function */
173 /* -------------------------------------------------------------------- */
174    double  xshift = 0, yshift = 0;  /* NO SHIFT */
175       
176 int main( int argc, char ** argv )
177 {
178
179 /* -------------------------------------------------------------------- */
180 /*      Check command line usage.                                       */
181 /* -------------------------------------------------------------------- */
182     if( argc < 2 ) error();
183     strcpy(infile, argv[1]);
184     if (argc > 2) {
185         strcpy(outfile,argv[2]);
186         if (strncasecmp2(outfile, "LIST",0) == 0) { ilist = TRUE; }
187         if (strncasecmp2(outfile, "ALL",0) == 0)  { iall  = TRUE; }
188     } 
189     if (ilist || iall || argc == 2 ) {
190         setext(infile, "shp");
191         printf("DESCRIBE: %s\n",infile);
192         strcpy(outfile,"");
193     }
194 /* -------------------------------------------------------------------- */
195 /*      Look for other functions on the command line. (SELECT, UNIT)    */
196 /* -------------------------------------------------------------------- */
197     for (i = 3; i < argc; i++)
198     {
199         if ((strncasecmp2(argv[i],  "SEL",3) == 0) ||
200             (strncasecmp2(argv[i],  "UNSEL",5) == 0))
201         {
202             if (strncasecmp2(argv[i],  "UNSEL",5) == 0) iunselect=TRUE;
203             i++;
204             if (i >= argc) error();
205             strcpy(selectitem,argv[i]);
206             i++;
207             if (i >= argc) error();
208             selcount=0;
209             strcpy(temp,argv[i]);
210             cpt=temp;
211             tj = atoi(cpt);
212             ti = 0;
213             while (tj>0) {
214                 selectvalues[selcount] = tj;
215                 while( *cpt >= '0' && *cpt <= '9')
216                     cpt++; 
217                 while( *cpt > '\0' && (*cpt < '0' || *cpt > '9') )
218                     cpt++; 
219                 tj=atoi(cpt);
220                 selcount++;
221             }
222             iselect=TRUE;
223         }  /*** End SEL & UNSEL ***/
224         else
225             if ((strncasecmp2(argv[i], "CLIP",4) == 0) ||
226                 (strncasecmp2(argv[i],  "ERASE",5) == 0))
227             {
228                 if (strncasecmp2(argv[i],  "ERASE",5) == 0) ierase=TRUE;
229                 i++;
230                 if (i >= argc) error();
231                 strcpy(clipfile,argv[i]);
232                 sscanf(argv[i],"%lf",&cxmin);
233                 i++;
234                 if (i >= argc) error();
235                 if (strncasecmp2(argv[i],  "BOUND",5) == 0) {
236                     setext(clipfile, "shp");
237                     hSHP = SHPOpen( clipfile, "rb" );
238                     if( hSHP == NULL )
239                     {
240                         printf( "ERROR: Unable to open the clip shape file:%s\n", clipfile );
241                         exit( 1 );
242                     }
243                     SHPGetInfo( hSHPappend, NULL, NULL,
244                                 adfBoundsMin, adfBoundsMax );
245                     cxmin = adfBoundsMin[0];
246                     cymin = adfBoundsMin[1];
247                     cxmax = adfBoundsMax[0];
248                     cymax = adfBoundsMax[1];
249                     printf("Theme Clip Boundary: (%lf,%lf) - (%lf,%lf)\n",
250                            cxmin, cymin, cxmax, cymax);
251                     ibound=TRUE;
252                 } else {  /*** xmin,ymin,xmax,ymax ***/
253                     sscanf(argv[i],"%lf",&cymin);
254                     i++;
255                     if (i >= argc) error();
256                     sscanf(argv[i],"%lf",&cxmax);
257                     i++;
258                     if (i >= argc) error();
259                     sscanf(argv[i],"%lf",&cymax);
260                     printf("Clip Box: (%lf,%lf) - (%lf,%lf)\n",cxmin, cymin, cxmax, cymax);
261                 }
262                 i++;
263                 if (i >= argc) error();
264                 if      (strncasecmp2(argv[i], "CUT",3) == 0)    icut=TRUE;
265                 else if (strncasecmp2(argv[i], "TOUCH",5) == 0)  itouch=TRUE;
266                 else if (strncasecmp2(argv[i], "INSIDE",6) == 0) iinside=TRUE;
267                 else error();
268                 iclip=TRUE;
269             } /*** End CLIP & ERASE ***/
270             else if (strncasecmp2(argv[i],  "FACTOR",0) == 0)
271             {
272                 i++;
273                 if (i >= argc) error();
274                 infactor=findunit(argv[i]);
275                 if (infactor == 0) error();
276                 iunit=TRUE;
277                 i++;
278                 if (i >= argc) error();
279                 outfactor=findunit(argv[i]);
280                 if (outfactor == 0)
281                 {
282                     sscanf(argv[i],"%lf",&factor);
283                     if (factor == 0) error();
284                 }
285                 if (factor == 0)
286                 {
287                     if (infactor ==0)
288                     { puts("ERROR: Input unit must be defined before output unit"); exit(1); }
289                     factor=infactor/outfactor;
290                 }
291                 printf("Output file coordinate values will be factored by %lg\n",factor);
292                 ifactor=(factor != 1); /* True if a valid factor */
293             } /*** End FACTOR ***/
294             else if (strncasecmp2(argv[i],"SHIFT",5) == 0)
295             {
296                 i++;
297                 if (i >= argc) error();
298                 sscanf(argv[i],"%lf",&xshift);
299                 i++;
300                 if (i >= argc) error();
301                 sscanf(argv[i],"%lf",&yshift);
302                 iunit=TRUE;
303                 printf("X Shift: %lg   Y Shift: %lg\n",xshift,yshift);
304             } /*** End SHIFT ***/
305             else {
306                 printf("ERROR: Unknown function %s\n",argv[i]);  error();
307             }
308     }
309 /* -------------------------------------------------------------------- */
310 /*      If there is no data in this file let the user know.             */
311 /* -------------------------------------------------------------------- */
312     openfiles();  /* Open the infile and the outfile for shape and dbf. */
313     if( DBFGetFieldCount(hDBF) == 0 )
314     {
315         puts( "There are no fields in this table!" );
316         exit( 1 );
317     }
318 /* -------------------------------------------------------------------- */
319 /*      Print out the file bounds.                                      */
320 /* -------------------------------------------------------------------- */
321     iRecord = DBFGetRecordCount( hDBF );
322     SHPGetInfo( hSHP, NULL, NULL, adfBoundsMin, adfBoundsMax );
323
324     printf( "Input Bounds:  (%lg,%lg) - (%lg,%lg)   Entities: %d   DBF: %d\n",
325             adfBoundsMin[0], adfBoundsMin[1],
326             adfBoundsMax[0], adfBoundsMax[1],
327             nEntities, iRecord );
328             
329     if (strcmp(outfile,"") == 0) /* Describe the shapefile; No other functions */
330     {
331         ti = DBFGetFieldCount( hDBF );
332         showitems();
333         exit(0);
334     }
335      
336     if (iclip) check_theme_bnd();
337     
338     jRecord = DBFGetRecordCount( hDBFappend );
339     SHPGetInfo( hSHPappend, NULL, NULL, adfBoundsMin, adfBoundsMax );
340     if (nEntitiesAppend == 0)
341         puts("New Output File\n");
342     else
343         printf( "Append Bounds: (%lg,%lg)-(%lg,%lg)   Entities: %d  DBF: %d\n",
344                 adfBoundsMin[0], adfBoundsMin[1],
345                 adfBoundsMax[0], adfBoundsMax[1],
346                 nEntitiesAppend, jRecord );
347     
348 /* -------------------------------------------------------------------- */
349 /*      Find matching fields in the append file or add new items.       */
350 /* -------------------------------------------------------------------- */
351     mergefields();
352 /* -------------------------------------------------------------------- */
353 /*      Find selection field if needed.                                 */
354 /* -------------------------------------------------------------------- */
355     if (iselect)    findselect();
356
357 /* -------------------------------------------------------------------- */
358 /*  Read all the records                                                */
359 /* -------------------------------------------------------------------- */
360     jRecord = DBFGetRecordCount( hDBFappend );
361     for( iRecord = 0; iRecord < nEntities; iRecord++)  /** DBFGetRecordCount(hDBF) **/
362     {
363 /* -------------------------------------------------------------------- */
364 /*      SELECT for values if needed. (Can the record be skipped.)       */
365 /* -------------------------------------------------------------------- */
366         if (iselect)
367             if (selectrec() == 0) goto SKIP_RECORD;   /** SKIP RECORD **/
368
369 /* -------------------------------------------------------------------- */
370 /*      Read a Shape record                                             */
371 /* -------------------------------------------------------------------- */
372         psCShape = SHPReadObject( hSHP, iRecord );
373
374 /* -------------------------------------------------------------------- */
375 /*      Clip coordinates of shapes if needed.                           */
376 /* -------------------------------------------------------------------- */
377         if (iclip)
378             if (clip_boundary() == 0) goto SKIP_RECORD; /** SKIP RECORD **/
379
380 /* -------------------------------------------------------------------- */
381 /*      Read a DBF record and copy each field.                          */
382 /* -------------------------------------------------------------------- */
383         for( i = 0; i < DBFGetFieldCount(hDBF); i++ )
384         {
385 /* -------------------------------------------------------------------- */
386 /*      Store the record according to the type and formatting           */
387 /*      information implicit in the DBF field description.              */
388 /* -------------------------------------------------------------------- */
389             if (pt[i] > -1)  /* if the current field exists in output file */
390             {
391                 switch( DBFGetFieldInfo( hDBF, i, NULL, &iWidth, &iDecimals ) )
392                 {
393                   case FTString:
394                   case FTLogical:
395                     DBFWriteStringAttribute(hDBFappend, jRecord, pt[i],
396                                             (DBFReadStringAttribute( hDBF, iRecord, i )) );
397                     break;
398
399                   case FTInteger:
400                     DBFWriteIntegerAttribute(hDBFappend, jRecord, pt[i],
401                                              (DBFReadIntegerAttribute( hDBF, iRecord, i )) );
402                     break;
403
404                   case FTDouble:
405                     DBFWriteDoubleAttribute(hDBFappend, jRecord, pt[i],
406                                             (DBFReadDoubleAttribute( hDBF, iRecord, i )) );
407                     break;
408
409                   case FTInvalid:
410                     break;
411                 }
412             }
413         }
414         jRecord++;
415 /* -------------------------------------------------------------------- */
416 /*      Change FACTOR and SHIFT coordinates of shapes if needed.        */
417 /* -------------------------------------------------------------------- */
418         if (iunit)
419         {
420             for( j = 0; j < psCShape->nVertices; j++ ) 
421             {
422                 psCShape->padfX[j] = psCShape->padfX[j] * factor + xshift;
423                 psCShape->padfY[j] = psCShape->padfY[j] * factor + yshift;
424             }
425         }
426         
427 /* -------------------------------------------------------------------- */
428 /*      Write the Shape record after recomputing current extents.       */
429 /* -------------------------------------------------------------------- */
430         SHPComputeExtents( psCShape );
431         SHPWriteObject( hSHPappend, -1, psCShape );
432
433       SKIP_RECORD:
434         SHPDestroyObject( psCShape );
435         psCShape = NULL;
436         j=0;
437     }
438
439 /* -------------------------------------------------------------------- */
440 /*      Print out the # of Entities and the file bounds.                */
441 /* -------------------------------------------------------------------- */
442     jRecord = DBFGetRecordCount( hDBFappend );
443     SHPGetInfo( hSHPappend, &nEntitiesAppend, &nShapeTypeAppend,
444                 adfBoundsMin, adfBoundsMax );
445     
446     printf( "Output Bounds: (%lg,%lg) - (%lg,%lg)   Entities: %d  DBF: %d\n\n",
447             adfBoundsMin[0], adfBoundsMin[1],
448             adfBoundsMax[0], adfBoundsMax[1],
449             nEntitiesAppend, jRecord );
450
451 /* -------------------------------------------------------------------- */
452 /*      Close the both shapefiles.                                      */
453 /* -------------------------------------------------------------------- */
454     SHPClose( hSHP );
455     SHPClose( hSHPappend );
456     DBFClose( hDBF );
457     DBFClose( hDBFappend );
458     if (nEntitiesAppend == 0) {
459         puts("Remove the output files.");
460         setext(outfile, "dbf");
461         remove(outfile);
462         setext(outfile, "shp");
463         remove(outfile);
464         setext(outfile, "shx");
465         remove(outfile);
466     }
467     return( 0 );
468 }
469
470
471 /************************************************************************/
472 /*                             openfiles()                              */
473 /************************************************************************/
474
475 void openfiles() {
476 /* -------------------------------------------------------------------- */
477 /*      Open the DBF file.                                              */
478 /* -------------------------------------------------------------------- */
479     setext(infile, "dbf");
480     hDBF = DBFOpen( infile, "rb" );
481     if( hDBF == NULL )
482     {
483         printf( "ERROR: Unable to open the input DBF:%s\n", infile );
484         exit( 1 );
485     }
486 /* -------------------------------------------------------------------- */
487 /*      Open the append DBF file.                                       */
488 /* -------------------------------------------------------------------- */
489     if (strcmp(outfile,"")) {
490         setext(outfile, "dbf");
491         hDBFappend = DBFOpen( outfile, "rb+" );
492         newdbf=0;
493         if( hDBFappend == NULL )
494         {
495             newdbf=1;
496             hDBFappend = DBFCreate( outfile );
497             if( hDBFappend == NULL )
498             {
499                 printf( "ERROR: Unable to open the append DBF:%s\n", outfile );
500                 exit( 1 );
501             }
502         }
503     }
504 /* -------------------------------------------------------------------- */
505 /*      Open the passed shapefile.                                      */
506 /* -------------------------------------------------------------------- */
507     setext(infile, "shp");
508     hSHP = SHPOpen( infile, "rb" );
509
510     if( hSHP == NULL )
511     {
512         printf( "ERROR: Unable to open the input shape file:%s\n", infile );
513         exit( 1 );
514     }
515
516     SHPGetInfo( hSHP, &nEntities, &nShapeType, NULL, NULL );
517
518 /* -------------------------------------------------------------------- */
519 /*      Open the passed append shapefile.                               */
520 /* -------------------------------------------------------------------- */
521     if (strcmp(outfile,"")) {
522         setext(outfile, "shp");
523         hSHPappend = SHPOpen( outfile, "rb+" );
524
525         if( hSHPappend == NULL )
526         {
527             hSHPappend = SHPCreate( outfile, nShapeType );
528             if( hSHPappend == NULL )
529             {
530                 printf( "ERROR: Unable to open the append shape file:%s\n",
531                         outfile );
532                 exit( 1 );
533             }
534         }
535         SHPGetInfo( hSHPappend, &nEntitiesAppend, &nShapeTypeAppend,
536                     NULL, NULL );
537
538         if (nShapeType != nShapeTypeAppend) 
539         {
540             puts( "ERROR: Input and Append shape files are of different types.");
541             exit( 1 );
542         }
543     }
544 }
545
546 /* -------------------------------------------------------------------- */
547 /*      Change the extension.  If there is any extension on the         */
548 /*      filename, strip it off and add the new extension                */
549 /* -------------------------------------------------------------------- */
550 void setext(char *pt, char *ext)
551 {
552 int i;
553     for( i = strlen(pt)-1; 
554          i > 0 && pt[i] != '.' && pt[i] != '/' && pt[i] != '\\';
555          i-- ) {}
556
557     if( pt[i] == '.' )
558         pt[i] = '\0';
559         
560     strcat(pt,".");
561     strcat(pt,ext);
562 }
563
564
565
566 /* -------------------------------------------------------------------- */
567 /*      Find matching fields in the append file.                        */
568 /*      Output file must have zero records to add any new fields.       */
569 /* -------------------------------------------------------------------- */
570 void mergefields()
571 {
572     int i,j;
573     ti = DBFGetFieldCount( hDBF );
574     tj = DBFGetFieldCount( hDBFappend );
575     /* Create a pointer array for the max # of fields in the output file */
576     pt = (int *) malloc( (ti+tj+1) * sizeof(int) ); 
577     
578     for( i = 0; i < ti; i++ )
579     {
580         pt[i]= -1;  /* Initial pt values to -1 */
581     }
582     /* DBF must be empty before adding items */
583     jRecord = DBFGetRecordCount( hDBFappend );
584     for( i = 0; i < ti; i++ )
585     {
586         iType = DBFGetFieldInfo( hDBF, i, iszTitle, &iWidth, &iDecimals );
587         found=FALSE;
588         {
589             for( j = 0; j < tj; j++ )   /* Search all field names for a match */
590             {
591                 jType = DBFGetFieldInfo( hDBFappend, j, jszTitle, &jWidth, &jDecimals );
592                 if (iType == jType && (strcmp(iszTitle, jszTitle) == 0) )
593                 {
594                     if (found || newdbf)
595                     {
596                         if (i == j)  pt[i]=j;
597                         printf("Warning: Duplicate field name found (%s)\n",iszTitle);
598                         /* Duplicate field name
599                            (Try to guess the correct field by position) */
600                     }
601                     else
602                     {
603                         pt[i]=j;  found=TRUE; 
604                     }
605                 }
606             }
607         }
608         
609         if (pt[i] == -1  && (! found) )  /* Try to force into an existing field */
610         {                                /* Ignore the field name, width, and decimal places */
611             jType = DBFGetFieldInfo( hDBFappend, j, jszTitle, &jWidth, &jDecimals );
612             if (iType == jType) 
613             {
614                 pt[i]=i;  found=1;
615             }
616         }
617         if ( (! found) &&  jRecord == 0)  /* Add missing field to the append table */
618         {                 /* The output DBF must be is empty */
619             pt[i]=tj;
620             tj++;
621             if( DBFAddField( hDBFappend, iszTitle, iType, iWidth, iDecimals )
622                 == -1 )
623             {
624                 printf( "Warning: DBFAddField(%s, TYPE:%d, WIDTH:%d  DEC:%d, ITEM#:%d of %d) failed.\n",
625                         iszTitle, iType, iWidth, iDecimals, (i+1), (ti+1) );
626                 pt[i]=-1;
627             }
628         }
629     }
630 }
631
632
633 void findselect()
634 {
635     /* Find the select field name */
636     iselectitem = -1;
637     for( i = 0; i < ti  &&  iselectitem < 0; i++ )
638     {
639         iType = DBFGetFieldInfo( hDBF, i, iszTitle, &iWidth, &iDecimals );
640         if (strncasecmp2(iszTitle, selectitem, 0) == 0) iselectitem = i;
641     }
642     if (iselectitem == -1) 
643     {
644         printf("Warning: Item not found for selection (%s)\n",selectitem);
645         iselect = FALSE;
646         iall = FALSE;
647         showitems();
648         printf("Continued... (Selecting entire file)\n");
649     }
650     /* Extract all of the select values (by field type) */
651     
652 }
653
654 void showitems()
655 {
656     char      stmp[40],slow[40],shigh[40];
657     double    dtmp,dlow,dhigh,dsum,mean;
658     long int  itmp,ilow,ihigh,isum;
659     long int  maxrec;
660     char      *pt;
661
662     printf("Available Items: (%d)",ti);
663     maxrec = DBFGetRecordCount(hDBF);
664     if (maxrec > 5000 && ! iall) 
665     { maxrec=5000; printf("  ** ESTIMATED RANGES (MEAN)  For more records use \"All\""); }
666     else  { printf("          RANGES (MEAN)"); }
667           
668     for( i = 0; i < ti; i++ )
669     {
670         switch( DBFGetFieldInfo( hDBF, i, iszTitle, &iWidth, &iDecimals ) )
671         {
672           case FTString:
673           case FTLogical:
674             strcpy(slow, "~");
675             strcpy(shigh,"\0");
676             printf("\n  String  %3d  %-16s",iWidth,iszTitle);
677             for( iRecord = 0; iRecord < maxrec; iRecord++ ) {
678                 strncpy(stmp,DBFReadStringAttribute( hDBF, iRecord, i ),39);
679                 if (strcmp(stmp,"!!") > 0) {
680                     if (strncasecmp2(stmp,slow,0)  < 0) strncpy(slow, stmp,39);
681                     if (strncasecmp2(stmp,shigh,0) > 0) strncpy(shigh,stmp,39);
682                 }
683             }
684             pt=slow+strlen(slow)-1; 
685             while(*pt == ' ') { *pt='\0'; pt--; }
686             pt=shigh+strlen(shigh)-1;
687             while(*pt == ' ') { *pt='\0'; pt--; }
688             if (strncasecmp2(slow,shigh,0) < 0)         printf("%s to %s",slow,shigh);
689             else if (strncasecmp2(slow,shigh,0) == 0)   printf("= %s",slow);
690             else        printf("No Values");
691             break;
692           case FTInteger:
693             printf("\n  Integer %3d  %-16s",iWidth,iszTitle);
694             ilow =  1999999999;
695             ihigh= -1999999999;
696             isum =  0;
697             for( iRecord = 0; iRecord < maxrec; iRecord++ ) {
698                 itmp = DBFReadIntegerAttribute( hDBF, iRecord, i );
699                 if (ilow > itmp)  ilow = itmp;
700                 if (ihigh < itmp) ihigh = itmp;
701                 isum = isum + itmp;
702             }
703             mean=isum/maxrec;
704             if (ilow < ihigh)       printf("%ld to %ld \t(%.1f)",ilow,ihigh,mean);
705             else if (ilow == ihigh) printf("= %ld",ilow);
706             else printf("No Values");
707             break;
708
709           case FTDouble:
710             printf("\n  Real  %3d,%d  %-16s",iWidth,iDecimals,iszTitle);
711             dlow =  999999999999999.0;
712             dhigh= -999999999999999.0;
713             dsum =  0;
714             for( iRecord = 0; iRecord < maxrec; iRecord++ ) {
715                 dtmp = DBFReadDoubleAttribute( hDBF, iRecord, i );
716                 if (dlow > dtmp) dlow = dtmp;
717                 if (dhigh < dtmp) dhigh = dtmp;
718                 dsum = dsum + dtmp;
719             }
720             mean=dsum/maxrec;
721             sprintf(stmp,"%%.%df to %%.%df \t(%%.%df)",iDecimals,iDecimals,iDecimals);
722             if (dlow < dhigh)       printf(stmp,dlow,dhigh,mean);
723             else if (dlow == dhigh) {
724                 sprintf(stmp,"= %%.%df",iDecimals);
725                 printf(stmp,dlow);
726             }
727             else printf("No Values");
728             break;
729
730           case FTInvalid:
731             break;
732
733         }
734
735     }
736     printf("\n");
737 }
738
739 int selectrec()
740 {
741     long int value, ty;
742
743     ty = DBFGetFieldInfo( hDBF, iselectitem, NULL, &iWidth, &iDecimals);
744     switch(ty)
745     {
746       case FTString:
747         puts("Invalid Item");
748         iselect=FALSE;
749         break;
750       case FTInteger:
751         value = DBFReadIntegerAttribute( hDBF, iRecord, iselectitem );
752         for (j = 0; j<selcount; j++)
753         {
754             if (selectvalues[j] == value)
755             {
756                 if (iunselect) return(0);  /* Keep this record */
757                 else  return(1);  /* Skip this record */
758             }
759         }
760         break;
761       case FTDouble:
762         puts("Invalid Item");
763         iselect=FALSE;
764         break;
765     }
766     if (iunselect) return(1);  /* Skip this record */
767     else  return(0);  /* Keep this record */
768 }
769
770
771 void check_theme_bnd()
772 {
773     if ( (adfBoundsMin[0] >= cxmin) && (adfBoundsMax[0] <= cxmax) &&
774          (adfBoundsMin[1] >= cymin) && (adfBoundsMax[1] <= cymax) )
775     {   /** Theme is totally inside clip area **/
776         if (ierase) nEntities=0; /** SKIP THEME  **/
777         else   iclip=FALSE; /** WRITE THEME (Clip not needed) **/
778     }
779             
780     if ( ( (adfBoundsMin[0] < cxmin) && (adfBoundsMax[0] < cxmin) ) ||
781          ( (adfBoundsMin[1] < cymin) && (adfBoundsMax[1] < cymin) ) ||
782          ( (adfBoundsMin[0] > cxmax) && (adfBoundsMax[0] > cxmax) ) ||
783          ( (adfBoundsMin[1] > cymax) && (adfBoundsMax[1] > cymax) ) )
784     {   /** Theme is totally outside clip area **/
785         if (ierase) iclip=FALSE; /** WRITE THEME (Clip not needed) **/
786              else   nEntities=0; /** SKIP THEME  **/
787     }
788             
789     if (nEntities == 0)
790         puts("WARNING: Theme is outside the clip area."); /** SKIP THEME  **/
791 }
792
793 int clip_boundary()
794 {
795     int  inside;
796     int  prev_outside;
797     int  i2;
798     int  j2;
799     
800     /*** FIRST check the boundary of the feature ***/
801     if ( ( (psCShape->dfXMin < cxmin) && (psCShape->dfXMax < cxmin) ) ||
802          ( (psCShape->dfYMin < cymin) && (psCShape->dfYMax < cymin) ) ||
803          ( (psCShape->dfXMin > cxmax) && (psCShape->dfXMax > cxmax) ) ||
804          ( (psCShape->dfYMin > cymax) && (psCShape->dfYMax > cymax) ) )
805     {   /** Feature is totally outside clip area **/
806         if (ierase) return(1); /** WRITE RECORD **/
807         else   return(0); /** SKIP  RECORD **/
808     }
809        
810     if ( (psCShape->dfXMin >= cxmin) && (psCShape->dfXMax <= cxmax) &&
811          (psCShape->dfYMin >= cymin) && (psCShape->dfYMax <= cymax) )
812     {   /** Feature is totally inside clip area **/
813         if (ierase) return(0); /** SKIP  RECORD **/
814         else   return(1); /** WRITE RECORD **/
815     }
816             
817     if (iinside) 
818     { /** INSIDE * Feature might touch the boundary or could be outside **/
819         if (ierase)  return(1); /** WRITE RECORD **/
820         else    return(0); /** SKIP  RECORD **/
821     }
822        
823     if (itouch)
824     {   /** TOUCH **/
825         if ( ( (psCShape->dfXMin <= cxmin) || (psCShape->dfXMax >= cxmax) ) && 
826              (psCShape->dfYMin >= cymin) && (psCShape->dfYMax <= cymax)    )
827         {   /** Feature intersects the clip boundary only on the X axis **/
828             if (ierase) return(0); /** SKIP  RECORD **/
829             else   return(1); /** WRITE RECORD **/
830         }
831
832         if (   (psCShape->dfXMin >= cxmin) && (psCShape->dfXMax <= cxmax)   && 
833                ( (psCShape->dfYMin <= cymin) || (psCShape->dfYMax >= cymax) )  )
834         {   /** Feature intersects the clip boundary only on the Y axis **/
835             if (ierase) return(0); /** SKIP  RECORD **/
836             else   return(1); /** WRITE RECORD **/
837         }
838                
839         for( j2 = 0; j2 < psCShape->nVertices; j2++ ) 
840         {   /** At least one vertex must be inside the clip boundary **/
841             if ( (psCShape->padfX[j2] >= cxmin  &&  psCShape->padfX[j2] <= cxmax) ||
842                  (psCShape->padfY[j2] >= cymin  &&  psCShape->padfY[j2] <= cymax)  )
843             {
844                 if (ierase) return(0); /** SKIP  RECORD **/
845                 else   return(1); /** WRITE RECORD **/
846             }
847         }
848                
849         /** All vertices are outside the clip boundary **/ 
850         if (ierase) return(1); /** WRITE RECORD **/
851         else   return(0); /** SKIP  RECORD **/
852     }   /** End TOUCH **/
853           
854     if (icut)
855     {   /** CUT **/
856         /*** Check each vertex in the feature with the Boundary and "CUT" ***/
857         /*** THIS CODE WAS NOT COMPLETED!  READ NOTE AT THE BOTTOM ***/
858         i2=0;
859         prev_outside=FALSE;
860         for( j2 = 0; j2 < psCShape->nVertices; j2++ ) 
861         {
862             inside = psCShape->padfX[j2] >= cxmin  &&  psCShape->padfX[j2] <= cxmax  &&
863                 psCShape->padfY[j2] >= cymin  &&  psCShape->padfY[j2] <= cymax ;
864                       
865             if (ierase) inside=(! inside);
866             if (inside)
867             {
868                 if (i2 != j2)
869                 {
870                     if (prev_outside)
871                     {
872                         /*** AddIntersection(i2); ***/  /*** Add intersection ***/
873                         prev_outside=FALSE;
874                     }
875                     psCShape->padfX[i2]=psCShape->padfX[j2];     /** move vertex **/
876                     psCShape->padfY[i2]=psCShape->padfY[j2];
877                 }
878                 i2++;
879             } else {
880                 if ( (! prev_outside) && (j2 > 0) )
881                 {
882                     /*** AddIntersection(i2); ***//*** Add intersection (Watch out for j2==i2-1) ***/
883                     /*** Also a polygon may overlap twice and will split into a several parts  ***/
884                     prev_outside=TRUE;
885                 }
886             }
887         }
888              
889         printf("Vertices:%d   OUT:%d   Number of Parts:%d\n",
890                psCShape->nVertices,i2, psCShape->nParts );
891                
892         psCShape->nVertices = i2;
893              
894         if (i2 < 2) return(0); /** SKIP RECORD **/
895         /*** (WE ARE NOT CREATING INTERESECTIONS and some lines could be reduced to one point) **/
896         
897         if (i2 == 0) return(0); /** SKIP  RECORD **/
898         else    return(1); /** WRITE RECORD **/
899     }  /** End CUT **/
900 }
901
902
903 /************************************************************************/
904 /*                            strncasecmp2()                            */
905 /*                                                                      */
906 /*      Compare two strings up to n characters                          */
907 /*      If n=0 then s1 and s2 must be an exact match                    */
908 /************************************************************************/
909
910 int strncasecmp2(char *s1, char *s2, int n)
911
912 {
913 int j,i;
914    if (n<1) n=strlen(s1)+1;
915    for (i=0; i<n; i++)
916    {
917       if (*s1 != *s2)
918       {
919          if (*s1 >= 'a' && *s1 <= 'z') {
920             j=*s1-32;
921             if (j != *s2) return(*s1-*s2);
922          } else {
923             if (*s1 >= 'A' && *s1 <= 'Z') { j=*s1+32; }
924                                    else   { j=*s1;    }
925             if (j != *s2) return(*s1-*s2); 
926          }
927       }
928       s1++;
929       s2++;
930    }
931    return(0);
932 }
933
934
935 #define  NKEYS (sizeof(unitkeytab) / sizeof(struct unitkey))
936 double findunit(char *unit)
937    {
938    struct unitkey {
939      char   *name;
940      double value;
941    } unitkeytab[] = {
942      "CM",            39.37,
943      "CENTIMETER",    39.37,
944      "CENTIMETERS",   39.37,  /** # of inches * 100 in unit **/
945      "METER",          3937,
946      "METERS",         3937,
947      "KM",          3937000,
948      "KILOMETER",   3937000, 
949      "KILOMETERS",  3937000,
950      "INCH",            100,
951      "INCHES",          100,
952      "FEET",           1200,
953      "FOOT",           1200,
954      "YARD",           3600,
955      "YARDS",          3600,       
956      "MILE",        6336000,
957      "MILES",       6336000  
958    };
959
960    double unitfactor=0;
961    for (j = 0; j < NKEYS; j++) {
962     if (strncasecmp2(unit, unitkeytab[j].name, 0) == 0) unitfactor=unitkeytab[j].value;
963    }
964    return(unitfactor);
965 }
966
967 /* -------------------------------------------------------------------- */
968 /*      Display a usage message.                                        */
969 /* -------------------------------------------------------------------- */
970 void error()
971 {       
972     puts( "The program will append to an existing shape file or it will" );
973     puts( "create a new file if needed." );
974     puts( "Only the items in the first output file will be preserved." );
975     puts( "When an item does not match with the append theme then the item");
976     puts( "might be placed to an existing item at the same position and type." );
977     puts( "  OTHER FUNCTIONS:" );
978     puts( "  - Describe all items in the dbase file (Use ALL for more than 5000 recs.)");
979     puts( "  - Select a group of shapes from a comma separated selection list.");
980     puts( "  - UnSelect a group of shapes from a comma separated selection list.");
981     puts( "  - Clip boundary extent or by theme boundary." );
982     puts( "      Touch writes all the shapes that touch the boundary.");
983     puts( "      Inside writes all the shapes that are completely within the boundary.");
984     puts( "      Boundary clips are only the min and max of a theme boundary." );
985     puts( "  - Erase boundary extent or by theme boundary." );
986     puts( "      Erase is the direct opposite of the Clip function." );
987     puts( "  - Change coordinate value units between meters and feet.");
988     puts( "      There is no way to determine the input unit of a shape file.");
989     puts( "      Skip this function if the shape file is already in the correct unit.");
990     puts( "      Clip and Erase will be done before the unit is changed.");
991     puts( "      A shift will be done after the unit is changed."); 
992     puts( "  - Shift X and Y coordinates.\n" );
993     puts( "Finally, There can only be one select or unselect in the command line.");
994     puts( "         There can only be one clip or erase in the command line.");
995     puts( "         There can only be one unit and only one shift in the command line.\n");
996     puts( "Ex: shputils in.shp out.shp   SELECT countycode 3,5,9,13,17,27");
997     puts( "    shputils in.shp out.shp   CLIP   10 10 90 90 Touch   FACTOR Meter Feet");
998     puts( "    shputils in.shp out.shp   FACTOR Meter 3.0");
999     puts( "    shputils in.shp out.shp   CLIP   clip.shp Boundary Touch   SHIFT 40 40");
1000     puts( "    shputils in.shp out.shp   SELECT co 112   CLIP clip.shp Boundary Touch\n");
1001     puts( "USAGE: shputils  <DescribeShape>   {ALL}");
1002     puts( "USAGE: shputils  <InputShape>  <OutShape|AppendShape>" );
1003     puts( "   { <FACTOR>       <FEET|MILES|METERS|KM> <FEET|MILES|METERS|KM|factor> }" );
1004     puts( "   { <SHIFT>        <xshift> <yshift> }" );
1005     puts( "   { <SELECT|UNSEL> <Item> <valuelist> }" );
1006     puts( "   { <CLIP|ERASE>   <xmin> <ymin> <xmax> <ymax> <TOUCH|INSIDE|CUT> }" );
1007     puts( "   { <CLIP|ERASE>   <theme>      <BOUNDARY>     <TOUCH|INSIDE|CUT> }" );
1008     puts( "     Note: CUT is not complete and does not create intersections.");
1009     puts( "           For more information read programmer comment.");
1010         
1011     /****   Clip functions for Polygon and Cut is not supported
1012             There are several web pages that describe methods of doing this function.
1013             It seem easy to impliment until you start writting code.  I don't have the
1014             time to add these functions but a did leave a simple cut routine in the 
1015             program that can be called by using CUT instead of TOUCH in the 
1016             CLIP or ERASE functions.  It does not add the intersection of the line and
1017             the clip box, so polygons could look incomplete and lines will come up short.
1018         
1019             Information about clipping lines with a box:
1020             http://www.csclub.uwaterloo.ca/u/mpslager/articles/sutherland/wr.html
1021             Information about finding the intersection of two lines:
1022             http://www.whisqu.se/per/docs/math28.htm
1023            
1024             THE CODE LOOKS LIKE THIS:
1025             ********************************************************      
1026             void Intersect_Lines(float x0,float y0,float x1,float y1,
1027             float x2,float y2,float x3,float y3,
1028             float *xi,float *yi)
1029             {
1030 //  this function computes the intersection of the sent lines
1031 //  and returns the intersection point, note that the function assumes
1032 //  the lines intersect. the function can handle vertical as well
1033 //  as horizontal lines. note the function isn't very clever, it simply
1034 //  applies the math, but we don't need speed since this is a
1035 //  pre-processing step
1036 //  The Intersect_lines program came from (http://www.whisqu.se/per/docs/math28.htm)
1037
1038 float a1,b1,c1, // constants of linear equations 
1039 a2,b2,c2,
1040 det_inv,  // the inverse of the determinant of the coefficientmatrix
1041 m1,m2;    // the slopes of each line
1042       
1043 // compute slopes, note the cludge for infinity, however, this will
1044 // be close enough
1045 if ((x1-x0)!=0)
1046 m1 = (y1-y0)/(x1-x0);
1047 else
1048 m1 = (float)1e+10;  // close enough to infinity
1049    
1050    
1051 if ((x3-x2)!=0) 
1052 m2 = (y3-y2)/(x3-x2);
1053 else
1054 m2 = (float)1e+10;  // close enough to infinity
1055    
1056 // compute constants
1057 a1 = m1;
1058 a2 = m2;
1059 b1 = -1;
1060 b2 = -1;
1061 c1 = (y0-m1*x0);
1062 c2 = (y2-m2*x2);
1063 // compute the inverse of the determinate
1064 det_inv = 1/(a1*b2 - a2*b1);
1065 // use Kramers rule to compute xi and yi
1066 *xi=((b1*c2 - b2*c1)*det_inv);
1067 *yi=((a2*c1 - a1*c2)*det_inv);
1068 } // end Intersect_Lines
1069     **********************************************************/
1070
1071     exit( 1 );
1072 }