1 : /******************************************************************************
2 : * $Id: rasterlitecreatecopy.cpp 22035 2011-03-25 23:57:59Z rouault $
3 : *
4 : * Project: GDAL Rasterlite driver
5 : * Purpose: Implement GDAL Rasterlite support using OGR SQLite driver
6 : * Author: Even Rouault, <even dot rouault at mines dash paris dot org>
7 : *
8 : **********************************************************************
9 : * Copyright (c) 2009, Even Rouault, <even dot rouault at mines dash paris dot org>
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 "cpl_string.h"
31 : #include "ogr_api.h"
32 : #include "ogr_srs_api.h"
33 :
34 : #include "rasterlitedataset.h"
35 :
36 : CPL_CVSID("$Id: rasterlitecreatecopy.cpp 22035 2011-03-25 23:57:59Z rouault $");
37 :
38 : /************************************************************************/
39 : /* RasterliteGetTileDriverOptions () */
40 : /************************************************************************/
41 :
42 144 : static char** RasterliteAddTileDriverOptionsForDriver(char** papszOptions,
43 : char** papszTileDriverOptions,
44 : const char* pszOptionName,
45 : const char* pszExpectedDriverName)
46 : {
47 144 : const char* pszVal = CSLFetchNameValue(papszOptions, pszOptionName);
48 144 : if (pszVal)
49 : {
50 : const char* pszDriverName =
51 0 : CSLFetchNameValueDef(papszOptions, "DRIVER", "GTiff");
52 0 : if (EQUAL(pszDriverName, pszExpectedDriverName))
53 : {
54 : papszTileDriverOptions =
55 0 : CSLSetNameValue(papszTileDriverOptions, pszOptionName, pszVal);
56 : }
57 : else
58 : {
59 : CPLError(CE_Warning, CPLE_NotSupported,
60 : "Unexpected option '%s' for driver '%s'",
61 0 : pszOptionName, pszDriverName);
62 : }
63 : }
64 144 : return papszTileDriverOptions;
65 : }
66 :
67 36 : static char** RasterliteGetTileDriverOptions(char** papszOptions)
68 : {
69 36 : char** papszTileDriverOptions = NULL;
70 :
71 : const char* pszDriverName =
72 36 : CSLFetchNameValueDef(papszOptions, "DRIVER", "GTiff");
73 :
74 36 : if (EQUAL(pszDriverName, "EPSILON"))
75 : {
76 : papszTileDriverOptions = CSLSetNameValue(papszTileDriverOptions,
77 4 : "RASTERLITE_OUTPUT", "YES");
78 : }
79 :
80 36 : const char* pszQuality = CSLFetchNameValue(papszOptions, "QUALITY");
81 36 : if (pszQuality)
82 : {
83 0 : if (EQUAL(pszDriverName, "GTiff"))
84 : {
85 : papszTileDriverOptions =
86 0 : CSLSetNameValue(papszTileDriverOptions, "JPEG_QUALITY", pszQuality);
87 : }
88 0 : else if (EQUAL(pszDriverName, "JPEG") || EQUAL(pszDriverName, "WEBP"))
89 : {
90 : papszTileDriverOptions =
91 0 : CSLSetNameValue(papszTileDriverOptions, "QUALITY", pszQuality);
92 : }
93 : else
94 : {
95 : CPLError(CE_Warning, CPLE_NotSupported,
96 : "Unexpected option '%s' for driver '%s'",
97 0 : "QUALITY", pszDriverName);
98 : }
99 : }
100 :
101 : papszTileDriverOptions = RasterliteAddTileDriverOptionsForDriver(
102 36 : papszOptions, papszTileDriverOptions, "COMPRESS", "GTiff");
103 : papszTileDriverOptions = RasterliteAddTileDriverOptionsForDriver(
104 36 : papszOptions, papszTileDriverOptions, "PHOTOMETRIC", "GTiff");
105 : papszTileDriverOptions = RasterliteAddTileDriverOptionsForDriver(
106 36 : papszOptions, papszTileDriverOptions, "TARGET", "EPSILON");
107 : papszTileDriverOptions = RasterliteAddTileDriverOptionsForDriver(
108 36 : papszOptions, papszTileDriverOptions, "FILTER", "EPSILON");
109 :
110 36 : return papszTileDriverOptions;
111 : }
112 :
113 : /************************************************************************/
114 : /* RasterliteInsertSRID () */
115 : /************************************************************************/
116 :
117 36 : static int RasterliteInsertSRID(OGRDataSourceH hDS, const char* pszWKT)
118 : {
119 36 : CPLString osSQL;
120 :
121 36 : int nAuthorityCode = 0;
122 36 : CPLString osAuthorityName, osProjCS, osProj4;
123 36 : if (pszWKT != NULL && strlen(pszWKT) != 0)
124 : {
125 36 : OGRSpatialReferenceH hSRS = OSRNewSpatialReference(pszWKT);
126 36 : if (hSRS)
127 : {
128 36 : const char* pszAuthorityName = OSRGetAuthorityName(hSRS, NULL);
129 36 : if (pszAuthorityName) osAuthorityName = pszAuthorityName;
130 :
131 36 : const char* pszProjCS = OSRGetAttrValue(hSRS, "PROJCS", 0);
132 36 : if (pszProjCS) osProjCS = pszProjCS;
133 :
134 36 : const char* pszAuthorityCode = OSRGetAuthorityCode(hSRS, NULL);
135 36 : if (pszAuthorityCode) nAuthorityCode = atoi(pszAuthorityCode);
136 :
137 36 : char *pszProj4 = NULL;
138 36 : if( OSRExportToProj4( hSRS, &pszProj4 ) != OGRERR_NONE )
139 0 : pszProj4 = CPLStrdup("");
140 36 : osProj4 = pszProj4;
141 36 : CPLFree(pszProj4);
142 : }
143 36 : OSRDestroySpatialReference(hSRS);
144 : }
145 :
146 36 : int nSRSId = -1;
147 36 : if (nAuthorityCode != 0 && osAuthorityName.size() != 0)
148 : {
149 36 : osSQL.Printf ("SELECT srid FROM spatial_ref_sys WHERE auth_srid = %d", nAuthorityCode);
150 36 : OGRLayerH hLyr = OGR_DS_ExecuteSQL(hDS, osSQL.c_str(), NULL, NULL);
151 36 : if (hLyr == NULL)
152 : {
153 0 : nSRSId = nAuthorityCode;
154 :
155 0 : if ( osProjCS.size() != 0 )
156 : osSQL.Printf(
157 : "INSERT INTO spatial_ref_sys "
158 : "(srid, auth_name, auth_srid, ref_sys_name, proj4text) "
159 : "VALUES (%d, '%s', '%d', '%s', '%s')",
160 : nSRSId, osAuthorityName.c_str(),
161 0 : nAuthorityCode, osProjCS.c_str(), osProj4.c_str() );
162 : else
163 : osSQL.Printf(
164 : "INSERT INTO spatial_ref_sys "
165 : "(srid, auth_name, auth_srid, proj4text) "
166 : "VALUES (%d, '%s', '%d', '%s')",
167 : nSRSId, osAuthorityName.c_str(),
168 0 : nAuthorityCode, osProj4.c_str() );
169 :
170 :
171 0 : OGR_DS_ExecuteSQL(hDS, osSQL.c_str(), NULL, NULL);
172 : }
173 : else
174 : {
175 36 : OGRFeatureH hFeat = OGR_L_GetNextFeature(hLyr);
176 36 : if (hFeat)
177 : {
178 36 : nSRSId = OGR_F_GetFieldAsInteger(hFeat, 0);
179 36 : OGR_F_Destroy(hFeat);
180 : }
181 36 : OGR_DS_ReleaseResultSet(hDS, hLyr);
182 : }
183 : }
184 :
185 36 : return nSRSId;
186 : }
187 :
188 : /************************************************************************/
189 : /* RasterliteCreateTables () */
190 : /************************************************************************/
191 :
192 36 : OGRDataSourceH RasterliteCreateTables(OGRDataSourceH hDS, const char* pszTableName,
193 : int nSRSId, int bWipeExistingData)
194 : {
195 36 : CPLString osSQL;
196 :
197 36 : CPLString osOldVal = CPLGetConfigOption("SQLITE_LIST_ALL_TABLES", "FALSE");
198 36 : CPLString osDBName = OGR_DS_GetName(hDS);
199 :
200 36 : CPLString osRasterLayer;
201 36 : osRasterLayer.Printf("%s_rasters", pszTableName);
202 :
203 36 : CPLString osMetatadataLayer;
204 36 : osMetatadataLayer.Printf("%s_metadata", pszTableName);
205 :
206 : OGRLayerH hLyr;
207 :
208 36 : if (OGR_DS_GetLayerByName(hDS, osRasterLayer.c_str()) == NULL)
209 : {
210 : /* -------------------------------------------------------------------- */
211 : /* The table don't exist. Create them */
212 : /* -------------------------------------------------------------------- */
213 :
214 : /* Create _rasters table */
215 : osSQL.Printf ("CREATE TABLE \"%s\" ("
216 : "id INTEGER NOT NULL PRIMARY KEY AUTOINCREMENT,"
217 36 : "raster BLOB NOT NULL)", osRasterLayer.c_str());
218 36 : OGR_DS_ExecuteSQL(hDS, osSQL.c_str(), NULL, NULL);
219 :
220 : /* Create _metadata table */
221 : osSQL.Printf ("CREATE TABLE \"%s\" ("
222 : "id INTEGER NOT NULL PRIMARY KEY,"
223 : "source_name TEXT NOT NULL,"
224 : "tile_id INTEGER NOT NULL,"
225 : "width INTEGER NOT NULL,"
226 : "height INTEGER NOT NULL,"
227 : "pixel_x_size DOUBLE NOT NULL,"
228 : "pixel_y_size DOUBLE NOT NULL)",
229 36 : osMetatadataLayer.c_str());
230 36 : OGR_DS_ExecuteSQL(hDS, osSQL.c_str(), NULL, NULL);
231 :
232 : /* Add geometry column to _metadata table */
233 : osSQL.Printf("SELECT AddGeometryColumn('%s', 'geometry', %d, 'POLYGON', 2)",
234 36 : osMetatadataLayer.c_str(), nSRSId);
235 36 : if ((hLyr = OGR_DS_ExecuteSQL(hDS, osSQL.c_str(), NULL, NULL)) == NULL)
236 : {
237 : CPLError(CE_Failure, CPLE_AppDefined,
238 0 : "Check that the OGR SQLite driver has Spatialite support");
239 0 : OGRReleaseDataSource(hDS);
240 0 : return NULL;
241 : }
242 36 : OGR_DS_ReleaseResultSet(hDS, hLyr);
243 :
244 : /* Create spatial index on _metadata table */
245 : osSQL.Printf("SELECT CreateSpatialIndex('%s', 'geometry')",
246 36 : osMetatadataLayer.c_str());
247 36 : if ((hLyr = OGR_DS_ExecuteSQL(hDS, osSQL.c_str(), NULL, NULL)) == NULL)
248 : {
249 0 : OGRReleaseDataSource(hDS);
250 0 : return NULL;
251 : }
252 36 : OGR_DS_ReleaseResultSet(hDS, hLyr);
253 :
254 : /* Re-open the DB to take into account the new tables*/
255 36 : OGRReleaseDataSource(hDS);
256 :
257 36 : CPLSetThreadLocalConfigOption("SQLITE_LIST_ALL_TABLES", "TRUE");
258 36 : hDS = OGROpen(osDBName.c_str(), TRUE, NULL);
259 36 : CPLSetThreadLocalConfigOption("SQLITE_LIST_ALL_TABLES", osOldVal.c_str());
260 : }
261 : else
262 : {
263 : /* Check that the existing SRS is consistent with the one of the new */
264 : /* data to be inserted */
265 : osSQL.Printf("SELECT srid FROM geometry_columns WHERE f_table_name = '%s'",
266 0 : osMetatadataLayer.c_str());
267 0 : hLyr = OGR_DS_ExecuteSQL(hDS, osSQL.c_str(), NULL, NULL);
268 0 : if (hLyr)
269 : {
270 0 : int nExistingSRID = -1;
271 0 : OGRFeatureH hFeat = OGR_L_GetNextFeature(hLyr);
272 0 : if (hFeat)
273 : {
274 0 : nExistingSRID = OGR_F_GetFieldAsInteger(hFeat, 0);
275 0 : OGR_F_Destroy(hFeat);
276 : }
277 0 : OGR_DS_ReleaseResultSet(hDS, hLyr);
278 :
279 0 : if (nExistingSRID != nSRSId)
280 : {
281 0 : if (bWipeExistingData)
282 : {
283 : osSQL.Printf("UPDATE geometry_columns SET srid = %d "
284 : "WHERE f_table_name = \"%s\"",
285 0 : nSRSId, osMetatadataLayer.c_str());
286 0 : OGR_DS_ExecuteSQL(hDS, osSQL.c_str(), NULL, NULL);
287 :
288 : /* Re-open the DB to take into account the change of SRS */
289 0 : OGRReleaseDataSource(hDS);
290 :
291 0 : CPLSetThreadLocalConfigOption("SQLITE_LIST_ALL_TABLES", "TRUE");
292 0 : hDS = OGROpen(osDBName.c_str(), TRUE, NULL);
293 0 : CPLSetThreadLocalConfigOption("SQLITE_LIST_ALL_TABLES", osOldVal.c_str());
294 : }
295 : else
296 : {
297 : CPLError(CE_Failure, CPLE_NotSupported,
298 0 : "New data has not the same SRS as existing data");
299 0 : OGRReleaseDataSource(hDS);
300 0 : return NULL;
301 : }
302 : }
303 : }
304 :
305 0 : if (bWipeExistingData)
306 : {
307 0 : osSQL.Printf("DELETE FROM \"%s\"", osRasterLayer.c_str());
308 0 : OGR_DS_ExecuteSQL(hDS, osSQL.c_str(), NULL, NULL);
309 :
310 0 : osSQL.Printf("DELETE FROM \"%s\"", osMetatadataLayer.c_str());
311 0 : OGR_DS_ExecuteSQL(hDS, osSQL.c_str(), NULL, NULL);
312 : }
313 : }
314 :
315 36 : return hDS;
316 : }
317 :
318 : /************************************************************************/
319 : /* RasterliteCreateCopy () */
320 : /************************************************************************/
321 :
322 : GDALDataset *
323 42 : RasterliteCreateCopy( const char * pszFilename, GDALDataset *poSrcDS,
324 : int bStrict, char ** papszOptions,
325 : GDALProgressFunc pfnProgress, void * pProgressData )
326 : {
327 42 : int nBands = poSrcDS->GetRasterCount();
328 42 : if (nBands == 0)
329 : {
330 2 : CPLError(CE_Failure, CPLE_NotSupported, "nBands == 0");
331 2 : return NULL;
332 : }
333 :
334 40 : const char* pszDriverName = CSLFetchNameValueDef(papszOptions, "DRIVER", "GTiff");
335 40 : if (EQUAL(pszDriverName, "MEM") || EQUAL(pszDriverName, "VRT"))
336 : {
337 : CPLError(CE_Failure, CPLE_AppDefined, "GDAL %s driver cannot be used as underlying driver",
338 0 : pszDriverName);
339 0 : return NULL;
340 : }
341 :
342 40 : GDALDriverH hTileDriver = GDALGetDriverByName(pszDriverName);
343 40 : if ( hTileDriver == NULL)
344 : {
345 0 : CPLError(CE_Failure, CPLE_AppDefined, "Cannot load GDAL %s driver", pszDriverName);
346 0 : return NULL;
347 : }
348 :
349 40 : GDALDriverH hMemDriver = GDALGetDriverByName("MEM");
350 40 : if (hMemDriver == NULL)
351 : {
352 0 : CPLError(CE_Failure, CPLE_AppDefined, "Cannot load GDAL MEM driver");
353 0 : return NULL;
354 : }
355 :
356 40 : int nXSize = GDALGetRasterXSize(poSrcDS);
357 40 : int nYSize = GDALGetRasterYSize(poSrcDS);
358 :
359 : double adfGeoTransform[6];
360 40 : if (poSrcDS->GetGeoTransform(adfGeoTransform) != CE_None)
361 : {
362 0 : adfGeoTransform[0] = 0;
363 0 : adfGeoTransform[1] = 1;
364 0 : adfGeoTransform[2] = 0;
365 0 : adfGeoTransform[3] = 0;
366 0 : adfGeoTransform[4] = 0;
367 0 : adfGeoTransform[5] = -1;
368 : }
369 40 : else if (adfGeoTransform[2] != 0.0 || adfGeoTransform[4] != 0.0)
370 : {
371 0 : CPLError(CE_Failure, CPLE_AppDefined, "Cannot use geotransform with rotational terms");
372 0 : return NULL;
373 : }
374 :
375 40 : int bTiled = CSLTestBoolean(CSLFetchNameValueDef(papszOptions, "TILED", "YES"));
376 : int nBlockXSize, nBlockYSize;
377 40 : if (bTiled)
378 : {
379 40 : nBlockXSize = atoi(CSLFetchNameValueDef(papszOptions, "BLOCKXSIZE", "256"));
380 40 : nBlockYSize = atoi(CSLFetchNameValueDef(papszOptions, "BLOCKYSIZE", "256"));
381 40 : if (nBlockXSize < 64) nBlockXSize = 64;
382 40 : else if (nBlockXSize > 4096) nBlockXSize = 4096;
383 40 : if (nBlockYSize < 64) nBlockYSize = 64;
384 40 : else if (nBlockYSize > 4096) nBlockYSize = 4096;
385 : }
386 : else
387 : {
388 0 : nBlockXSize = nXSize;
389 0 : nBlockYSize = nYSize;
390 : }
391 :
392 : /* -------------------------------------------------------------------- */
393 : /* Analyze arguments */
394 : /* -------------------------------------------------------------------- */
395 :
396 40 : CPLString osDBName;
397 40 : CPLString osTableName;
398 : VSIStatBuf sBuf;
399 : int bExists;
400 :
401 : /* Skip optionnal RASTERLITE: prefix */
402 40 : const char* pszFilenameWithoutPrefix = pszFilename;
403 40 : if (EQUALN(pszFilename, "RASTERLITE:", 11))
404 2 : pszFilenameWithoutPrefix += 11;
405 :
406 : char** papszTokens = CSLTokenizeStringComplex(
407 40 : pszFilenameWithoutPrefix, ", ", FALSE, FALSE );
408 40 : int nTokens = CSLCount(papszTokens);
409 40 : if (nTokens == 0)
410 : {
411 0 : osDBName = pszFilenameWithoutPrefix;
412 0 : osTableName = CPLGetBasename(pszFilenameWithoutPrefix);
413 : }
414 : else
415 : {
416 40 : osDBName = papszTokens[0];
417 :
418 : int i;
419 42 : for(i=1;i<nTokens;i++)
420 : {
421 2 : if (EQUALN(papszTokens[i], "table=", 6))
422 2 : osTableName = papszTokens[i] + 6;
423 : else
424 : {
425 : CPLError(CE_Warning, CPLE_AppDefined,
426 0 : "Invalid option : %s", papszTokens[i]);
427 : }
428 : }
429 : }
430 :
431 40 : CSLDestroy(papszTokens);
432 40 : papszTokens = NULL;
433 :
434 40 : bExists = (VSIStat(osDBName.c_str(), &sBuf) == 0);
435 :
436 40 : if (osTableName.size() == 0)
437 : {
438 38 : if (bExists)
439 : {
440 : CPLError(CE_Failure, CPLE_AppDefined,
441 0 : "Database already exists. Explicit table name must be specified");
442 0 : return NULL;
443 : }
444 38 : osTableName = CPLGetBasename(osDBName.c_str());
445 : }
446 :
447 40 : CPLString osRasterLayer;
448 40 : osRasterLayer.Printf("%s_rasters", osTableName.c_str());
449 :
450 40 : CPLString osMetatadataLayer;
451 40 : osMetatadataLayer.Printf("%s_metadata", osTableName.c_str());
452 :
453 : /* -------------------------------------------------------------------- */
454 : /* Create or open the SQLite DB */
455 : /* -------------------------------------------------------------------- */
456 :
457 40 : if (OGRGetDriverCount() == 0)
458 2 : OGRRegisterAll();
459 :
460 40 : OGRSFDriverH hSQLiteDriver = OGRGetDriverByName("SQLite");
461 40 : if (hSQLiteDriver == NULL)
462 : {
463 0 : CPLError(CE_Failure, CPLE_AppDefined, "Cannot load OGR SQLite driver");
464 0 : return NULL;
465 : }
466 :
467 : OGRDataSourceH hDS;
468 :
469 : CPLString osOldVal =
470 40 : CPLGetConfigOption("SQLITE_LIST_ALL_TABLES", "FALSE");
471 40 : CPLSetThreadLocalConfigOption("SQLITE_LIST_ALL_TABLES", "TRUE");
472 40 : if (!bExists)
473 : {
474 40 : char** papszOGROptions = CSLAddString(NULL, "SPATIALITE=YES");
475 : hDS = OGR_Dr_CreateDataSource(hSQLiteDriver,
476 40 : osDBName.c_str(), papszOGROptions);
477 40 : CSLDestroy(papszOGROptions);
478 : }
479 : else
480 : {
481 0 : hDS = OGROpen(osDBName.c_str(), TRUE, NULL);
482 : }
483 40 : CPLSetThreadLocalConfigOption("SQLITE_LIST_ALL_TABLES", osOldVal.c_str());
484 :
485 40 : if (hDS == NULL)
486 : {
487 : CPLError(CE_Failure, CPLE_AppDefined,
488 4 : "Cannot load or create SQLite database");
489 4 : return NULL;
490 : }
491 :
492 36 : CPLString osSQL;
493 :
494 : /* -------------------------------------------------------------------- */
495 : /* Get the SRID for the SRS */
496 : /* -------------------------------------------------------------------- */
497 36 : int nSRSId = RasterliteInsertSRID(hDS, poSrcDS->GetProjectionRef());
498 :
499 : /* -------------------------------------------------------------------- */
500 : /* Create or wipe existing tables */
501 : /* -------------------------------------------------------------------- */
502 : int bWipeExistingData =
503 36 : CSLTestBoolean(CSLFetchNameValueDef(papszOptions, "WIPE", "NO"));
504 :
505 : hDS = RasterliteCreateTables(hDS, osTableName.c_str(),
506 36 : nSRSId, bWipeExistingData);
507 36 : if (hDS == NULL)
508 0 : return NULL;
509 :
510 36 : OGRLayerH hRasterLayer = OGR_DS_GetLayerByName(hDS, osRasterLayer.c_str());
511 36 : OGRLayerH hMetadataLayer = OGR_DS_GetLayerByName(hDS, osMetatadataLayer.c_str());
512 36 : if (hRasterLayer == NULL || hMetadataLayer == NULL)
513 : {
514 : CPLError(CE_Failure, CPLE_AppDefined,
515 0 : "Cannot find metadata and/or raster tables");
516 0 : OGRReleaseDataSource(hDS);
517 0 : return NULL;
518 : }
519 :
520 : /* -------------------------------------------------------------------- */
521 : /* Check if there is overlapping data and warn the user */
522 : /* -------------------------------------------------------------------- */
523 36 : double minx = adfGeoTransform[0];
524 36 : double maxx = adfGeoTransform[0] + nXSize * adfGeoTransform[1];
525 36 : double maxy = adfGeoTransform[3];
526 36 : double miny = adfGeoTransform[3] + nYSize * adfGeoTransform[5];
527 :
528 : osSQL.Printf("SELECT COUNT(geometry) FROM \"%s\" "
529 : "WHERE rowid IN "
530 : "(SELECT pkid FROM \"idx_%s_metadata_geometry\" "
531 : "WHERE xmin < %.15f AND xmax > %.15f "
532 : "AND ymin < %.15f AND ymax > %.15f) "
533 : "AND pixel_x_size >= %.15f AND pixel_x_size <= %.15f AND "
534 : "pixel_y_size >= %.15f AND pixel_y_size <= %.15f",
535 : osMetatadataLayer.c_str(),
536 : osTableName.c_str(),
537 : maxx, minx, maxy, miny,
538 72 : adfGeoTransform[1] - 1e-15, adfGeoTransform[1] + 1e-15,
539 108 : - adfGeoTransform[5] - 1e-15, - adfGeoTransform[5] + 1e-15);
540 :
541 36 : int nOverlappingGeoms = 0;
542 36 : OGRLayerH hCountLyr = OGR_DS_ExecuteSQL(hDS, osSQL.c_str(), NULL, NULL);
543 36 : if (hCountLyr)
544 : {
545 36 : OGRFeatureH hFeat = OGR_L_GetNextFeature(hCountLyr);
546 36 : if (hFeat)
547 : {
548 36 : nOverlappingGeoms = OGR_F_GetFieldAsInteger(hFeat, 0);
549 36 : OGR_F_Destroy(hFeat);
550 : }
551 36 : OGR_DS_ReleaseResultSet(hDS, hCountLyr);
552 : }
553 :
554 36 : if (nOverlappingGeoms != 0)
555 : {
556 : CPLError(CE_Warning, CPLE_AppDefined,
557 : "Raster tiles already exist in the %s table within "
558 : "the extent of the data to be inserted in",
559 0 : osTableName.c_str());
560 : }
561 :
562 : /* -------------------------------------------------------------------- */
563 : /* Iterate over blocks to add data into raster and metadata tables */
564 : /* -------------------------------------------------------------------- */
565 36 : int nXBlocks = (nXSize + nBlockXSize - 1) / nBlockXSize;
566 36 : int nYBlocks = (nYSize + nBlockYSize - 1) / nBlockYSize;
567 :
568 36 : GDALDataType eDataType = poSrcDS->GetRasterBand(1)->GetRasterDataType();
569 36 : int nDataTypeSize = GDALGetDataTypeSize(eDataType) / 8;
570 : GByte* pabyMEMDSBuffer =
571 36 : (GByte*)VSIMalloc3(nBlockXSize, nBlockYSize, nBands * nDataTypeSize);
572 36 : if (pabyMEMDSBuffer == NULL)
573 : {
574 0 : OGRReleaseDataSource(hDS);
575 0 : return NULL;
576 : }
577 :
578 36 : CPLString osTempFileName;
579 36 : osTempFileName.Printf("/vsimem/%p", hDS);
580 :
581 36 : int nTileId = 0;
582 36 : int nBlocks = 0;
583 36 : int nTotalBlocks = nXBlocks * nYBlocks;
584 :
585 36 : char** papszTileDriverOptions = RasterliteGetTileDriverOptions(papszOptions);
586 :
587 36 : OGR_DS_ExecuteSQL(hDS, "BEGIN", NULL, NULL);
588 :
589 36 : CPLErr eErr = CE_None;
590 : int nBlockXOff, nBlockYOff;
591 72 : for(nBlockYOff=0;eErr == CE_None && nBlockYOff<nYBlocks;nBlockYOff++)
592 : {
593 72 : for(nBlockXOff=0;eErr == CE_None && nBlockXOff<nXBlocks;nBlockXOff++)
594 : {
595 : /* -------------------------------------------------------------------- */
596 : /* Create in-memory tile */
597 : /* -------------------------------------------------------------------- */
598 36 : int nReqXSize = nBlockXSize, nReqYSize = nBlockYSize;
599 36 : if ((nBlockXOff+1) * nBlockXSize > nXSize)
600 36 : nReqXSize = nXSize - nBlockXOff * nBlockXSize;
601 36 : if ((nBlockYOff+1) * nBlockYSize > nYSize)
602 36 : nReqYSize = nYSize - nBlockYOff * nBlockYSize;
603 :
604 : eErr = poSrcDS->RasterIO(GF_Read,
605 : nBlockXOff * nBlockXSize,
606 : nBlockYOff * nBlockYSize,
607 : nReqXSize, nReqYSize,
608 : pabyMEMDSBuffer, nReqXSize, nReqYSize,
609 : eDataType, nBands, NULL,
610 36 : 0, 0, 0);
611 36 : if (eErr != CE_None)
612 : {
613 0 : break;
614 : }
615 :
616 : GDALDatasetH hMemDS = GDALCreate(hMemDriver, "MEM:::",
617 : nReqXSize, nReqYSize, 0,
618 36 : eDataType, NULL);
619 36 : if (hMemDS == NULL)
620 : {
621 0 : eErr = CE_Failure;
622 0 : break;
623 : }
624 :
625 : int iBand;
626 96 : for(iBand = 0; iBand < nBands; iBand ++)
627 : {
628 60 : char** papszMEMDSOptions = NULL;
629 : char szTmp[64];
630 60 : memset(szTmp, 0, sizeof(szTmp));
631 : CPLPrintPointer(szTmp,
632 : pabyMEMDSBuffer + iBand * nDataTypeSize *
633 60 : nReqXSize * nReqYSize, sizeof(szTmp));
634 60 : papszMEMDSOptions = CSLSetNameValue(papszMEMDSOptions, "DATAPOINTER", szTmp);
635 60 : GDALAddBand(hMemDS, eDataType, papszMEMDSOptions);
636 60 : CSLDestroy(papszMEMDSOptions);
637 : }
638 :
639 : GDALDatasetH hOutDS = GDALCreateCopy(hTileDriver,
640 : osTempFileName.c_str(), hMemDS, FALSE,
641 36 : papszTileDriverOptions, NULL, NULL);
642 :
643 36 : GDALClose(hMemDS);
644 36 : if (hOutDS)
645 36 : GDALClose(hOutDS);
646 : else
647 : {
648 0 : eErr = CE_Failure;
649 0 : break;
650 : }
651 :
652 : /* -------------------------------------------------------------------- */
653 : /* Insert new entry into raster table */
654 : /* -------------------------------------------------------------------- */
655 :
656 36 : vsi_l_offset nDataLength = 0;
657 : GByte *pabyData = VSIGetMemFileBuffer( osTempFileName.c_str(),
658 72 : &nDataLength, FALSE);
659 :
660 36 : OGRFeatureH hFeat = OGR_F_Create( OGR_L_GetLayerDefn(hRasterLayer) );
661 36 : OGR_F_SetFieldBinary(hFeat, 0, (int)nDataLength, pabyData);
662 :
663 36 : OGR_L_CreateFeature(hRasterLayer, hFeat);
664 : /* Query raster ID to set it as the ID of the associated metadata */
665 36 : int nRasterID = (int)OGR_F_GetFID(hFeat);
666 :
667 36 : OGR_F_Destroy(hFeat);
668 :
669 36 : VSIUnlink(osTempFileName.c_str());
670 :
671 : /* -------------------------------------------------------------------- */
672 : /* Insert new entry into metadata table */
673 : /* -------------------------------------------------------------------- */
674 :
675 36 : hFeat = OGR_F_Create( OGR_L_GetLayerDefn(hMetadataLayer) );
676 36 : OGR_F_SetFID(hFeat, nRasterID);
677 36 : OGR_F_SetFieldString(hFeat, 0, GDALGetDescription(poSrcDS));
678 36 : OGR_F_SetFieldInteger(hFeat, 1, nTileId ++);
679 36 : OGR_F_SetFieldInteger(hFeat, 2, nReqXSize);
680 36 : OGR_F_SetFieldInteger(hFeat, 3, nReqYSize);
681 36 : OGR_F_SetFieldDouble(hFeat, 4, adfGeoTransform[1]);
682 36 : OGR_F_SetFieldDouble(hFeat, 5, -adfGeoTransform[5]);
683 :
684 36 : minx = adfGeoTransform[0] +
685 36 : (nBlockXSize * nBlockXOff) * adfGeoTransform[1];
686 36 : maxx = adfGeoTransform[0] +
687 36 : (nBlockXSize * nBlockXOff + nReqXSize) * adfGeoTransform[1];
688 36 : maxy = adfGeoTransform[3] +
689 36 : (nBlockYSize * nBlockYOff) * adfGeoTransform[5];
690 36 : miny = adfGeoTransform[3] +
691 36 : (nBlockYSize * nBlockYOff + nReqYSize) * adfGeoTransform[5];
692 :
693 36 : OGRGeometryH hRectangle = OGR_G_CreateGeometry(wkbPolygon);
694 36 : OGRGeometryH hLinearRing = OGR_G_CreateGeometry(wkbLinearRing);
695 36 : OGR_G_AddPoint_2D(hLinearRing, minx, miny);
696 36 : OGR_G_AddPoint_2D(hLinearRing, minx, maxy);
697 36 : OGR_G_AddPoint_2D(hLinearRing, maxx, maxy);
698 36 : OGR_G_AddPoint_2D(hLinearRing, maxx, miny);
699 36 : OGR_G_AddPoint_2D(hLinearRing, minx, miny);
700 36 : OGR_G_AddGeometryDirectly(hRectangle, hLinearRing);
701 :
702 36 : OGR_F_SetGeometryDirectly(hFeat, hRectangle);
703 :
704 36 : OGR_L_CreateFeature(hMetadataLayer, hFeat);
705 36 : OGR_F_Destroy(hFeat);
706 :
707 36 : nBlocks++;
708 36 : if (pfnProgress && !pfnProgress(1.0 * nBlocks / nTotalBlocks,
709 : NULL, pProgressData))
710 0 : eErr = CE_Failure;
711 : }
712 : }
713 :
714 36 : if (eErr == CE_None)
715 36 : OGR_DS_ExecuteSQL(hDS, "COMMIT", NULL, NULL);
716 : else
717 0 : OGR_DS_ExecuteSQL(hDS, "ROLLBACK", NULL, NULL);
718 :
719 36 : CSLDestroy(papszTileDriverOptions);
720 :
721 36 : VSIFree(pabyMEMDSBuffer);
722 :
723 36 : OGRReleaseDataSource(hDS);
724 :
725 36 : return (GDALDataset*) GDALOpen(pszFilename, GA_Update);
726 : }
727 :
728 : /************************************************************************/
729 : /* RasterliteDelete () */
730 : /************************************************************************/
731 :
732 6 : CPLErr RasterliteDelete(const char* pszFilename)
733 : {
734 6 : return CE_None;
735 : }
|