1 : /******************************************************************************
2 : * $Id: rasterliteoverviews.cpp 24801 2012-08-18 14:29:41Z 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: rasterliteoverviews.cpp 24801 2012-08-18 14:29:41Z rouault $");
37 :
38 : /************************************************************************/
39 : /* ReloadOverviews() */
40 : /************************************************************************/
41 :
42 4 : CPLErr RasterliteDataset::ReloadOverviews()
43 : {
44 4 : if (nLevel != 0)
45 0 : return CE_Failure;
46 :
47 : /* -------------------------------------------------------------------- */
48 : /* Fetch resolutions */
49 : /* -------------------------------------------------------------------- */
50 :
51 4 : CPLString osSQL;
52 4 : OGRLayerH hRasterPyramidsLyr = OGR_DS_GetLayerByName(hDS, "raster_pyramids");
53 4 : if (hRasterPyramidsLyr)
54 : {
55 : osSQL.Printf("SELECT pixel_x_size, pixel_y_size "
56 : "FROM raster_pyramids WHERE table_prefix = '%s' "
57 : "ORDER BY pixel_x_size ASC",
58 4 : osTableName.c_str());
59 : }
60 : else
61 : {
62 : osSQL.Printf("SELECT DISTINCT(pixel_x_size), pixel_y_size "
63 : "FROM \"%s_metadata\" WHERE pixel_x_size != 0 "
64 : "ORDER BY pixel_x_size ASC",
65 0 : osTableName.c_str());
66 : }
67 :
68 4 : OGRLayerH hSQLLyr = OGR_DS_ExecuteSQL(hDS, osSQL.c_str(), NULL, NULL);
69 4 : if (hSQLLyr == NULL)
70 : {
71 0 : if (hRasterPyramidsLyr == NULL)
72 0 : return CE_Failure;
73 :
74 : osSQL.Printf("SELECT DISTINCT(pixel_x_size), pixel_y_size "
75 : "FROM \"%s_metadata\" WHERE pixel_x_size != 0 "
76 : "ORDER BY pixel_x_size ASC",
77 0 : osTableName.c_str());
78 :
79 0 : hSQLLyr = OGR_DS_ExecuteSQL(hDS, osSQL.c_str(), NULL, NULL);
80 0 : if (hSQLLyr == NULL)
81 0 : return CE_Failure;
82 : }
83 :
84 : /* -------------------------------------------------------------------- */
85 : /* Cleanup */
86 : /* -------------------------------------------------------------------- */
87 :
88 : int i;
89 6 : for(i=1;i<nResolutions;i++)
90 2 : delete papoOverviews[i-1];
91 4 : CPLFree(papoOverviews);
92 4 : papoOverviews = NULL;
93 4 : CPLFree(padfXResolutions);
94 4 : padfXResolutions = NULL;
95 4 : CPLFree(padfYResolutions);
96 4 : padfYResolutions = NULL;
97 :
98 : /* -------------------------------------------------------------------- */
99 : /* Rebuild arrays */
100 : /* -------------------------------------------------------------------- */
101 :
102 4 : nResolutions = OGR_L_GetFeatureCount(hSQLLyr, TRUE);
103 :
104 : padfXResolutions =
105 4 : (double*)CPLMalloc(sizeof(double) * nResolutions);
106 : padfYResolutions =
107 4 : (double*)CPLMalloc(sizeof(double) * nResolutions);
108 :
109 4 : i = 0;
110 : OGRFeatureH hFeat;
111 18 : while((hFeat = OGR_L_GetNextFeature(hSQLLyr)) != NULL)
112 : {
113 10 : padfXResolutions[i] = OGR_F_GetFieldAsDouble(hFeat, 0);
114 10 : padfYResolutions[i] = OGR_F_GetFieldAsDouble(hFeat, 1);
115 :
116 10 : OGR_F_Destroy(hFeat);
117 :
118 10 : i ++;
119 : }
120 :
121 4 : OGR_DS_ReleaseResultSet(hDS, hSQLLyr);
122 4 : hSQLLyr = NULL;
123 :
124 : /* -------------------------------------------------------------------- */
125 : /* Add overview levels as internal datasets */
126 : /* -------------------------------------------------------------------- */
127 4 : if (nResolutions > 1)
128 : {
129 4 : CPLString osRasterTableName = osTableName;
130 4 : osRasterTableName += "_rasters";
131 :
132 4 : OGRLayerH hRasterLyr = OGR_DS_GetLayerByName(hDS, osRasterTableName.c_str());
133 :
134 : papoOverviews = (RasterliteDataset**)
135 4 : CPLCalloc(nResolutions - 1, sizeof(RasterliteDataset*));
136 : int nLev;
137 10 : for(nLev=1;nLev<nResolutions;nLev++)
138 : {
139 : int nOvrBands;
140 : GDALDataType eOvrDataType;
141 : int nBlockXSize, nBlockYSize;
142 6 : if (GetBlockParams(hRasterLyr, nLev, &nOvrBands, &eOvrDataType,
143 : &nBlockXSize, &nBlockYSize))
144 : {
145 6 : if (eOvrDataType == GDT_Byte && nOvrBands == 1 && nBands == 3)
146 0 : nOvrBands = 3;
147 :
148 6 : papoOverviews[nLev-1] = new RasterliteDataset(this, nLev);
149 :
150 : int iBand;
151 12 : for(iBand=0;iBand<nBands;iBand++)
152 : {
153 6 : papoOverviews[nLev-1]->SetBand(iBand+1,
154 : new RasterliteBand(papoOverviews[nLev-1], iBand+1, eOvrDataType,
155 12 : nBlockXSize, nBlockYSize));
156 : }
157 : }
158 : else
159 : {
160 : CPLError(CE_Failure, CPLE_AppDefined,
161 0 : "Cannot find block characteristics for overview %d", nLev);
162 0 : papoOverviews[nLev-1] = NULL;
163 : }
164 4 : }
165 : }
166 :
167 4 : return CE_None;
168 : }
169 :
170 : /************************************************************************/
171 : /* CleanOverviews() */
172 : /************************************************************************/
173 :
174 2 : CPLErr RasterliteDataset::CleanOverviews()
175 : {
176 2 : CPLString osSQL;
177 :
178 2 : if (nLevel != 0)
179 0 : return CE_Failure;
180 :
181 2 : osSQL.Printf("BEGIN");
182 2 : OGR_DS_ExecuteSQL(hDS, osSQL.c_str(), NULL, NULL);
183 :
184 2 : CPLString osResolutionCond;
185 : osResolutionCond.Printf(
186 : "(pixel_x_size < %.15f OR pixel_x_size > %.15f) AND "
187 : "(pixel_y_size < %.15f OR pixel_y_size > %.15f)",
188 4 : padfXResolutions[0] - 1e-15, padfXResolutions[0] + 1e-15,
189 6 : padfYResolutions[0] - 1e-15, padfYResolutions[0] + 1e-15);
190 :
191 : osSQL.Printf("DELETE FROM \"%s_rasters\" WHERE id "
192 : "IN(SELECT id FROM \"%s_metadata\" WHERE %s)",
193 : osTableName.c_str(), osTableName.c_str(),
194 2 : osResolutionCond.c_str());
195 2 : OGR_DS_ExecuteSQL(hDS, osSQL.c_str(), NULL, NULL);
196 :
197 : osSQL.Printf("DELETE FROM \"%s_metadata\" WHERE %s",
198 2 : osTableName.c_str(), osResolutionCond.c_str());
199 2 : OGR_DS_ExecuteSQL(hDS, osSQL.c_str(), NULL, NULL);
200 :
201 2 : OGRLayerH hRasterPyramidsLyr = OGR_DS_GetLayerByName(hDS, "raster_pyramids");
202 2 : if (hRasterPyramidsLyr)
203 : {
204 : osSQL.Printf("DELETE FROM raster_pyramids WHERE table_prefix = '%s' AND %s",
205 2 : osTableName.c_str(), osResolutionCond.c_str());
206 2 : OGR_DS_ExecuteSQL(hDS, osSQL.c_str(), NULL, NULL);
207 : }
208 :
209 2 : osSQL.Printf("COMMIT");
210 2 : OGR_DS_ExecuteSQL(hDS, osSQL.c_str(), NULL, NULL);
211 :
212 : int i;
213 4 : for(i=1;i<nResolutions;i++)
214 2 : delete papoOverviews[i-1];
215 2 : CPLFree(papoOverviews);
216 2 : papoOverviews = NULL;
217 2 : nResolutions = 1;
218 :
219 2 : return CE_None;
220 : }
221 :
222 : /************************************************************************/
223 : /* CleanOverviewLevel() */
224 : /************************************************************************/
225 :
226 4 : CPLErr RasterliteDataset::CleanOverviewLevel(int nOvrFactor)
227 : {
228 4 : CPLString osSQL;
229 :
230 4 : if (nLevel != 0)
231 0 : return CE_Failure;
232 :
233 : /* -------------------------------------------------------------------- */
234 : /* Find the index of the overview matching the overview factor */
235 : /* -------------------------------------------------------------------- */
236 : int iLev;
237 6 : for(iLev=1;iLev<nResolutions;iLev++)
238 : {
239 2 : if (fabs(padfXResolutions[0] * nOvrFactor - padfXResolutions[iLev]) < 1e-15 &&
240 0 : fabs(padfYResolutions[0] * nOvrFactor - padfYResolutions[iLev]) < 1e-15)
241 0 : break;
242 : }
243 :
244 4 : if (iLev == nResolutions)
245 4 : return CE_None;
246 :
247 : /* -------------------------------------------------------------------- */
248 : /* Now clean existing overviews at that resolution */
249 : /* -------------------------------------------------------------------- */
250 :
251 0 : osSQL.Printf("BEGIN");
252 0 : OGR_DS_ExecuteSQL(hDS, osSQL.c_str(), NULL, NULL);
253 :
254 0 : CPLString osResolutionCond;
255 : osResolutionCond.Printf(
256 : "pixel_x_size >= %.15f AND pixel_x_size <= %.15f AND "
257 : "pixel_y_size >= %.15f AND pixel_y_size <= %.15f",
258 0 : padfXResolutions[iLev] - 1e-15, padfXResolutions[iLev] + 1e-15,
259 0 : padfYResolutions[iLev] - 1e-15, padfYResolutions[iLev] + 1e-15);
260 :
261 : osSQL.Printf("DELETE FROM \"%s_rasters\" WHERE id "
262 : "IN(SELECT id FROM \"%s_metadata\" WHERE %s)",
263 : osTableName.c_str(), osTableName.c_str(),
264 0 : osResolutionCond.c_str());
265 0 : OGR_DS_ExecuteSQL(hDS, osSQL.c_str(), NULL, NULL);
266 :
267 : osSQL.Printf("DELETE FROM \"%s_metadata\" WHERE %s",
268 0 : osTableName.c_str(), osResolutionCond.c_str());
269 0 : OGR_DS_ExecuteSQL(hDS, osSQL.c_str(), NULL, NULL);
270 :
271 0 : OGRLayerH hRasterPyramidsLyr = OGR_DS_GetLayerByName(hDS, "raster_pyramids");
272 0 : if (hRasterPyramidsLyr)
273 : {
274 : osSQL.Printf("DELETE FROM raster_pyramids WHERE table_prefix = '%s' AND %s",
275 0 : osTableName.c_str(), osResolutionCond.c_str());
276 0 : OGR_DS_ExecuteSQL(hDS, osSQL.c_str(), NULL, NULL);
277 : }
278 :
279 0 : osSQL.Printf("COMMIT");
280 0 : OGR_DS_ExecuteSQL(hDS, osSQL.c_str(), NULL, NULL);
281 :
282 0 : return CE_None;
283 : }
284 :
285 : /************************************************************************/
286 : /* CleanOverviewLevel() */
287 : /************************************************************************/
288 :
289 4 : CPLErr RasterliteDataset::CreateOverviewLevel(const char * pszResampling,
290 : int nOvrFactor,
291 : char** papszOptions,
292 : GDALProgressFunc pfnProgress,
293 : void * pProgressData)
294 : {
295 :
296 4 : double dfXResolution = padfXResolutions[0] * nOvrFactor;
297 4 : double dfYResolution = padfXResolutions[0] * nOvrFactor;
298 :
299 4 : CPLString osSQL;
300 :
301 4 : int nOvrXSize = nRasterXSize / nOvrFactor;
302 4 : int nOvrYSize = nRasterYSize / nOvrFactor;
303 :
304 4 : if (nOvrXSize == 0 || nOvrYSize == 0)
305 0 : return CE_Failure;
306 :
307 4 : int bTiled = CSLTestBoolean(CSLFetchNameValueDef(papszOptions, "TILED", "YES"));
308 : int nBlockXSize, nBlockYSize;
309 4 : if (bTiled)
310 : {
311 4 : nBlockXSize = atoi(CSLFetchNameValueDef(papszOptions, "BLOCKXSIZE", "256"));
312 4 : nBlockYSize = atoi(CSLFetchNameValueDef(papszOptions, "BLOCKYSIZE", "256"));
313 4 : if (nBlockXSize < 64) nBlockXSize = 64;
314 4 : else if (nBlockXSize > 4096) nBlockXSize = 4096;
315 4 : if (nBlockYSize < 64) nBlockYSize = 64;
316 4 : else if (nBlockYSize > 4096) nBlockYSize = 4096;
317 : }
318 : else
319 : {
320 0 : nBlockXSize = nOvrXSize;
321 0 : nBlockYSize = nOvrYSize;
322 : }
323 :
324 4 : int nXBlocks = (nOvrXSize + nBlockXSize - 1) / nBlockXSize;
325 4 : int nYBlocks = (nOvrYSize + nBlockYSize - 1) / nBlockYSize;
326 :
327 4 : const char* pszDriverName = CSLFetchNameValueDef(papszOptions, "DRIVER", "GTiff");
328 4 : if (EQUAL(pszDriverName, "MEM") || EQUAL(pszDriverName, "VRT"))
329 : {
330 : CPLError(CE_Failure, CPLE_AppDefined, "GDAL %s driver cannot be used as underlying driver",
331 0 : pszDriverName);
332 0 : return CE_Failure;
333 : }
334 4 : GDALDriverH hTileDriver = GDALGetDriverByName(pszDriverName);
335 4 : if (hTileDriver == NULL)
336 : {
337 0 : CPLError(CE_Failure, CPLE_AppDefined, "Cannot load GDAL %s driver", pszDriverName);
338 0 : return CE_Failure;
339 : }
340 :
341 4 : GDALDriverH hMemDriver = GDALGetDriverByName("MEM");
342 4 : if (hMemDriver == NULL)
343 : {
344 0 : CPLError(CE_Failure, CPLE_AppDefined, "Cannot load GDAL MEM driver");
345 0 : return CE_Failure;
346 : }
347 :
348 4 : GDALDataType eDataType = GetRasterBand(1)->GetRasterDataType();
349 4 : int nDataTypeSize = GDALGetDataTypeSize(eDataType) / 8;
350 : GByte* pabyMEMDSBuffer =
351 4 : (GByte*)VSIMalloc3(nBlockXSize, nBlockYSize, nBands * nDataTypeSize);
352 4 : if (pabyMEMDSBuffer == NULL)
353 : {
354 0 : return CE_Failure;
355 : }
356 :
357 4 : CPLString osTempFileName;
358 4 : osTempFileName.Printf("/vsimem/%p", hDS);
359 :
360 4 : int nTileId = 0;
361 4 : int nBlocks = 0;
362 4 : int nTotalBlocks = nXBlocks * nYBlocks;
363 :
364 4 : CPLString osRasterLayer;
365 4 : osRasterLayer.Printf("%s_rasters", osTableName.c_str());
366 :
367 4 : CPLString osMetatadataLayer;
368 4 : osMetatadataLayer.Printf("%s_metadata", osTableName.c_str());
369 :
370 4 : OGRLayerH hRasterLayer = OGR_DS_GetLayerByName(hDS, osRasterLayer.c_str());
371 4 : OGRLayerH hMetadataLayer = OGR_DS_GetLayerByName(hDS, osMetatadataLayer.c_str());
372 :
373 4 : CPLString osSourceName = "unknown";
374 :
375 : osSQL.Printf("SELECT source_name FROM \"%s\" WHERE "
376 : "pixel_x_size >= %.15f AND pixel_x_size <= %.15f AND "
377 : "pixel_y_size >= %.15f AND pixel_y_size <= %.15f LIMIT 1",
378 : osMetatadataLayer.c_str(),
379 8 : padfXResolutions[0] - 1e-15, padfXResolutions[0] + 1e-15,
380 12 : padfYResolutions[0] - 1e-15, padfYResolutions[0] + 1e-15);
381 4 : OGRLayerH hSQLLyr = OGR_DS_ExecuteSQL(hDS, osSQL.c_str(), NULL, NULL);
382 4 : if (hSQLLyr)
383 : {
384 4 : OGRFeatureH hFeat = OGR_L_GetNextFeature(hSQLLyr);
385 4 : if (hFeat)
386 : {
387 4 : const char* pszVal = OGR_F_GetFieldAsString(hFeat, 0);
388 4 : if (pszVal)
389 4 : osSourceName = pszVal;
390 4 : OGR_F_Destroy(hFeat);
391 : }
392 4 : OGR_DS_ReleaseResultSet(hDS, hSQLLyr);
393 : }
394 :
395 : /* -------------------------------------------------------------------- */
396 : /* Compute up to which existing overview level we can use for */
397 : /* computing the requested overview */
398 : /* -------------------------------------------------------------------- */
399 : int iLev;
400 4 : nLimitOvrCount = 0;
401 6 : for(iLev=1;iLev<nResolutions;iLev++)
402 : {
403 4 : if (!(padfXResolutions[iLev] < dfXResolution - 1e-10 &&
404 2 : padfYResolutions[iLev] < dfYResolution - 1e-10))
405 : {
406 0 : break;
407 : }
408 2 : nLimitOvrCount++;
409 : }
410 : /* -------------------------------------------------------------------- */
411 : /* Allocate buffer for tile of previous overview level */
412 : /* -------------------------------------------------------------------- */
413 :
414 : GDALDataset* poPrevOvrLevel =
415 2 : (papoOverviews != NULL && iLev >= 2 && iLev <= nResolutions && papoOverviews[iLev-2]) ?
416 6 : papoOverviews[iLev-2] : this;
417 4 : double dfRatioPrevOvr = poPrevOvrLevel->GetRasterBand(1)->GetXSize() / nOvrXSize;
418 4 : int nPrevOvrBlockXSize = (int)(nBlockXSize * dfRatioPrevOvr + 0.5);
419 4 : int nPrevOvrBlockYSize = (int)(nBlockYSize * dfRatioPrevOvr + 0.5);
420 4 : GByte* pabyPrevOvrMEMDSBuffer = NULL;
421 :
422 4 : if( !EQUALN(pszResampling, "NEAR", 4))
423 : {
424 : pabyPrevOvrMEMDSBuffer =
425 2 : (GByte*)VSIMalloc3(nPrevOvrBlockXSize, nPrevOvrBlockYSize, nBands * nDataTypeSize);
426 2 : if (pabyPrevOvrMEMDSBuffer == NULL)
427 : {
428 0 : VSIFree(pabyMEMDSBuffer);
429 0 : return CE_Failure;
430 : }
431 : }
432 :
433 : /* -------------------------------------------------------------------- */
434 : /* Iterate over blocks to add data into raster and metadata tables */
435 : /* -------------------------------------------------------------------- */
436 :
437 4 : char** papszTileDriverOptions = RasterliteGetTileDriverOptions(papszOptions);
438 :
439 4 : OGR_DS_ExecuteSQL(hDS, "BEGIN", NULL, NULL);
440 :
441 4 : CPLErr eErr = CE_None;
442 : int nBlockXOff, nBlockYOff;
443 8 : for(nBlockYOff=0;eErr == CE_None && nBlockYOff<nYBlocks;nBlockYOff++)
444 : {
445 8 : for(nBlockXOff=0;eErr == CE_None && nBlockXOff<nXBlocks;nBlockXOff++)
446 : {
447 4 : GDALDatasetH hPrevOvrMemDS = NULL;
448 :
449 : /* -------------------------------------------------------------------- */
450 : /* Create in-memory tile */
451 : /* -------------------------------------------------------------------- */
452 4 : int nReqXSize = nBlockXSize, nReqYSize = nBlockYSize;
453 4 : if ((nBlockXOff+1) * nBlockXSize > nOvrXSize)
454 4 : nReqXSize = nOvrXSize - nBlockXOff * nBlockXSize;
455 4 : if ((nBlockYOff+1) * nBlockYSize > nOvrYSize)
456 4 : nReqYSize = nOvrYSize - nBlockYOff * nBlockYSize;
457 :
458 4 : if( pabyPrevOvrMEMDSBuffer != NULL )
459 : {
460 : int nPrevOvrReqXSize =
461 2 : (int)(nReqXSize * dfRatioPrevOvr + 0.5);
462 : int nPrevOvrReqYSize =
463 2 : (int)(nReqYSize * dfRatioPrevOvr + 0.5);
464 :
465 : eErr = RasterIO(GF_Read,
466 : nBlockXOff * nBlockXSize * nOvrFactor,
467 : nBlockYOff * nBlockYSize * nOvrFactor,
468 : nReqXSize * nOvrFactor, nReqYSize * nOvrFactor,
469 : pabyPrevOvrMEMDSBuffer, nPrevOvrReqXSize, nPrevOvrReqYSize,
470 : eDataType, nBands, NULL,
471 2 : 0, 0, 0);
472 :
473 2 : if (eErr != CE_None)
474 : {
475 0 : break;
476 : }
477 :
478 : hPrevOvrMemDS = GDALCreate(hMemDriver, "MEM:::",
479 : nPrevOvrReqXSize, nPrevOvrReqYSize, 0,
480 2 : eDataType, NULL);
481 :
482 2 : if (hPrevOvrMemDS == NULL)
483 : {
484 0 : eErr = CE_Failure;
485 0 : break;
486 : }
487 :
488 : int iBand;
489 4 : for(iBand = 0; iBand < nBands; iBand ++)
490 : {
491 2 : char** papszOptions = NULL;
492 : char szTmp[64];
493 2 : memset(szTmp, 0, sizeof(szTmp));
494 : CPLPrintPointer(szTmp,
495 : pabyPrevOvrMEMDSBuffer + iBand * nDataTypeSize *
496 2 : nPrevOvrReqXSize * nPrevOvrReqYSize, sizeof(szTmp));
497 2 : papszOptions = CSLSetNameValue(papszOptions, "DATAPOINTER", szTmp);
498 2 : GDALAddBand(hPrevOvrMemDS, eDataType, papszOptions);
499 2 : CSLDestroy(papszOptions);
500 : }
501 : }
502 : else
503 : {
504 : eErr = RasterIO(GF_Read,
505 : nBlockXOff * nBlockXSize * nOvrFactor,
506 : nBlockYOff * nBlockYSize * nOvrFactor,
507 : nReqXSize * nOvrFactor, nReqYSize * nOvrFactor,
508 : pabyMEMDSBuffer, nReqXSize, nReqYSize,
509 : eDataType, nBands, NULL,
510 2 : 0, 0, 0);
511 2 : if (eErr != CE_None)
512 : {
513 0 : break;
514 : }
515 : }
516 :
517 : GDALDatasetH hMemDS = GDALCreate(hMemDriver, "MEM:::",
518 : nReqXSize, nReqYSize, 0,
519 4 : eDataType, NULL);
520 4 : if (hMemDS == NULL)
521 : {
522 0 : eErr = CE_Failure;
523 0 : break;
524 : }
525 :
526 : int iBand;
527 8 : for(iBand = 0; iBand < nBands; iBand ++)
528 : {
529 4 : char** papszOptions = NULL;
530 : char szTmp[64];
531 4 : memset(szTmp, 0, sizeof(szTmp));
532 : CPLPrintPointer(szTmp,
533 : pabyMEMDSBuffer + iBand * nDataTypeSize *
534 4 : nReqXSize * nReqYSize, sizeof(szTmp));
535 4 : papszOptions = CSLSetNameValue(papszOptions, "DATAPOINTER", szTmp);
536 4 : GDALAddBand(hMemDS, eDataType, papszOptions);
537 4 : CSLDestroy(papszOptions);
538 : }
539 :
540 4 : if( hPrevOvrMemDS != NULL )
541 : {
542 4 : for(iBand = 0; iBand < nBands; iBand ++)
543 : {
544 2 : GDALRasterBandH hDstOvrBand = GDALGetRasterBand(hMemDS, iBand+1);
545 :
546 : eErr = GDALRegenerateOverviews( GDALGetRasterBand(hPrevOvrMemDS, iBand+1),
547 : 1, &hDstOvrBand,
548 : pszResampling,
549 2 : NULL, NULL );
550 2 : if( eErr != CE_None )
551 0 : break;
552 : }
553 :
554 2 : GDALClose(hPrevOvrMemDS);
555 : }
556 :
557 : GDALDatasetH hOutDS = GDALCreateCopy(hTileDriver,
558 : osTempFileName.c_str(), hMemDS, FALSE,
559 4 : papszTileDriverOptions, NULL, NULL);
560 :
561 4 : GDALClose(hMemDS);
562 4 : if (hOutDS)
563 4 : GDALClose(hOutDS);
564 : else
565 : {
566 0 : eErr = CE_Failure;
567 0 : break;
568 : }
569 :
570 : /* -------------------------------------------------------------------- */
571 : /* Insert new entry into raster table */
572 : /* -------------------------------------------------------------------- */
573 :
574 : vsi_l_offset nDataLength;
575 : GByte *pabyData = VSIGetMemFileBuffer( osTempFileName.c_str(),
576 8 : &nDataLength, FALSE);
577 :
578 4 : OGRFeatureH hFeat = OGR_F_Create( OGR_L_GetLayerDefn(hRasterLayer) );
579 4 : OGR_F_SetFieldBinary(hFeat, 0, (int)nDataLength, pabyData);
580 :
581 4 : OGR_L_CreateFeature(hRasterLayer, hFeat);
582 : /* Query raster ID to set it as the ID of the associated metadata */
583 4 : int nRasterID = (int)OGR_F_GetFID(hFeat);
584 :
585 4 : OGR_F_Destroy(hFeat);
586 :
587 4 : VSIUnlink(osTempFileName.c_str());
588 :
589 : /* -------------------------------------------------------------------- */
590 : /* Insert new entry into metadata table */
591 : /* -------------------------------------------------------------------- */
592 :
593 4 : hFeat = OGR_F_Create( OGR_L_GetLayerDefn(hMetadataLayer) );
594 4 : OGR_F_SetFID(hFeat, nRasterID);
595 4 : OGR_F_SetFieldString(hFeat, 0, osSourceName);
596 4 : OGR_F_SetFieldInteger(hFeat, 1, nTileId ++);
597 4 : OGR_F_SetFieldInteger(hFeat, 2, nReqXSize);
598 4 : OGR_F_SetFieldInteger(hFeat, 3, nReqYSize);
599 4 : OGR_F_SetFieldDouble(hFeat, 4, dfXResolution);
600 4 : OGR_F_SetFieldDouble(hFeat, 5, dfYResolution);
601 :
602 : double minx, maxx, maxy, miny;
603 4 : minx = adfGeoTransform[0] +
604 4 : (nBlockXSize * nBlockXOff) * dfXResolution;
605 4 : maxx = adfGeoTransform[0] +
606 4 : (nBlockXSize * nBlockXOff + nReqXSize) * dfXResolution;
607 4 : maxy = adfGeoTransform[3] +
608 4 : (nBlockYSize * nBlockYOff) * (-dfYResolution);
609 4 : miny = adfGeoTransform[3] +
610 4 : (nBlockYSize * nBlockYOff + nReqYSize) * (-dfYResolution);
611 :
612 4 : OGRGeometryH hRectangle = OGR_G_CreateGeometry(wkbPolygon);
613 4 : OGRGeometryH hLinearRing = OGR_G_CreateGeometry(wkbLinearRing);
614 4 : OGR_G_AddPoint_2D(hLinearRing, minx, miny);
615 4 : OGR_G_AddPoint_2D(hLinearRing, minx, maxy);
616 4 : OGR_G_AddPoint_2D(hLinearRing, maxx, maxy);
617 4 : OGR_G_AddPoint_2D(hLinearRing, maxx, miny);
618 4 : OGR_G_AddPoint_2D(hLinearRing, minx, miny);
619 4 : OGR_G_AddGeometryDirectly(hRectangle, hLinearRing);
620 :
621 4 : OGR_F_SetGeometryDirectly(hFeat, hRectangle);
622 :
623 4 : OGR_L_CreateFeature(hMetadataLayer, hFeat);
624 4 : OGR_F_Destroy(hFeat);
625 :
626 4 : nBlocks++;
627 4 : if (pfnProgress && !pfnProgress(1.0 * nBlocks / nTotalBlocks,
628 : NULL, pProgressData))
629 0 : eErr = CE_Failure;
630 : }
631 : }
632 :
633 4 : nLimitOvrCount = -1;
634 :
635 4 : if (eErr == CE_None)
636 4 : OGR_DS_ExecuteSQL(hDS, "COMMIT", NULL, NULL);
637 : else
638 0 : OGR_DS_ExecuteSQL(hDS, "ROLLBACK", NULL, NULL);
639 :
640 4 : VSIFree(pabyMEMDSBuffer);
641 4 : VSIFree(pabyPrevOvrMEMDSBuffer);
642 :
643 4 : CSLDestroy(papszTileDriverOptions);
644 4 : papszTileDriverOptions = NULL;
645 :
646 : /* -------------------------------------------------------------------- */
647 : /* Update raster_pyramids table */
648 : /* -------------------------------------------------------------------- */
649 4 : if (eErr == CE_None)
650 : {
651 4 : OGRLayerH hRasterPyramidsLyr = OGR_DS_GetLayerByName(hDS, "raster_pyramids");
652 4 : if (hRasterPyramidsLyr == NULL)
653 : {
654 : osSQL.Printf ("CREATE TABLE raster_pyramids ("
655 : "table_prefix TEXT NOT NULL,"
656 : "pixel_x_size DOUBLE NOT NULL,"
657 : "pixel_y_size DOUBLE NOT NULL,"
658 1 : "tile_count INTEGER NOT NULL)");
659 1 : OGR_DS_ExecuteSQL(hDS, osSQL.c_str(), NULL, NULL);
660 :
661 : /* Re-open the DB to take into account the new tables*/
662 1 : OGRReleaseDataSource(hDS);
663 :
664 1 : CPLString osOldVal = CPLGetConfigOption("SQLITE_LIST_ALL_TABLES", "FALSE");
665 1 : CPLSetThreadLocalConfigOption("SQLITE_LIST_ALL_TABLES", "TRUE");
666 1 : hDS = OGROpen(osFileName.c_str(), TRUE, NULL);
667 1 : CPLSetThreadLocalConfigOption("SQLITE_LIST_ALL_TABLES", osOldVal.c_str());
668 :
669 1 : hRasterPyramidsLyr = OGR_DS_GetLayerByName(hDS, "raster_pyramids");
670 1 : if (hRasterPyramidsLyr == NULL)
671 0 : return CE_Failure;
672 : }
673 4 : OGRFeatureDefnH hFDefn = OGR_L_GetLayerDefn(hRasterPyramidsLyr);
674 :
675 : /* Insert base resolution into raster_pyramids if not already done */
676 4 : int bHasBaseResolution = FALSE;
677 : osSQL.Printf("SELECT * FROM raster_pyramids WHERE "
678 : "table_prefix = '%s' AND pixel_x_size >= %.15f AND pixel_x_size <= %.15f AND "
679 : "pixel_y_size >= %.15f AND pixel_y_size <= %.15f",
680 : osTableName.c_str(),
681 8 : padfXResolutions[0] - 1e-15, padfXResolutions[0] + 1e-15,
682 12 : padfYResolutions[0] - 1e-15, padfYResolutions[0] + 1e-15);
683 4 : hSQLLyr = OGR_DS_ExecuteSQL(hDS, osSQL.c_str(), NULL, NULL);
684 4 : if (hSQLLyr)
685 : {
686 4 : OGRFeatureH hFeat = OGR_L_GetNextFeature(hSQLLyr);
687 4 : if (hFeat)
688 : {
689 3 : bHasBaseResolution = TRUE;
690 3 : OGR_F_Destroy(hFeat);
691 : }
692 4 : OGR_DS_ReleaseResultSet(hDS, hSQLLyr);
693 : }
694 :
695 4 : if (!bHasBaseResolution)
696 : {
697 : osSQL.Printf("SELECT COUNT(*) FROM \"%s\" WHERE "
698 : "pixel_x_size >= %.15f AND pixel_x_size <= %.15f AND "
699 : "pixel_y_size >= %.15f AND pixel_y_size <= %.15f",
700 : osMetatadataLayer.c_str(),
701 2 : padfXResolutions[0] - 1e-15, padfXResolutions[0] + 1e-15,
702 3 : padfYResolutions[0] - 1e-15, padfYResolutions[0] + 1e-15);
703 :
704 1 : int nBlocksMainRes = 0;
705 :
706 1 : hSQLLyr = OGR_DS_ExecuteSQL(hDS, osSQL.c_str(), NULL, NULL);
707 1 : if (hSQLLyr)
708 : {
709 1 : OGRFeatureH hFeat = OGR_L_GetNextFeature(hSQLLyr);
710 1 : if (hFeat)
711 : {
712 1 : nBlocksMainRes = OGR_F_GetFieldAsInteger(hFeat, 0);
713 1 : OGR_F_Destroy(hFeat);
714 : }
715 1 : OGR_DS_ReleaseResultSet(hDS, hSQLLyr);
716 : }
717 :
718 1 : OGRFeatureH hFeat = OGR_F_Create( hFDefn );
719 1 : OGR_F_SetFieldString(hFeat, OGR_FD_GetFieldIndex(hFDefn, "table_prefix"), osTableName.c_str());
720 1 : OGR_F_SetFieldDouble(hFeat, OGR_FD_GetFieldIndex(hFDefn, "pixel_x_size"), padfXResolutions[0]);
721 1 : OGR_F_SetFieldDouble(hFeat, OGR_FD_GetFieldIndex(hFDefn, "pixel_y_size"), padfYResolutions[0]);
722 1 : OGR_F_SetFieldInteger(hFeat, OGR_FD_GetFieldIndex(hFDefn, "tile_count"), nBlocksMainRes);
723 1 : OGR_L_CreateFeature(hRasterPyramidsLyr, hFeat);
724 1 : OGR_F_Destroy(hFeat);
725 : }
726 :
727 4 : OGRFeatureH hFeat = OGR_F_Create( hFDefn );
728 4 : OGR_F_SetFieldString(hFeat, OGR_FD_GetFieldIndex(hFDefn, "table_prefix"), osTableName.c_str());
729 4 : OGR_F_SetFieldDouble(hFeat, OGR_FD_GetFieldIndex(hFDefn, "pixel_x_size"), dfXResolution);
730 4 : OGR_F_SetFieldDouble(hFeat, OGR_FD_GetFieldIndex(hFDefn, "pixel_y_size"), dfYResolution);
731 4 : OGR_F_SetFieldInteger(hFeat, OGR_FD_GetFieldIndex(hFDefn, "tile_count"), nTotalBlocks);
732 4 : OGR_L_CreateFeature(hRasterPyramidsLyr, hFeat);
733 4 : OGR_F_Destroy(hFeat);
734 : }
735 :
736 4 : return eErr;
737 : }
738 :
739 : /************************************************************************/
740 : /* IBuildOverviews() */
741 : /************************************************************************/
742 :
743 4 : CPLErr RasterliteDataset::IBuildOverviews( const char * pszResampling,
744 : int nOverviews, int * panOverviewList,
745 : int nBands, int * panBandList,
746 : GDALProgressFunc pfnProgress,
747 : void * pProgressData )
748 : {
749 4 : CPLErr eErr = CE_None;
750 :
751 4 : if (nLevel != 0)
752 : {
753 : CPLError(CE_Failure, CPLE_AppDefined,
754 0 : "Overviews can only be computed on the base dataset");
755 0 : return CE_Failure;
756 : }
757 :
758 4 : if (osTableName.size() == 0)
759 0 : return CE_Failure;
760 :
761 : /* -------------------------------------------------------------------- */
762 : /* If we don't have read access, then create the overviews */
763 : /* externally. */
764 : /* -------------------------------------------------------------------- */
765 4 : if( GetAccess() != GA_Update )
766 : {
767 : CPLDebug( "Rasterlite",
768 : "File open for read-only accessing, "
769 0 : "creating overviews externally." );
770 :
771 0 : if (nResolutions != 1)
772 : {
773 : CPLError(CE_Failure, CPLE_NotSupported,
774 : "Cannot add external overviews to a "
775 0 : "dataset with internal overviews");
776 0 : return CE_Failure;
777 : }
778 :
779 0 : bCheckForExistingOverview = FALSE;
780 : eErr = GDALDataset::IBuildOverviews(
781 : pszResampling, nOverviews, panOverviewList,
782 0 : nBands, panBandList, pfnProgress, pProgressData );
783 0 : bCheckForExistingOverview = TRUE;
784 0 : return eErr;
785 : }
786 :
787 : /* -------------------------------------------------------------------- */
788 : /* If zero overviews were requested, we need to clear all */
789 : /* existing overviews. */
790 : /* -------------------------------------------------------------------- */
791 4 : if (nOverviews == 0)
792 : {
793 2 : return CleanOverviews();
794 : }
795 :
796 2 : if( nBands != GetRasterCount() )
797 : {
798 : CPLError( CE_Failure, CPLE_NotSupported,
799 : "Generation of overviews in RASTERLITE only"
800 : " supported when operating on all bands.\n"
801 0 : "Operation failed.\n" );
802 0 : return CE_Failure;
803 : }
804 :
805 2 : const char* pszOvrOptions = CPLGetConfigOption("RASTERLITE_OVR_OPTIONS", NULL);
806 2 : char** papszOptions = (pszOvrOptions) ? CSLTokenizeString2( pszOvrOptions, ",", 0) : NULL;
807 2 : GDALValidateCreationOptions( GetDriver(), papszOptions);
808 :
809 : int i;
810 6 : for(i=0;i<nOverviews && eErr == CE_None;i++)
811 : {
812 4 : if (panOverviewList[i] <= 1)
813 0 : continue;
814 :
815 4 : eErr = CleanOverviewLevel(panOverviewList[i]);
816 4 : if (eErr == CE_None)
817 4 : eErr = CreateOverviewLevel(pszResampling, panOverviewList[i], papszOptions, pfnProgress, pProgressData);
818 :
819 4 : ReloadOverviews();
820 : }
821 :
822 2 : CSLDestroy(papszOptions);
823 :
824 2 : return eErr;
825 : }
|