1 : /******************************************************************************
2 : * $Id: gdaltindex.c 18344 2009-12-19 10:42:32Z rouault $
3 : *
4 : * Project: MapServer
5 : * Purpose: Commandline App to build tile index for raster files.
6 : * Author: Frank Warmerdam, warmerdam@pobox.com
7 : *
8 : ******************************************************************************
9 : * Copyright (c) 2001, Frank Warmerdam, DM Solutions Group Inc
10 : *
11 : * Permission is hereby granted, free of charge, to any person obtaining a
12 : * copy of this software and associated documentation files (the "Software"),
13 : * to deal in the Software without restriction, including without limitation
14 : * the rights to use, copy, modify, merge, publish, distribute, sublicense,
15 : * and/or sell copies of the Software, and to permit persons to whom the
16 : * Software is furnished to do so, subject to the following conditions:
17 : *
18 : * The above copyright notice and this permission notice shall be included
19 : * in all copies or substantial portions of the Software.
20 : *
21 : * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
22 : * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
23 : * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
24 : * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
25 : * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
26 : * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
27 : * DEALINGS IN THE SOFTWARE.
28 : ****************************************************************************/
29 :
30 : #include "ogr_api.h"
31 : #include "ogr_srs_api.h"
32 : #include "gdal.h"
33 : #include "cpl_port.h"
34 : #include "cpl_conv.h"
35 :
36 : CPL_CVSID("$Id: gdaltindex.c 18344 2009-12-19 10:42:32Z rouault $");
37 :
38 : /************************************************************************/
39 : /* Usage() */
40 : /************************************************************************/
41 :
42 0 : static void Usage()
43 :
44 : {
45 0 : fprintf(stdout, "%s",
46 : "\n"
47 : "Usage: gdaltindex [-tileindex field_name] [-write_absolute_path] \n"
48 : " [-skip_different_projection] index_file [gdal_file]*\n"
49 : "\n"
50 : "eg.\n"
51 : " % gdaltindex doq_index.shp doq/*.tif\n"
52 : "\n"
53 : "NOTES:\n"
54 : " o The shapefile (index_file) will be created if it doesn't already exist.\n"
55 : " o The default tile index field is 'location'.\n"
56 : " o Raster filenames will be put in the file exactly as they are specified\n"
57 : " on the commandline unless the option -write_absolute_path is used.\n"
58 : " o If -skip_different_projection is specified, only files with same projection ref\n"
59 : " as files already inserted in the tileindex will be inserted.\n"
60 : " o Simple rectangular polygons are generated in the same\n"
61 : " coordinate system as the rasters.\n" );
62 0 : exit(1);
63 : }
64 :
65 : /************************************************************************/
66 : /* main() */
67 : /************************************************************************/
68 :
69 6 : int main(int argc, char *argv[])
70 : {
71 6 : const char *index_filename = NULL;
72 6 : const char *tile_index = "location";
73 : int i_arg, ti_field;
74 : OGRDataSourceH hTileIndexDS;
75 6 : OGRLayerH hLayer = NULL;
76 : OGRFeatureDefnH hFDefn;
77 6 : int write_absolute_path = FALSE;
78 6 : char* current_path = NULL;
79 : int i;
80 : int nExistingFiles;
81 6 : int skip_different_projection = FALSE;
82 6 : char** existingFilesTab = NULL;
83 6 : int alreadyExistingProjectionRefValid = FALSE;
84 6 : char* alreadyExistingProjectionRef = NULL;
85 : char* index_filename_mod;
86 : int bExists;
87 : VSIStatBuf sStatBuf;
88 :
89 : /* Check that we are running against at least GDAL 1.4 */
90 : /* Note to developers : if we use newer API, please change the requirement */
91 6 : if (atoi(GDALVersionInfo("VERSION_NUM")) < 1400)
92 : {
93 0 : fprintf(stderr, "At least, GDAL >= 1.4.0 is required for this version of %s, "
94 : "which was compiled against GDAL %s\n", argv[0], GDAL_RELEASE_NAME);
95 0 : exit(1);
96 : }
97 :
98 6 : GDALAllRegister();
99 6 : OGRRegisterAll();
100 :
101 6 : argc = GDALGeneralCmdLineProcessor( argc, &argv, 0 );
102 6 : if( argc < 1 )
103 0 : exit( -argc );
104 :
105 : /* -------------------------------------------------------------------- */
106 : /* Get commandline arguments other than the GDAL raster filenames. */
107 : /* -------------------------------------------------------------------- */
108 7 : for( i_arg = 1; i_arg < argc; i_arg++ )
109 : {
110 7 : if( EQUAL(argv[i_arg], "--utility_version") )
111 : {
112 1 : printf("%s was compiled against GDAL %s and is running against GDAL %s\n",
113 : argv[0], GDAL_RELEASE_NAME, GDALVersionInfo("RELEASE_NAME"));
114 1 : return 0;
115 : }
116 6 : else if( strcmp(argv[i_arg],"-tileindex") == 0 )
117 : {
118 0 : tile_index = argv[++i_arg];
119 : }
120 6 : else if ( strcmp(argv[i_arg],"-write_absolute_path") == 0 )
121 : {
122 0 : write_absolute_path = TRUE;
123 : }
124 6 : else if ( strcmp(argv[i_arg],"-skip_different_projection") == 0 )
125 : {
126 1 : skip_different_projection = TRUE;
127 : }
128 5 : else if( argv[i_arg][0] == '-' )
129 0 : Usage();
130 5 : else if( index_filename == NULL )
131 : {
132 5 : index_filename = argv[i_arg];
133 5 : i_arg++;
134 5 : break;
135 : }
136 : }
137 :
138 5 : if( index_filename == NULL || i_arg == argc )
139 0 : Usage();
140 :
141 : /* -------------------------------------------------------------------- */
142 : /* Open or create the target shapefile and DBF file. */
143 : /* -------------------------------------------------------------------- */
144 5 : index_filename_mod = CPLStrdup(CPLResetExtension(index_filename, "shp"));
145 :
146 5 : bExists = (VSIStat(index_filename_mod, &sStatBuf) == 0);
147 5 : if (!bExists)
148 : {
149 2 : CPLFree(index_filename_mod);
150 2 : index_filename_mod = CPLStrdup(CPLResetExtension(index_filename, "SHP"));
151 2 : bExists = (VSIStat(index_filename_mod, &sStatBuf) == 0);
152 : }
153 5 : CPLFree(index_filename_mod);
154 :
155 5 : if (bExists)
156 : {
157 3 : hTileIndexDS = OGROpen( index_filename, TRUE, NULL );
158 3 : if (hTileIndexDS != NULL)
159 : {
160 3 : hLayer = OGR_DS_GetLayer(hTileIndexDS, 0);
161 : }
162 : }
163 : else
164 : {
165 : OGRSFDriverH hDriver;
166 2 : const char* pszDriverName = "ESRI Shapefile";
167 :
168 2 : printf( "Creating new index file...\n" );
169 2 : hDriver = OGRGetDriverByName( pszDriverName );
170 2 : if( hDriver == NULL )
171 : {
172 0 : printf( "%s driver not available.\n", pszDriverName );
173 0 : exit( 1 );
174 : }
175 :
176 2 : hTileIndexDS = OGR_Dr_CreateDataSource( hDriver, index_filename, NULL );
177 2 : if (hTileIndexDS)
178 : {
179 2 : char* pszLayerName = CPLStrdup(CPLGetBasename(index_filename));
180 2 : OGRSpatialReferenceH hSpatialRef = NULL;
181 2 : GDALDatasetH hDS = GDALOpen( argv[i_arg], GA_ReadOnly );
182 2 : if (hDS)
183 : {
184 2 : const char* pszWKT = GDALGetProjectionRef(hDS);
185 2 : if (pszWKT != NULL && pszWKT[0] != '\0')
186 : {
187 2 : hSpatialRef = OSRNewSpatialReference(pszWKT);
188 : }
189 2 : GDALClose(hDS);
190 : }
191 :
192 2 : hLayer = OGR_DS_CreateLayer( hTileIndexDS, pszLayerName, hSpatialRef, wkbPolygon, NULL );
193 2 : CPLFree(pszLayerName);
194 2 : if (hSpatialRef)
195 2 : OSRRelease(hSpatialRef);
196 :
197 2 : if (hLayer)
198 : {
199 2 : OGRFieldDefnH hFieldDefn = OGR_Fld_Create( tile_index, OFTString );
200 2 : OGR_Fld_SetWidth( hFieldDefn, 255);
201 2 : OGR_L_CreateField( hLayer, hFieldDefn, TRUE );
202 2 : OGR_Fld_Destroy(hFieldDefn);
203 : }
204 : }
205 : }
206 :
207 5 : if( hTileIndexDS == NULL || hLayer == NULL )
208 : {
209 0 : fprintf( stderr, "Unable to open/create shapefile `%s'.\n",
210 : index_filename );
211 0 : exit(2);
212 : }
213 :
214 5 : hFDefn = OGR_L_GetLayerDefn(hLayer);
215 :
216 5 : for( ti_field = 0; ti_field < OGR_FD_GetFieldCount(hFDefn); ti_field++ )
217 : {
218 5 : OGRFieldDefnH hFieldDefn = OGR_FD_GetFieldDefn( hFDefn, ti_field );
219 5 : if( strcmp(OGR_Fld_GetNameRef(hFieldDefn), tile_index) == 0 )
220 5 : break;
221 : }
222 :
223 5 : if( ti_field == OGR_FD_GetFieldCount(hFDefn) )
224 : {
225 0 : fprintf( stderr, "Unable to find field `%s' in DBF file `%s'.\n",
226 : tile_index, index_filename );
227 0 : exit(2);
228 : }
229 :
230 : /* Load in memory existing file names in SHP */
231 5 : nExistingFiles = OGR_L_GetFeatureCount(hLayer, FALSE);
232 5 : if (nExistingFiles)
233 : {
234 : OGRFeatureH hFeature;
235 3 : existingFilesTab = (char**)CPLMalloc(nExistingFiles * sizeof(char*));
236 13 : for(i=0;i<nExistingFiles;i++)
237 : {
238 10 : hFeature = OGR_L_GetNextFeature(hLayer);
239 10 : existingFilesTab[i] = CPLStrdup(OGR_F_GetFieldAsString( hFeature, ti_field ));
240 10 : if (i == 0)
241 : {
242 3 : GDALDatasetH hDS = GDALOpen(existingFilesTab[i], GA_ReadOnly );
243 3 : if (hDS)
244 : {
245 3 : alreadyExistingProjectionRefValid = TRUE;
246 3 : alreadyExistingProjectionRef = CPLStrdup(GDALGetProjectionRef(hDS));
247 3 : GDALClose(hDS);
248 : }
249 : }
250 10 : OGR_F_Destroy( hFeature );
251 : }
252 : }
253 :
254 5 : if (write_absolute_path)
255 : {
256 0 : current_path = CPLGetCurrentDir();
257 0 : if (current_path == NULL)
258 : {
259 0 : fprintf( stderr, "This system does not support the CPLGetCurrentDir call. "
260 : "The option -write_absolute_path will have no effect\n");
261 0 : write_absolute_path = FALSE;
262 : }
263 : }
264 :
265 : /* -------------------------------------------------------------------- */
266 : /* loop over GDAL files, processing. */
267 : /* -------------------------------------------------------------------- */
268 18 : for( ; i_arg < argc; i_arg++ )
269 : {
270 : GDALDatasetH hDS;
271 : double adfGeoTransform[6];
272 : double adfX[5], adfY[5];
273 : int nXSize, nYSize;
274 : char* fileNameToWrite;
275 : const char* projectionRef;
276 : VSIStatBuf sStatBuf;
277 : int k;
278 : OGRFeatureH hFeature;
279 : OGRGeometryH hPoly, hRing;
280 :
281 : /* Make sure it is a file before building absolute path name */
282 13 : if (write_absolute_path && CPLIsFilenameRelative( argv[i_arg] ) &&
283 0 : VSIStat( argv[i_arg], &sStatBuf ) == 0)
284 : {
285 0 : fileNameToWrite = CPLStrdup(CPLProjectRelativeFilename(current_path, argv[i_arg]));
286 : }
287 : else
288 : {
289 13 : fileNameToWrite = CPLStrdup(argv[i_arg]);
290 : }
291 :
292 : /* Checks that file is not already in tileindex */
293 27 : for(i=0;i<nExistingFiles;i++)
294 : {
295 18 : if (EQUAL(fileNameToWrite, existingFilesTab[i]))
296 : {
297 4 : fprintf(stderr, "File %s is already in tileindex. Skipping it.\n",
298 : fileNameToWrite);
299 4 : break;
300 : }
301 : }
302 13 : if (i != nExistingFiles)
303 : {
304 4 : CPLFree(fileNameToWrite);
305 4 : continue;
306 : }
307 :
308 9 : hDS = GDALOpen( argv[i_arg], GA_ReadOnly );
309 9 : if( hDS == NULL )
310 : {
311 0 : fprintf( stderr, "Unable to open %s, skipping.\n",
312 0 : argv[i_arg] );
313 0 : CPLFree(fileNameToWrite);
314 0 : continue;
315 : }
316 :
317 9 : GDALGetGeoTransform( hDS, adfGeoTransform );
318 9 : if( adfGeoTransform[0] == 0.0
319 0 : && adfGeoTransform[1] == 1.0
320 0 : && adfGeoTransform[3] == 0.0
321 0 : && ABS(adfGeoTransform[5]) == 1.0 )
322 : {
323 0 : fprintf( stderr,
324 : "It appears no georeferencing is available for\n"
325 : "`%s', skipping.\n",
326 0 : argv[i_arg] );
327 0 : GDALClose( hDS );
328 0 : CPLFree(fileNameToWrite);
329 0 : continue;
330 : }
331 :
332 9 : projectionRef = GDALGetProjectionRef(hDS);
333 9 : if (alreadyExistingProjectionRefValid)
334 : {
335 : int projectionRefNotNull, alreadyExistingProjectionRefNotNull;
336 7 : projectionRefNotNull = projectionRef && projectionRef[0];
337 7 : alreadyExistingProjectionRefNotNull = alreadyExistingProjectionRef && alreadyExistingProjectionRef[0];
338 14 : if ((projectionRefNotNull &&
339 : alreadyExistingProjectionRefNotNull &&
340 7 : EQUAL(projectionRef, alreadyExistingProjectionRef) == 0) ||
341 : (projectionRefNotNull != alreadyExistingProjectionRefNotNull))
342 : {
343 2 : fprintf(stderr, "Warning : %s is not using the same projection system as "
344 : "other files in the tileindex. This may cause problems when "
345 1 : "using it in MapServer for example.%s\n", argv[i_arg],
346 : (skip_different_projection) ? " Skipping it" : "");
347 1 : if (skip_different_projection)
348 : {
349 1 : CPLFree(fileNameToWrite);
350 1 : GDALClose( hDS );
351 1 : continue;
352 : }
353 : }
354 : }
355 : else
356 : {
357 2 : alreadyExistingProjectionRefValid = TRUE;
358 2 : alreadyExistingProjectionRef = CPLStrdup(projectionRef);
359 : }
360 :
361 8 : nXSize = GDALGetRasterXSize( hDS );
362 8 : nYSize = GDALGetRasterYSize( hDS );
363 :
364 16 : adfX[0] = adfGeoTransform[0]
365 8 : + 0 * adfGeoTransform[1]
366 8 : + 0 * adfGeoTransform[2];
367 16 : adfY[0] = adfGeoTransform[3]
368 8 : + 0 * adfGeoTransform[4]
369 8 : + 0 * adfGeoTransform[5];
370 :
371 16 : adfX[1] = adfGeoTransform[0]
372 8 : + nXSize * adfGeoTransform[1]
373 16 : + 0 * adfGeoTransform[2];
374 16 : adfY[1] = adfGeoTransform[3]
375 8 : + nXSize * adfGeoTransform[4]
376 16 : + 0 * adfGeoTransform[5];
377 :
378 16 : adfX[2] = adfGeoTransform[0]
379 8 : + nXSize * adfGeoTransform[1]
380 16 : + nYSize * adfGeoTransform[2];
381 16 : adfY[2] = adfGeoTransform[3]
382 8 : + nXSize * adfGeoTransform[4]
383 16 : + nYSize * adfGeoTransform[5];
384 :
385 16 : adfX[3] = adfGeoTransform[0]
386 8 : + 0 * adfGeoTransform[1]
387 8 : + nYSize * adfGeoTransform[2];
388 16 : adfY[3] = adfGeoTransform[3]
389 8 : + 0 * adfGeoTransform[4]
390 8 : + nYSize * adfGeoTransform[5];
391 :
392 16 : adfX[4] = adfGeoTransform[0]
393 8 : + 0 * adfGeoTransform[1]
394 8 : + 0 * adfGeoTransform[2];
395 16 : adfY[4] = adfGeoTransform[3]
396 8 : + 0 * adfGeoTransform[4]
397 8 : + 0 * adfGeoTransform[5];
398 :
399 8 : hFeature = OGR_F_Create( OGR_L_GetLayerDefn( hLayer ) );
400 8 : OGR_F_SetFieldString( hFeature, ti_field, fileNameToWrite );
401 :
402 8 : hPoly = OGR_G_CreateGeometry(wkbPolygon);
403 8 : hRing = OGR_G_CreateGeometry(wkbLinearRing);
404 48 : for(k=0;k<5;k++)
405 40 : OGR_G_SetPoint_2D(hRing, k, adfX[k], adfY[k]);
406 8 : OGR_G_AddGeometryDirectly( hPoly, hRing );
407 8 : OGR_F_SetGeometryDirectly( hFeature, hPoly );
408 :
409 8 : if( OGR_L_CreateFeature( hLayer, hFeature ) != OGRERR_NONE )
410 : {
411 0 : printf( "Failed to create feature in shapefile.\n" );
412 0 : break;
413 : }
414 :
415 8 : OGR_F_Destroy( hFeature );
416 :
417 :
418 8 : CPLFree(fileNameToWrite);
419 :
420 8 : GDALClose( hDS );
421 : }
422 :
423 5 : CPLFree(current_path);
424 :
425 5 : if (nExistingFiles)
426 : {
427 13 : for(i=0;i<nExistingFiles;i++)
428 : {
429 10 : CPLFree(existingFilesTab[i]);
430 : }
431 3 : CPLFree(existingFilesTab);
432 : }
433 5 : CPLFree(alreadyExistingProjectionRef);
434 :
435 5 : OGR_DS_Destroy( hTileIndexDS );
436 :
437 5 : GDALDestroyDriverManager();
438 :
439 5 : exit( 0 );
440 : }
|
