1 : /*************************************************************************
2 : * File : postgisrasterdataset.cpp
3 : * Project: PostGIS Raster driver
4 : * Purpose: GDAL Dataset implementation for PostGIS Raster driver
5 : * Author: Jorge Arevalo, jorge.arevalo@deimos-space.com
6 : *
7 : * Last changes:
8 : * $Id:$
9 : *
10 : ************************************************************************
11 : * Copyright (c) 2009 - 2011, Jorge Arevalo, jorge.arevalo@deimos-space.com
12 : *
13 : * Permission is hereby granted, free of charge, to any person obtaining a
14 : * copy of this software and associated documentation files (the
15 : * "Software"), to deal in the Software without restriction, including
16 : * without limitation the rights to use, copy, modify, merge, publish,
17 : * distribute, sublicense, and/or sell copies of the Software, and to
18 : * permit persons to whom the Software is furnished to do so, subject to
19 : * the following conditions:
20 : *
21 : * The above copyright notice and this permission notice shall be included
22 : * in all copies or substantial portions of the Software.
23 : *
24 : * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
25 : * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
26 : * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT.
27 : * IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY
28 : * CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT,
29 : * TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE
30 : * SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
31 : ************************************************************************/
32 :
33 : #include "postgisraster.h"
34 : #include <stdlib.h>
35 : #include "ogr_api.h"
36 : #include "ogr_geometry.h"
37 : #include "gdal.h"
38 : #include "cpl_conv.h"
39 : #include "cpl_string.h"
40 : #include "gdal_priv.h"
41 : #include <math.h>
42 : #include "cpl_error.h"
43 : #include "ogr_core.h"
44 :
45 : #ifdef _WIN32
46 : #define rint(x) floor((x) + 0.5)
47 : #endif
48 :
49 : CPL_C_START
50 : void GDALRegister_PostGISRaster(void);
51 :
52 : CPL_C_END
53 :
54 :
55 :
56 : /************************
57 : * \brief Constructor
58 : ************************/
59 0 : PostGISRasterDataset::PostGISRasterDataset() {
60 0 : papszSubdatasets = NULL;
61 0 : nSrid = -1;
62 0 : poConn = NULL;
63 0 : bRegularBlocking = false;
64 0 : bRegisteredInRasterColumns = false;
65 0 : pszSchema = NULL;
66 0 : pszTable = NULL;
67 0 : pszColumn = NULL;
68 0 : pszWhere = NULL;
69 0 : pszProjection = NULL;
70 0 : nMode = NO_MODE;
71 0 : poDriver = NULL;
72 0 : nBlockXSize = 0;
73 0 : nBlockYSize = 0;
74 0 : adfGeoTransform[0] = 0.0; /* X Origin (top left corner) */
75 0 : adfGeoTransform[1] = 1.0; /* X Pixel size */
76 0 : adfGeoTransform[2] = 0.0;
77 0 : adfGeoTransform[3] = 0.0; /* Y Origin (top left corner) */
78 0 : adfGeoTransform[4] = 0.0;
79 0 : adfGeoTransform[5] = 1.0; /* Y Pixel Size */
80 0 : bBlocksCached = false;
81 0 : bRegularBlocking = false;
82 0 : bAllTilesSnapToSameGrid = false;
83 :
84 : /**
85 : * TODO: Parametrize bAllTilesSnapToSameGrid. It controls if all the
86 : * raster rows, in ONE_RASTER_PER_TABLE mode, must be checked to test if
87 : * they snap to the same grid and have the same srid. It can be the user
88 : * decission, if he/she's sure all the rows pass the test and want more
89 : * speed.
90 : **/
91 :
92 0 : }
93 :
94 : /************************
95 : * \brief Constructor
96 : ************************/
97 0 : PostGISRasterDataset::~PostGISRasterDataset() {
98 :
99 0 : if (pszSchema)
100 0 : CPLFree(pszSchema);
101 0 : if (pszTable)
102 0 : CPLFree(pszTable);
103 0 : if (pszColumn)
104 0 : CPLFree(pszColumn);
105 0 : if (pszWhere)
106 0 : CPLFree(pszWhere);
107 0 : if (pszProjection)
108 0 : CPLFree(pszProjection);
109 :
110 0 : if (papszSubdatasets)
111 0 : CSLDestroy(papszSubdatasets);
112 :
113 0 : PQfinish(poConn);
114 :
115 : //delete poDriver;
116 0 : }
117 :
118 : /**************************************************************
119 : * \brief Replace the single quotes by " in the input string
120 : *
121 : * Needed before tokenize function
122 : *************************************************************/
123 : static
124 0 : char * ReplaceSingleQuotes(const char * pszInput, int nLength) {
125 : int i;
126 0 : char* pszOutput = NULL;
127 :
128 0 : if (nLength == -1)
129 0 : nLength = strlen(pszInput);
130 :
131 0 : pszOutput = (char*) CPLCalloc(nLength + 1, sizeof (char));
132 :
133 0 : for (i = 0; i < nLength; i++) {
134 0 : if (pszInput[i] == '\'')
135 0 : pszOutput[i] = '"';
136 : else
137 0 : pszOutput[i] = pszInput[i];
138 :
139 : }
140 :
141 0 : return pszOutput;
142 : }
143 :
144 : /**************************************************************
145 : * \brief Replace the quotes by single quotes in the input string
146 : *
147 : * Needed in the 'where' part of the input string
148 : *************************************************************/
149 : static
150 0 : char * ReplaceQuotes(const char * pszInput, int nLength) {
151 : int i;
152 0 : char * pszOutput = NULL;
153 :
154 0 : if (nLength == -1)
155 0 : nLength = strlen(pszInput);
156 :
157 0 : pszOutput = (char*) CPLCalloc(nLength + 1, sizeof (char));
158 :
159 0 : for (i = 0; i < nLength; i++) {
160 0 : if (pszInput[i] == '"')
161 0 : pszOutput[i] = '\'';
162 : else
163 0 : pszOutput[i] = pszInput[i];
164 : }
165 :
166 0 : return pszOutput;
167 : }
168 :
169 : /*****************************************************************************
170 : * \brief Split connection string into user, password, host, database...
171 : *
172 : * The parameters separated by spaces are return as a list of strings. The
173 : * function accepts all the PostgreSQL recognized parameter key words.
174 : *
175 : * The returned list must be freed with CSLDestroy when no longer needed
176 : *
177 : *****************************************************************************/
178 : static
179 0 : char** ParseConnectionString(const char * pszConnectionString) {
180 0 : char * pszEscapedConnectionString = NULL;
181 :
182 : /* Escape string following SQL scheme */
183 0 : pszEscapedConnectionString = ReplaceSingleQuotes(pszConnectionString, -1);
184 :
185 : /* Avoid PG: part */
186 0 : char* pszStartPos = (char*) strstr(pszEscapedConnectionString, ":") + 1;
187 :
188 : /* Tokenize */
189 : char** papszParams = CSLTokenizeString2(pszStartPos, " ",
190 0 : CSLT_HONOURSTRINGS);
191 :
192 : /* Free */
193 0 : CPLFree(pszEscapedConnectionString);
194 :
195 0 : return papszParams;
196 :
197 : }
198 :
199 : /**************************************************************************
200 : * \brief Look for raster tables in database and store them as subdatasets
201 : *
202 : * If no table is provided in connection string, the driver looks for the
203 : * existent raster tables in the schema given as argument. This argument,
204 : * however, is optional. If a NULL value is provided, the driver looks for
205 : * all raster tables in all schemas of the user-provided database.
206 : *
207 : * NOTE: Permissions are managed by libpq. The driver only returns an error
208 : * if an error is returned when trying to access to tables not allowed to
209 : * the current user.
210 : **************************************************************************/
211 0 : GBool PostGISRasterDataset::BrowseDatabase(const char* pszCurrentSchema,
212 : char* pszValidConnectionString) {
213 : /* Be careful! These 3 vars override the class ones! */
214 0 : char* pszSchema = NULL;
215 0 : char* pszTable = NULL;
216 0 : char* pszColumn = NULL;
217 0 : int i = 0;
218 0 : int nTuples = 0;
219 0 : PGresult * poResult = NULL;
220 0 : CPLString osCommand;
221 :
222 :
223 : /*************************************************************
224 : * Fetch all the raster tables and store them as subdatasets
225 : *************************************************************/
226 0 : if (pszCurrentSchema == NULL) {
227 : osCommand.Printf("select pg_namespace.nspname as schema, pg_class.relname as \
228 : table, pg_attribute.attname as column from pg_class, \
229 : pg_namespace,pg_attribute, pg_type where \
230 : pg_class.relnamespace = pg_namespace.oid and pg_class.oid = \
231 : pg_attribute.attrelid and pg_attribute.atttypid = pg_type.oid \
232 0 : and pg_type.typname = 'raster'");
233 :
234 0 : poResult = PQexec(poConn, osCommand.c_str());
235 0 : if (
236 : poResult == NULL ||
237 : PQresultStatus(poResult) != PGRES_TUPLES_OK ||
238 : PQntuples(poResult) <= 0
239 : ) {
240 : CPLError(CE_Failure, CPLE_AppDefined,
241 0 : "Error browsing database for PostGIS Raster tables: %s", PQerrorMessage(poConn));
242 0 : if (poResult != NULL)
243 0 : PQclear(poResult);
244 :
245 0 : return false;
246 : }
247 :
248 :
249 0 : nTuples = PQntuples(poResult);
250 0 : for (i = 0; i < nTuples; i++) {
251 0 : pszSchema = PQgetvalue(poResult, i, 0);
252 0 : pszTable = PQgetvalue(poResult, i, 1);
253 0 : pszColumn = PQgetvalue(poResult, i, 2);
254 :
255 : papszSubdatasets = CSLSetNameValue(papszSubdatasets,
256 : CPLSPrintf("SUBDATASET_%d_NAME", (i + 1)),
257 : CPLSPrintf("PG:%s schema=%s table=%s column=%s",
258 0 : pszValidConnectionString, pszSchema, pszTable, pszColumn));
259 :
260 : papszSubdatasets = CSLSetNameValue(papszSubdatasets,
261 : CPLSPrintf("SUBDATASET_%d_DESC", (i + 1)),
262 0 : CPLSPrintf("PostGIS Raster table at %s.%s (%s)", pszSchema, pszTable, pszColumn));
263 : }
264 :
265 0 : PQclear(poResult);
266 :
267 : } /**********************************************************************
268 : * Fetch all the schema's raster tables and store them as subdatasets
269 : **********************************************************************/
270 : else {
271 : osCommand.Printf("select pg_class.relname as table, pg_attribute.attname \
272 : as column from pg_class, pg_namespace,pg_attribute, pg_type where \
273 : pg_class.relnamespace = pg_namespace.oid and pg_class.oid = \
274 : pg_attribute.attrelid and pg_attribute.atttypid = pg_type.oid \
275 : and pg_type.typname = 'raster' and pg_namespace.nspname = '%s'",
276 0 : pszCurrentSchema);
277 :
278 0 : poResult = PQexec(poConn, osCommand.c_str());
279 0 : if (
280 : poResult == NULL ||
281 : PQresultStatus(poResult) != PGRES_TUPLES_OK ||
282 : PQntuples(poResult) <= 0
283 : ) {
284 : CPLError(CE_Failure, CPLE_AppDefined,
285 0 : "Error browsing database for PostGIS Raster tables: %s", PQerrorMessage(poConn));
286 0 : if (poResult != NULL)
287 0 : PQclear(poResult);
288 :
289 0 : return false;
290 : }
291 :
292 :
293 0 : nTuples = PQntuples(poResult);
294 0 : for (i = 0; i < nTuples; i++) {
295 0 : pszTable = PQgetvalue(poResult, i, 0);
296 0 : pszColumn = PQgetvalue(poResult, i, 1);
297 :
298 : papszSubdatasets = CSLSetNameValue(papszSubdatasets,
299 : CPLSPrintf("SUBDATASET_%d_NAME", (i + 1)),
300 : CPLSPrintf("PG:%s schema=%s table=%s column=%s",
301 0 : pszValidConnectionString, pszCurrentSchema, pszTable, pszColumn));
302 :
303 : papszSubdatasets = CSLSetNameValue(papszSubdatasets,
304 : CPLSPrintf("SUBDATASET_%d_DESC", (i + 1)),
305 : CPLSPrintf("PostGIS Raster table at %s.%s (%s)", pszCurrentSchema,
306 0 : pszTable, pszColumn));
307 : }
308 :
309 0 : PQclear(poResult);
310 : }
311 :
312 0 : return true;
313 : }
314 :
315 : /*************************************************************************
316 : * \brief Set the general raster properties.
317 : *
318 : * We must distinguish between tiled and untiled raster coverages. In
319 : * PostGIS Raster, there's no real difference between 'tile' and 'raster'.
320 : * There's only 'raster objects'. Each record of a raster table is a
321 : * raster object, and has its own georeference information, whether if
322 : * the record is a tile of a bigger raster coverage or is a complete
323 : * raster. So, <b>there's no a way of knowing if the rows of a raster
324 : * table are related or not</b>. It's user's responsibility. The only
325 : * thing driver can do is to suppose all the rows of a table are from
326 : * the same raster coverage if the user has queried for one table, without
327 : * specifying a where clause.
328 : *
329 : * The user is responsible to ensure that the raster layer meets the minimum
330 : * topological requirements for analysis. The ideal case is when all the raster
331 : * tiles of a continuous layer are the same size, snap to the same grid and do
332 : * not overlap.
333 : *
334 : * So, when we query for a raster table, we have 3 different cases:
335 : * - If the result is only one row, we can gather the raster properties
336 : * from the returned object, regardless is a tile or a whole raster
337 : * - If the result are several rows of a table, and the working mode is
338 : * ONE_RASTER_PER_TABLE, we assume all the rows are from the same raster
339 : * coverage. The rows are ordered by upper left y, upper left x, growing way,
340 : * and we can get raster size from the first and last elements.
341 : * - If the result are several rows of a table, and the working mode is
342 : * ONE_RASTER_PER_ROW, we assume each row is a different raster object,
343 : * and is reported as a subdataset. If you want only one of the raster rows,
344 : * you must specify a where clause to restrict the number of rows returned.
345 : **************************************************************************/
346 0 : GBool PostGISRasterDataset::SetRasterProperties
347 : (const char * pszValidConnectionString)
348 : {
349 0 : PGresult* poResult = NULL;
350 0 : CPLString osCommand;
351 0 : CPLString osCommand2;
352 0 : PGresult* poResult2 = NULL;
353 0 : int i = 0;
354 0 : int nTuples = 0;
355 0 : GBool bRetValue = false;
356 0 : OGRSpatialReference * poSR = NULL;
357 0 : OGRGeometry* poGeom = NULL;
358 0 : OGRErr OgrErr = OGRERR_NONE;
359 0 : OGREnvelope* poE = NULL;
360 : char* pszExtent;
361 : char* pszProjectionRef;
362 0 : int nTmpSrid = -1;
363 0 : double dfTmpScaleX = 0.0;
364 0 : double dfTmpScaleY = 0.0;
365 0 : double dfTmpSkewX = 0.0;
366 0 : double dfTmpSkewY = 0.0;
367 0 : int nWidth = 0;
368 0 : int nHeight = 0;
369 0 : int nTmpWidth = 0;
370 0 : int nTmpHeight = 0;
371 :
372 :
373 : /* Execute the query to fetch raster data from db */
374 0 : if (pszWhere == NULL) {
375 : osCommand.Printf("select (foo.md).*, foo.rid from (select rid, st_metadata(%s) as md \
376 0 : from %s.%s) as foo", pszColumn, pszSchema, pszTable);
377 : } else {
378 : osCommand.Printf("select (foo.md).*, foo.rid from (select rid, st_metadata(%s) as md \
379 0 : from %s.%s where %s) as foo", pszColumn, pszSchema, pszTable, pszWhere);
380 : }
381 :
382 : CPLDebug("PostGIS_Raster", "PostGISRasterDataset::SetRasterProperties(): "
383 0 : "Query: %s", osCommand.c_str());
384 0 : poResult = PQexec(poConn, osCommand.c_str());
385 0 : if (
386 : poResult == NULL ||
387 : PQresultStatus(poResult) != PGRES_TUPLES_OK ||
388 : PQntuples(poResult) <= 0
389 : ) {
390 : CPLError(CE_Failure, CPLE_AppDefined,
391 0 : "Error browsing database for PostGIS Raster properties");
392 0 : if (poResult != NULL)
393 0 : PQclear(poResult);
394 :
395 0 : return false;
396 : }
397 :
398 0 : nTuples = PQntuples(poResult);
399 :
400 :
401 : /******************************************
402 : * Easier case. Only one raster to fetch
403 : ******************************************/
404 0 : if (nTuples == 1) {
405 0 : nSrid = atoi(PQgetvalue(poResult, 0, 8));
406 0 : nBands = atoi(PQgetvalue(poResult, 0, 9));
407 0 : adfGeoTransform[0] = atof(PQgetvalue(poResult, 0, 0)); //upperleft x
408 0 : adfGeoTransform[1] = atof(PQgetvalue(poResult, 0, 4)); //pixelsize x
409 0 : adfGeoTransform[2] = atof(PQgetvalue(poResult, 0, 6)); //skew x
410 0 : adfGeoTransform[3] = atof(PQgetvalue(poResult, 0, 1)); //upperleft y
411 0 : adfGeoTransform[4] = atof(PQgetvalue(poResult, 0, 7)); //skew y
412 0 : adfGeoTransform[5] = atof(PQgetvalue(poResult, 0, 5)); //pixelsize y
413 :
414 0 : nRasterXSize = atoi(PQgetvalue(poResult, 0, 2));
415 0 : nRasterYSize = atoi(PQgetvalue(poResult, 0, 3));
416 :
417 : /**
418 : * Not tiled dataset: The whole raster.
419 : * TODO: 'invent' a good block size.
420 : */
421 0 : nBlockXSize = nRasterXSize;
422 0 : nBlockYSize = nRasterYSize;
423 :
424 0 : bRetValue = true;
425 : } else {
426 0 : switch (nMode) {
427 : /**
428 : * Each row is a different raster. Create subdatasets, one per row
429 : **/
430 :
431 : case ONE_RASTER_PER_ROW:
432 : {
433 :
434 0 : for (i = 0; i < nTuples; i++) {
435 0 : nSrid = atoi(PQgetvalue(poResult, i, 10));
436 :
437 : papszSubdatasets = CSLSetNameValue(papszSubdatasets,
438 : CPLSPrintf("SUBDATASET_%d_NAME", (i + 1)),
439 : CPLSPrintf("PG:%s schema=%s table=%s column=%s where='rid = %d'",
440 0 : pszValidConnectionString, pszSchema, pszTable, pszColumn, nSrid));
441 :
442 : papszSubdatasets = CSLSetNameValue(papszSubdatasets,
443 : CPLSPrintf("SUBDATASET_%d_DESC", (i + 1)),
444 : CPLSPrintf("PostGIS Raster at %s.%s (%s), rid = %d", pszSchema,
445 0 : pszTable, pszColumn, nSrid));
446 : }
447 :
448 : /* Not a single raster fetched */
449 0 : nRasterXSize = 0;
450 0 : nRasterYSize = 0;
451 :
452 0 : bRetValue = true;
453 :
454 : }
455 0 : break;
456 :
457 : /************************************************************
458 : * All rows form a whole raster coverage
459 : ************************************************************/
460 : case ONE_RASTER_PER_TABLE:
461 : {
462 :
463 : /**
464 : * Get the rest of raster properties from this object
465 : */
466 0 : nSrid = atoi(PQgetvalue(poResult, 0, 8));
467 :
468 0 : nBands = atoi(PQgetvalue(poResult, 0, 9));
469 0 : adfGeoTransform[0] = atof(PQgetvalue(poResult, 0, 0)); //upperleft x
470 0 : adfGeoTransform[1] = atof(PQgetvalue(poResult, 0, 4)); //pixelsize x
471 0 : adfGeoTransform[2] = atof(PQgetvalue(poResult, 0, 6)); //skew x
472 0 : adfGeoTransform[3] = atof(PQgetvalue(poResult, 0, 1)); //upperleft y
473 0 : adfGeoTransform[4] = atof(PQgetvalue(poResult, 0, 7)); //skew y
474 0 : adfGeoTransform[5] = atof(PQgetvalue(poResult, 0, 5)); //pixelsize y
475 0 : nWidth = atoi(PQgetvalue(poResult, 0, 2));
476 0 : nHeight = atoi(PQgetvalue(poResult, 0, 3));
477 :
478 : /**
479 : * Now check if all tiles have the same dimensions.
480 : *
481 : * NOTE: If bRegularBlocking is 'true', this is not checked.
482 : * It's user responsibility
483 : *
484 : * TODO: Find a good block size, that works in any situation.
485 : **/
486 0 : if (!bRegularBlocking)
487 : {
488 0 : for(i = 1; i < nTuples; i++)
489 : {
490 0 : nTmpWidth = atoi(PQgetvalue(poResult, i, 2));
491 0 : nTmpHeight = atoi(PQgetvalue(poResult, i, 3));
492 :
493 0 : if (nWidth != nTmpWidth || nHeight != nTmpHeight)
494 : {
495 : // Not supported until the above TODO is implemented
496 : CPLError(CE_Failure, CPLE_AppDefined,
497 : "Error, the table %s.%s contains tiles with "
498 : "different size, and irregular blocking is "
499 0 : "not supported yet", pszSchema, pszTable);
500 :
501 0 : PQclear(poResult);
502 0 : return false;
503 : }
504 : }
505 :
506 : // Now, we can ensure this
507 0 : bRegularBlocking = true;
508 0 : nBlockXSize = nWidth;
509 0 : nBlockYSize = nHeight;
510 : }
511 :
512 : // The user ensures this...
513 : else
514 : {
515 0 : nBlockXSize = nWidth;
516 0 : nBlockYSize = nHeight;
517 : }
518 :
519 : /**
520 : * Check all the raster tiles have the same srid and snap to the
521 : * same grid. If not, return an error
522 : *
523 : * NOTE: If bAllTilesSnapToSameGrid is 'true', this is not
524 : * checked. It's user responsibility.
525 : *
526 : * TODO: Work even if this requisites are not complained. For
527 : * example, by:
528 : * - Resampling all the rows to the grid of the first one
529 : * - Providing a new grid alignment for all the rows, with a
530 : * maximum of 6 parameters: ulx, uly, pixelsizex, pixelsizey,
531 : * skewx, skewy or a minimum of 3 parameters: ulx, uly,
532 : * pixelsize (x and y pixel sizes are equal and both skew are
533 : * 0).
534 : **/
535 0 : if (!bAllTilesSnapToSameGrid)
536 : {
537 0 : for(i = 1; i < nTuples; i++)
538 : {
539 0 : nTmpSrid = atoi(PQgetvalue(poResult, i, 8));
540 0 : dfTmpScaleX = atof(PQgetvalue(poResult, i, 4));
541 0 : dfTmpScaleY = atof(PQgetvalue(poResult, i, 5));
542 0 : dfTmpSkewX = atof(PQgetvalue(poResult, i, 6));
543 0 : dfTmpSkewY = atof(PQgetvalue(poResult, i, 7));
544 :
545 0 : if (nTmpSrid != nSrid ||
546 0 : FLT_NEQ(dfTmpScaleX, adfGeoTransform[1]) ||
547 0 : FLT_NEQ(dfTmpScaleY, adfGeoTransform[5]) ||
548 0 : FLT_NEQ(dfTmpSkewX, adfGeoTransform[2]) ||
549 0 : FLT_NEQ(dfTmpSkewY, adfGeoTransform[4]))
550 : {
551 : /**
552 : * In this mode, it is not allowed this situation,
553 : * unless while the above TODO is not implemented
554 : **/
555 : CPLError(CE_Failure, CPLE_AppDefined,
556 : "Error, the table %s.%s contains tiles with "
557 : "different SRID or snapping to different grids",
558 0 : pszSchema, pszTable);
559 :
560 0 : PQclear(poResult);
561 0 : return false;
562 : }
563 : }
564 :
565 : // Now, we can ensure this
566 0 : bAllTilesSnapToSameGrid = true;
567 : }
568 :
569 : /**
570 : * Now, if there's irregular blocking and/or the blocks don't
571 : * snap to the same grid or don't have the same srid, we should
572 : * fix these situations. Assuming that we don't return an error
573 : * in that cases, of course.
574 : **/
575 :
576 :
577 :
578 : /**
579 : * Get whole raster extent
580 : **/
581 0 : if (pszWhere == NULL)
582 : osCommand2.Printf("select st_astext(st_setsrid(st_extent(%s::geometry),%d)) from %s.%s",
583 0 : pszColumn, nSrid, pszSchema, pszTable);
584 : else
585 : osCommand2.Printf("select st_astext(st_setsrid(st_extent(%s::geometry),%d)) from %s.%s where %s",
586 0 : pszColumn, nSrid, pszSchema, pszTable, pszWhere);
587 :
588 :
589 0 : poResult2 = PQexec(poConn, osCommand2.c_str());
590 0 : if (poResult2 == NULL ||
591 : PQresultStatus(poResult2) != PGRES_TUPLES_OK ||
592 : PQntuples(poResult2) <= 0) {
593 : CPLError(CE_Failure, CPLE_AppDefined,
594 : "Error calculating whole raster extent: %s",
595 0 : PQerrorMessage(poConn));
596 :
597 0 : if (poResult2 != NULL)
598 0 : PQclear(poResult2);
599 :
600 0 : if (poResult != NULL)
601 0 : PQclear(poResult);
602 :
603 0 : return false;
604 : }
605 :
606 : /* Construct an OGR object with the raster extent */
607 0 : pszExtent = PQgetvalue(poResult2, 0, 0);
608 :
609 0 : pszProjectionRef = (char*) GetProjectionRef();
610 0 : poSR = new OGRSpatialReference(pszProjectionRef);
611 : OgrErr = OGRGeometryFactory::createFromWkt(&pszExtent,
612 0 : poSR, &poGeom);
613 0 : if (OgrErr != OGRERR_NONE) {
614 : CPLError(CE_Failure, CPLE_AppDefined,
615 0 : "Couldn't calculate raster extent");
616 :
617 0 : if (poResult2)
618 0 : PQclear(poResult2);
619 :
620 0 : if (poResult != NULL)
621 0 : PQclear(poResult);
622 :
623 0 : return false;
624 : }
625 :
626 0 : poE = new OGREnvelope();
627 0 : poGeom->getEnvelope(poE);
628 :
629 : /* Correction for upper left y coord*/
630 :
631 : /**
632 : * TODO: Review this. Is a good algorithm?
633 : * If the pixel size Y is negative, we can assume the raster's
634 : * reference system uses cartesian coordinates, in which the
635 : * origin is in lower-left corner, while the origin in an image
636 : * is un upper-left corner. In this case, the upper left Y value
637 : * will be MaxY from the envelope. Otherwise, it will be MinY.
638 : **/
639 : /*
640 : adfGeoTransform[0] = poE->MinX;
641 : if (adfGeoTransform[5] >= 0.0)
642 : adfGeoTransform[3] = poE->MinY;
643 : else
644 : adfGeoTransform[3] = poE->MaxY;
645 : */
646 :
647 : /**
648 : * The raster size is the extent covered for all the raster's
649 : * columns
650 : **/
651 : nRasterXSize = (int)
652 0 : fabs(rint((poE->MaxX - poE->MinX) / adfGeoTransform[1]));
653 : nRasterYSize = (int)
654 0 : fabs(rint((poE->MaxY - poE->MinY) / adfGeoTransform[5]));
655 :
656 :
657 : /* Free resources */
658 0 : OGRGeometryFactory::destroyGeometry(poGeom);
659 0 : delete poE;
660 0 : delete poSR;
661 0 : PQclear(poResult2);
662 :
663 0 : bRetValue = true;
664 :
665 : }
666 0 : break;
667 :
668 : /* TODO: take into account more cases, if applies */
669 : default:
670 : {
671 : CPLError(CE_Failure, CPLE_AppDefined,
672 0 : "Error, incorrect working mode");
673 :
674 0 : bRetValue = false;
675 : }
676 : }
677 : }
678 :
679 : CPLDebug("PostGIS_Raster", "PostGISRasterDataset::SetRasterProperties(): "
680 : "adfGeoTransform = {%f, %f, %f, %f, %f,%f}", adfGeoTransform[0],
681 : adfGeoTransform[1], adfGeoTransform[2], adfGeoTransform[3],
682 0 : adfGeoTransform[4], adfGeoTransform[5]);
683 :
684 : CPLDebug("PostGIS_Raster", "PostGISRasterDataset::SetRasterProperties(): "
685 0 : "Raster size = (%d, %d)", nRasterXSize, nRasterYSize);
686 :
687 :
688 : CPLDebug("PostGIS_Raster", "PostGISRasterDataset::SetRasterProperties(): "
689 0 : "Block dimensions = (%d x %d)", nBlockXSize, nBlockYSize);
690 :
691 0 : PQclear(poResult);
692 0 : return bRetValue;
693 : }
694 :
695 : /*********************************************
696 : * \brief Set raster bands for this dataset
697 : *********************************************/
698 0 : GBool PostGISRasterDataset::SetRasterBands() {
699 0 : GBool bSignedByte = false;
700 0 : int nBitDepth = 8;
701 0 : char* pszDataType = NULL;
702 0 : int iBand = 0;
703 0 : PGresult * poResult = NULL;
704 0 : CPLString osCommand;
705 0 : double dfNodata = 0.0;
706 0 : GDALDataType hDataType = GDT_Byte;
707 0 : int nTuples = 0;
708 0 : GBool bIsOffline = false;
709 :
710 : /* Create each PostGISRasterRasterBand using the band metadata */
711 0 : for (iBand = 0; iBand < nBands; iBand++) {
712 : /* Create query to fetch metadata from db */
713 0 : if (pszWhere == NULL) {
714 : osCommand.Printf("select (foo.md).* from (select"
715 : " distinct st_bandmetadata( %s, %d) as md from %s. %s) as foo",
716 0 : pszColumn, iBand + 1, pszSchema, pszTable);
717 : } else {
718 :
719 : osCommand.Printf("select (foo.md).* from (select"
720 : " distinct st_bandmetadata( %s, %d) as md from %s. %s where %s) as foo",
721 0 : pszColumn, iBand + 1, pszSchema, pszTable, pszWhere);
722 : }
723 :
724 0 : poResult = PQexec(poConn, osCommand.c_str());
725 0 : nTuples = PQntuples(poResult);
726 :
727 : /* Error getting info from database */
728 0 : if (poResult == NULL || PQresultStatus(poResult) != PGRES_TUPLES_OK ||
729 : nTuples <= 0) {
730 : CPLError(CE_Failure, CPLE_AppDefined,
731 : "Error getting band metadata: %s",
732 0 : PQerrorMessage(poConn));
733 0 : if (poResult)
734 0 : PQclear(poResult);
735 :
736 0 : return false;
737 : }
738 :
739 : /**
740 : * If we have more than one record here is because there are several
741 : * rows, belonging to the same raster coverage, with different band
742 : * metadata values. An error must be raised.
743 : *
744 : * TODO: Is there any way to fix this problem?
745 : *
746 : * TODO: Even when the difference between metadata values are only a
747 : * few decimal numbers (for example: 3.0000000 and 3.0000001) they're
748 : * different tuples. And in that case, they must be the same
749 : **/
750 : /*
751 : if (nTuples > 1) {
752 : CPLError(CE_Failure, CPLE_AppDefined, "Error, the \
753 : ONE_RASTER_PER_TABLE mode can't be applied if the raster \
754 : rows don't have the same metadata for band %d",
755 : iBand + 1);
756 : PQclear(poResult);
757 : return false;
758 : }
759 : */
760 :
761 :
762 : /* Get metadata and create raster band objects */
763 0 : pszDataType = CPLStrdup(PQgetvalue(poResult, 0, 0));
764 0 : dfNodata = atof(PQgetvalue(poResult, 0, 2));
765 0 : bIsOffline = EQUALN(PQgetvalue(poResult, 0, 3), "t", sizeof (char));
766 :
767 0 : if (EQUALN(pszDataType, "1BB", 3 * sizeof (char))) {
768 0 : hDataType = GDT_Byte;
769 0 : nBitDepth = 1;
770 0 : } else if (EQUALN(pszDataType, "2BUI", 4 * sizeof (char))) {
771 0 : hDataType = GDT_Byte;
772 0 : nBitDepth = 2;
773 0 : } else if (EQUALN(pszDataType, "4BUI", 4 * sizeof (char))) {
774 0 : hDataType = GDT_Byte;
775 0 : nBitDepth = 4;
776 0 : } else if (EQUALN(pszDataType, "8BUI", 4 * sizeof (char))) {
777 0 : hDataType = GDT_Byte;
778 0 : nBitDepth = 8;
779 0 : } else if (EQUALN(pszDataType, "8BSI", 4 * sizeof (char))) {
780 0 : hDataType = GDT_Byte;
781 : /**
782 : * To indicate the unsigned byte values between 128 and 255
783 : * should be interpreted as being values between -128 and -1 for
784 : * applications that recognise the SIGNEDBYTE type.
785 : **/
786 0 : bSignedByte = true;
787 0 : nBitDepth = 8;
788 0 : } else if (EQUALN(pszDataType, "16BSI", 5 * sizeof (char))) {
789 0 : hDataType = GDT_Int16;
790 0 : nBitDepth = 16;
791 0 : } else if (EQUALN(pszDataType, "16BUI", 5 * sizeof (char))) {
792 0 : hDataType = GDT_UInt16;
793 0 : nBitDepth = 16;
794 0 : } else if (EQUALN(pszDataType, "32BSI", 5 * sizeof (char))) {
795 0 : hDataType = GDT_Int32;
796 0 : nBitDepth = 32;
797 0 : } else if (EQUALN(pszDataType, "32BUI", 5 * sizeof (char))) {
798 0 : hDataType = GDT_UInt32;
799 0 : nBitDepth = 32;
800 0 : } else if (EQUALN(pszDataType, "32BF", 4 * sizeof (char))) {
801 0 : hDataType = GDT_Float32;
802 0 : nBitDepth = 32;
803 0 : } else if (EQUALN(pszDataType, "64BF", 4 * sizeof (char))) {
804 0 : hDataType = GDT_Float64;
805 0 : nBitDepth = 64;
806 : } else {
807 0 : hDataType = GDT_Byte;
808 0 : nBitDepth = 8;
809 : }
810 :
811 : /* Create raster band object */
812 : SetBand(iBand + 1, new PostGISRasterRasterBand(this, iBand + 1, hDataType,
813 0 : dfNodata, bSignedByte, nBitDepth, 0, bIsOffline));
814 :
815 0 : CPLFree(pszDataType);
816 0 : PQclear(poResult);
817 : }
818 :
819 0 : return true;
820 : }
821 :
822 : /**
823 : * Read/write a region of image data from multiple bands.
824 : *
825 : * This method allows reading a region of one or more PostGISRasterBands from
826 : * this dataset into a buffer. The write support is still under development
827 : *
828 : * The function fetches all the raster data that intersects with the region
829 : * provided, and store the data in the GDAL cache.
830 : *
831 : * TODO: This only works in case of regular blocking rasters. A more
832 : * general approach to allow non-regular blocking rasters is under development.
833 : *
834 : * It automatically takes care of data type translation if the data type
835 : * (eBufType) of the buffer is different than that of the
836 : * PostGISRasterRasterBand.
837 : *
838 : * TODO: The method should take care of image decimation / replication if the
839 : * buffer size (nBufXSize x nBufYSize) is different than the size of the region
840 : * being accessed (nXSize x nYSize).
841 : *
842 : * The nPixelSpace, nLineSpace and nBandSpace parameters allow reading into or
843 : * writing from various organization of buffers.
844 : *
845 : * @param eRWFlag Either GF_Read to read a region of data, or GF_Write to write
846 : * a region of data.
847 : *
848 : * @param nXOff The pixel offset to the top left corner of the region of the
849 : * band to be accessed. This would be zero to start from the left side.
850 : *
851 : * @param nYOff The line offset to the top left corner of the region of the band
852 : * to be accessed. This would be zero to start from the top.
853 : *
854 : * @param nXSize The width of the region of the band to be accessed in pixels.
855 : *
856 : * @param nYSize The height of the region of the band to be accessed in lines.
857 : *
858 : * @param pData The buffer into which the data should be read, or from which it
859 : * should be written. This buffer must contain at least
860 : * nBufXSize * nBufYSize * nBandCount words of type eBufType. It is organized in
861 : * left to right,top to bottom pixel order. Spacing is controlled by the
862 : * nPixelSpace, and nLineSpace parameters.
863 : *
864 : * @param nBufXSize the width of the buffer image into which the desired region
865 : * is to be read, or from which it is to be written.
866 : *
867 : * @param nBufYSize the height of the buffer image into which the desired region
868 : * is to be read, or from which it is to be written.
869 : *
870 : * @param eBufType the type of the pixel values in the pData data buffer. The
871 : * pixel values will automatically be translated to/from the
872 : * PostGISRasterRasterBand data type as needed.
873 : *
874 : * @param nBandCount the number of bands being read or written.
875 : *
876 : * @param panBandMap the list of nBandCount band numbers being read/written.
877 : * Note band numbers are 1 based. This may be NULL to select the first
878 : * nBandCount bands.
879 : *
880 : * @param nPixelSpace The byte offset from the start of one pixel value in pData
881 : * to the start of the next pixel value within a scanline. If defaulted (0) the
882 : * size of the datatype eBufType is used.
883 : *
884 : * @param nLineSpace The byte offset from the start of one scanline in pData to
885 : * the start of the next. If defaulted (0) the size of the datatype
886 : * eBufType * nBufXSize is used.
887 : *
888 : * @param nBandSpace the byte offset from the start of one bands data to the
889 : * start of the next. If defaulted (0) the value will be nLineSpace * nBufYSize
890 : * implying band sequential organization of the data buffer.
891 : *
892 : * @return CE_Failure if the access fails, otherwise CE_None.
893 : */
894 0 : CPLErr PostGISRasterDataset::IRasterIO(GDALRWFlag eRWFlag,
895 : int nXOff, int nYOff, int nXSize, int nYSize,
896 : void * pData, int nBufXSize, int nBufYSize,
897 : GDALDataType eBufType,
898 : int nBandCount, int *panBandMap,
899 : int nPixelSpace, int nLineSpace, int nBandSpace)
900 : {
901 : double adfTransform[6];
902 : double adfProjWin[8];
903 : int ulx, uly, lrx, lry;
904 0 : CPLString osCommand;
905 0 : PGresult* poResult = NULL;
906 0 : int nTuples = 0;
907 0 : int iTuplesIndex = 0;
908 0 : GByte* pbyData = NULL;
909 0 : int nWKBLength = 0;
910 : int iBandIndex;
911 0 : GDALRasterBlock * poBlock = NULL;
912 : int iBlockXOff, iBlockYOff;
913 : int nBandDataSize, nBandDataLength;
914 0 : char * pBandData = NULL;
915 0 : PostGISRasterRasterBand * poBand = NULL;
916 0 : GByte * pabySrcBlock = NULL;
917 : int nBlocksPerRow, nBlocksPerColumn;
918 : char orderByY[4];
919 : char orderByX[3];
920 : int rid;
921 :
922 :
923 : /**
924 : * TODO: Write support not implemented yet
925 : **/
926 0 : if (eRWFlag == GF_Write)
927 : {
928 : CPLError(CE_Failure, CPLE_NotSupported,
929 0 : "PostGIS Raster does not support writing");
930 0 : return CE_Failure;
931 : }
932 :
933 : /**
934 : * TODO: Data decimation / replication needed
935 : */
936 0 : if (nBufXSize != nXSize || nBufYSize != nYSize)
937 : {
938 : /**
939 : * This will cause individual IReadBlock calls
940 : *
941 : */
942 : return GDALDataset::IRasterIO(eRWFlag, nXOff, nYOff, nXSize, nYSize,
943 : pData, nBufXSize, nBufYSize, eBufType, nBandCount,
944 0 : panBandMap, nPixelSpace, nLineSpace, nBandSpace);
945 : }
946 :
947 : CPLDebug("PostGIS_Raster", "PostGISRasterDataset::IRasterIO: "
948 : "nBandSpace = %d, nLineSpace = %d, nPixelSpace = %d",
949 0 : nBandSpace, nLineSpace, nPixelSpace);
950 :
951 : /**************************************************************************
952 : * In the first call, we fetch the data from database and store it as
953 : * 'blocks' in the cache.
954 : *
955 : * TODO: If the data is not cached, we must 'invent' a good block size, and
956 : * divide the data in blocks. To get a proper block size, we should rewrite
957 : * the GetBlockSize function at band level.
958 : ***************************************************************************/
959 0 : if (!bBlocksCached)
960 : {
961 : CPLDebug("PostGIS_Raster", "PostGISRasterDataset::IRasterIO: "
962 : "Buffer size = (%d, %d), Region size = (%d, %d)",
963 0 : nBufXSize, nBufYSize, nXSize, nYSize);
964 :
965 : /*******************************************************************
966 : * Construct a projected window to intersect the band data
967 : *******************************************************************/
968 0 : GetGeoTransform(adfTransform);
969 0 : ulx = nXOff;
970 0 : uly = nYOff;
971 0 : lrx = nXOff + nXSize;
972 0 : lry = nYOff + nYSize;
973 :
974 : /* Calculate the intersection polygon */
975 0 : adfProjWin[0] = adfTransform[0] +
976 0 : ulx * adfTransform[1] +
977 0 : uly * adfTransform[2];
978 :
979 0 : adfProjWin[1] = adfTransform[3] +
980 0 : ulx * adfTransform[4] +
981 0 : uly * adfTransform[5];
982 :
983 0 : adfProjWin[2] = adfTransform[0] +
984 0 : lrx * adfTransform[1] +
985 0 : uly * adfTransform[2];
986 :
987 0 : adfProjWin[3] = adfTransform[3] +
988 0 : lrx * adfTransform[4] +
989 0 : uly * adfTransform[5];
990 :
991 0 : adfProjWin[4] = adfTransform[0] +
992 0 : lrx * adfTransform[1] +
993 0 : lry * adfTransform[2];
994 :
995 0 : adfProjWin[5] = adfTransform[3] +
996 0 : lrx * adfTransform[4] +
997 0 : lry * adfTransform[5];
998 :
999 0 : adfProjWin[6] = adfTransform[0] +
1000 0 : ulx * adfTransform[1] +
1001 0 : lry * adfTransform[2];
1002 :
1003 0 : adfProjWin[7] = adfTransform[3] +
1004 0 : ulx * adfTransform[4] +
1005 0 : lry * adfTransform[5];
1006 :
1007 :
1008 : /* Construct order by for the query */
1009 0 : memset(orderByX, 0, 3);
1010 0 : memset(orderByY, 0, 4);
1011 :
1012 0 : strcpy(orderByX, "asc");
1013 0 : if (nSrid == -1)
1014 0 : strcpy(orderByY, "asc"); // Y starts at 0 and grows
1015 : else
1016 0 : strcpy(orderByY, "desc");// Y starts at max and decreases
1017 :
1018 :
1019 : /*********************************************************************
1020 : * We first get the data from database (ordered from upper left pixel
1021 : * to lower right one)
1022 : *********************************************************************/
1023 0 : if (pszWhere == NULL)
1024 : {
1025 : osCommand.Printf("SELECT rid, %s, ST_ScaleX(%s), ST_SkewY(%s), "
1026 : "ST_SkewX(%s), ST_ScaleY(%s), ST_UpperLeftX(%s), "
1027 : "ST_UpperLeftY(%s), ST_Width(%s), ST_Height(%s) FROM %s.%s WHERE "
1028 : "ST_Intersects(%s, ST_PolygonFromText('POLYGON((%f %f, %f %f, %f %f"
1029 : ", %f %f, %f %f))', %d)) ORDER BY ST_UpperLeftY(%s) %s, "
1030 : "ST_UpperLeftX(%s) %s", pszColumn, pszColumn,
1031 : pszColumn, pszColumn, pszColumn, pszColumn, pszColumn, pszColumn,
1032 : pszColumn, pszSchema, pszTable, pszColumn, adfProjWin[0],
1033 : adfProjWin[1], adfProjWin[2], adfProjWin[3], adfProjWin[4],
1034 : adfProjWin[5], adfProjWin[6], adfProjWin[7], adfProjWin[0],
1035 0 : adfProjWin[1], nSrid, pszColumn, orderByY, pszColumn, orderByX);
1036 : }
1037 :
1038 :
1039 : else
1040 : {
1041 : osCommand.Printf("SELECT rid, %s ST_ScaleX(%s), ST_SkewY(%s), "
1042 : "ST_SkewX(%s), ST_ScaleY(%s), ST_UpperLeftX(%s), "
1043 : "ST_UpperLeftY(%s), ST_Width(%s), ST_Height(%s) FROM %s.%s WHERE %s AND "
1044 : "ST_Intersects(%s, ST_PolygonFromText('POLYGON((%f %f, %f %f, %f %f"
1045 : ", %f %f, %f %f))', %d)) ORDER BY ST_UpperLeftY(%s) %s, "
1046 : "ST_UpperLeftX(%s) %s", pszColumn, pszColumn,
1047 : pszColumn, pszColumn, pszColumn, pszColumn, pszColumn, pszColumn,
1048 : pszColumn,pszSchema, pszTable, pszWhere, pszColumn, adfProjWin[0],
1049 : adfProjWin[1], adfProjWin[2], adfProjWin[3], adfProjWin[4],
1050 : adfProjWin[5], adfProjWin[6], adfProjWin[7], adfProjWin[0],
1051 0 : adfProjWin[1], nSrid, pszColumn, orderByY, pszColumn, orderByX);
1052 : }
1053 :
1054 : CPLDebug("PostGIS_Raster", "PostGISRasterDataset::IRasterIO(): Query = %s",
1055 0 : osCommand.c_str());
1056 :
1057 0 : poResult = PQexec(poConn, osCommand.c_str());
1058 0 : if (poResult == NULL || PQresultStatus(poResult) != PGRES_TUPLES_OK ||
1059 : PQntuples(poResult) <= 0)
1060 : {
1061 0 : if (poResult)
1062 0 : PQclear(poResult);
1063 :
1064 : /**
1065 : * This will cause individual IReadBlock calls
1066 : *
1067 : */
1068 : return GDALDataset::IRasterIO(eRWFlag, nXOff, nYOff, nXSize, nYSize,
1069 : pData, nBufXSize, nBufYSize, eBufType, nBandCount,
1070 0 : panBandMap, nPixelSpace, nLineSpace, nBandSpace);
1071 : }
1072 :
1073 : /**
1074 : * NOTE: In case of any of the raster columns have different SRID, the
1075 : * query will fail. So, if we don't fail, we can assume all the rows
1076 : * have the same SRID. We don't need to check it
1077 : **/
1078 :
1079 :
1080 0 : nTuples = PQntuples(poResult);
1081 : CPLDebug("PostGIS_Raster", "PostGISRasterDataset::IRasterIO(): nTuples = %d",
1082 0 : nTuples);
1083 :
1084 : CPLDebug("PostGIS_Raster", "PostGISRasterDataset::IRasterIO(): "
1085 0 : "Raster size = (%d, %d)", nRasterXSize, nRasterYSize);
1086 :
1087 :
1088 :
1089 : /**************************************************************************
1090 : * This is the simplest case: all the rows have the same dimensions
1091 : * (regularly blocked raster)
1092 : *
1093 : * Now we'll cache each tuple as a data block. More accurately, we must
1094 : * access each tuple, get the band data, and store this data as a block. So,
1095 : * each tuple contains the data for nBands blocks (and nBandCount <=
1096 : * nBands)
1097 : *************************************************************************/
1098 0 : for(iBandIndex = 0; iBandIndex < nBandCount; iBandIndex++)
1099 : {
1100 0 : poBand = (PostGISRasterRasterBand *)GetRasterBand(iBandIndex + 1);
1101 :
1102 0 : nBandDataSize = GDALGetDataTypeSize(poBand->eDataType) / 8;
1103 :
1104 : nBandDataLength = poBand->nBlockXSize * poBand->nBlockYSize *
1105 0 : nBandDataSize;
1106 :
1107 : CPLDebug("PostGIS_Raster", "PostGISRasterDataset::IRasterIO(): "
1108 : "Block size (%d, %d) for band %d", poBand->nBlockXSize,
1109 0 : poBand->nBlockYSize, poBand->nBand);
1110 :
1111 : /* Enables block caching, if it wasn't enabled */
1112 0 : if (!poBand->InitBlockInfo())
1113 0 : continue;
1114 :
1115 : /**
1116 : * This can be different from poBand->nBlocksPerRow and
1117 : * poBand->nBlocksPerColumn, if the region size is different than
1118 : * the raster size. So, we calculate these values for this case.
1119 : */
1120 : nBlocksPerRow =
1121 0 : (nXSize + poBand->nBlockXSize - 1) / poBand->nBlockXSize;
1122 :
1123 : nBlocksPerColumn =
1124 0 : (nYSize + poBand->nBlockYSize - 1) / poBand->nBlockYSize;
1125 :
1126 : CPLDebug("PostGIS_Raster", "PostGISRasterDataset::IRasterIO(): "
1127 0 : "Number of blocks: %dx%d", nBlocksPerRow, nBlocksPerColumn);
1128 :
1129 0 : for(iBlockYOff = 0; iBlockYOff < nBlocksPerColumn;
1130 : iBlockYOff++)
1131 : {
1132 0 : for(iBlockXOff = 0; iBlockXOff < nBlocksPerRow;
1133 : iBlockXOff++)
1134 : {
1135 : iTuplesIndex = (iBlockYOff * nBlocksPerRow) +
1136 0 : iBlockXOff;
1137 :
1138 : CPLDebug("PostGIS_Raster", "PostGISRasterDataset::IRasterIO(): "
1139 : "iBlockXOff = %d, iBlockYOff = %d, "
1140 : "iTuplesIndex = %d", iBlockXOff, iBlockYOff,
1141 0 : iTuplesIndex);
1142 :
1143 : pbyData = CPLHexToBinary(PQgetvalue(poResult, iTuplesIndex,
1144 0 : 1), &nWKBLength);
1145 :
1146 : pBandData = (char *)GET_BAND_DATA(pbyData, poBand->nBand,
1147 0 : nBandDataSize, nBandDataLength);
1148 :
1149 : CPLDebug("PostGIS_Raster", "PostGISRasterDataset::IRasterIO(): "
1150 : "Block data length for band %d: %d", poBand->nBand,
1151 0 : nBandDataLength);
1152 :
1153 : CPLDebug("PostGIS_Raster", "PostGISRasterDataset::IRasterIO(): "
1154 0 : "Block (%d, %d)", iBlockXOff, iBlockYOff);
1155 :
1156 : /* Create a new block */
1157 : poBlock = new GDALRasterBlock(poBand, iBlockXOff,
1158 0 : iBlockYOff);
1159 :
1160 0 : poBlock->AddLock();
1161 :
1162 : /* Allocate data space */
1163 0 : if (poBlock->Internalize() != CE_None)
1164 : {
1165 0 : poBlock->DropLock();
1166 0 : delete poBlock;
1167 0 : continue;
1168 : }
1169 :
1170 : /* Add the block to the block matrix */
1171 0 : if (poBand->AdoptBlock(iBlockXOff, iBlockYOff, poBlock) !=
1172 : CE_None)
1173 : {
1174 0 : poBlock->DropLock();
1175 0 : delete poBlock;
1176 0 : continue;
1177 : }
1178 :
1179 : /**
1180 : * Copy data to block
1181 : *
1182 : * TODO: Enable write mode too (mark the block as dirty and
1183 : * create IWriteBlock in PostGISRasterRasterBand)
1184 : */
1185 0 : pabySrcBlock = (GByte *)poBlock->GetDataRef();
1186 :
1187 0 : if (poBand->eDataType == eBufType)
1188 : {
1189 0 : memcpy(pabySrcBlock, pBandData, nBandDataLength);
1190 : }
1191 :
1192 : /**
1193 : * As in GDALDataset class... expensive way of handling
1194 : * single words
1195 : */
1196 : else
1197 : {
1198 : GDALCopyWords(pBandData, poBand->eDataType, 0,
1199 0 : pabySrcBlock, eBufType, 0, 1);
1200 : }
1201 :
1202 0 : poBlock->DropLock();
1203 :
1204 0 : CPLFree(pbyData);
1205 0 : pbyData = NULL;
1206 : }
1207 :
1208 : }
1209 :
1210 : }
1211 :
1212 0 : PQclear(poResult);
1213 0 : bBlocksCached = true;
1214 : }
1215 :
1216 : /* Once the blocks are cached, we delegate in GDAL I/O system */
1217 : return GDALDataset::IRasterIO(eRWFlag, nXOff, nYOff, nXSize, nYSize, pData,
1218 : nBufXSize, nBufYSize, eBufType, nBandCount, panBandMap, nPixelSpace,
1219 0 : nLineSpace, nBandSpace);
1220 : }
1221 :
1222 : /******************************************************************************
1223 : * \brief Open a connection with PostgreSQL. The connection string will have
1224 : * the PostgreSQL accepted format, plus the next key=value pairs:
1225 : * schema = <schema_name>
1226 : * table = <table_name>
1227 : * column = <column_name>
1228 : * where = <SQL where>
1229 : * mode = <working mode> (1 or 2)
1230 : *
1231 : * These pairs are used for selecting the right raster table.
1232 : *****************************************************************************/
1233 10648 : GDALDataset* PostGISRasterDataset::Open(GDALOpenInfo* poOpenInfo) {
1234 10648 : char** papszParams = NULL;
1235 10648 : char* pszConnectionString = NULL;
1236 10648 : char* pszValidConnectionString = NULL;
1237 10648 : char* pszSchema = NULL;
1238 10648 : char* pszTable = NULL;
1239 10648 : char* pszColumn = NULL;
1240 10648 : char* pszWhere = NULL;
1241 10648 : char* pszTmp = NULL;
1242 10648 : int nMode = -1;
1243 10648 : int nPos = -1;
1244 10648 : int i = 0;
1245 10648 : PGconn * poConn = NULL;
1246 10648 : PostGISRasterDataset* poDS = NULL;
1247 10648 : GBool bBrowseDatabase = false;
1248 10648 : PGresult * poResult = NULL;
1249 10648 : CPLString osCommand;
1250 :
1251 : /**************************
1252 : * Check input parameter
1253 : **************************/
1254 10648 : if (poOpenInfo->pszFilename == NULL ||
1255 : poOpenInfo->fp != NULL ||
1256 : !EQUALN(poOpenInfo->pszFilename, "PG:", 3) ||
1257 : (papszParams = ParseConnectionString(poOpenInfo->pszFilename)) == NULL) {
1258 : /**
1259 : * Drivers must quietly return NULL if the passed file is not of
1260 : * their format. They should only produce an error if the file
1261 : * does appear to be of their supported format, but for some
1262 : * reason, unsupported or corrupt
1263 : */
1264 10648 : return NULL;
1265 : }
1266 :
1267 : /**************************
1268 : * Create driver instance
1269 : **************************/
1270 : /*
1271 : poDriver = (PostGISRasterDriver*) GDALGetDriverByName("PostGISRaster");
1272 : if (poDriver == NULL) {
1273 : CPLError(CE_Failure, CPLE_AppDefined,
1274 : "PostGIS Raster driver not registered");
1275 : return NULL;
1276 : }
1277 : */
1278 :
1279 : /***************************************************************
1280 : * Now check existence of db, schema, table, column and where.
1281 : * If there's no enough data for querying one single table,
1282 : * we'll create a GDAL metadata's list of sub-datasets, with all
1283 : * the table with columns of type 'rasters'
1284 : ***************************************************************/
1285 :
1286 :
1287 : /**************************************************************************
1288 : * Get mode:
1289 : * - 1. ONE_RASTER_PER_ROW: Each row is considered as a separate raster
1290 : * - 2. ONE_RASTER_PER_TABLE: All the table rows are considered as a whole
1291 : * raster coverage
1292 : **************************************************************************/
1293 0 : nPos = CSLFindName(papszParams, "mode");
1294 0 : if (nPos != -1) {
1295 0 : nMode = atoi(CPLParseNameValue(papszParams[nPos], NULL));
1296 :
1297 0 : if (nMode != ONE_RASTER_PER_ROW && nMode != ONE_RASTER_PER_TABLE) {
1298 : /* Unrecognized mode, using default one */
1299 : /*
1300 : CPLError(CE_Warning, CPLE_AppDefined, "Undefined working mode (%d)."
1301 : " Valid working modes are 1 (ONE_RASTER_PER_ROW) and 2"
1302 : " (ONE_RASTER_PER_TABLE). Using ONE_RASTER_PER_TABLE"
1303 : " by default", nMode);
1304 : */
1305 0 : nMode = ONE_RASTER_PER_ROW;
1306 : }
1307 :
1308 : /* Remove the mode from connection string */
1309 0 : papszParams = CSLRemoveStrings(papszParams, nPos, 1, NULL);
1310 : }
1311 : /* Default mode */
1312 : else
1313 0 : nMode = ONE_RASTER_PER_ROW;
1314 :
1315 : /**
1316 : * Case 1: There's no database name: Error, you need, at least,
1317 : * specify a database name (NOTE: insensitive search)
1318 : **/
1319 0 : nPos = CSLFindName(papszParams, "dbname");
1320 0 : if (nPos == -1) {
1321 : CPLError(CE_Failure, CPLE_AppDefined,
1322 0 : "You must specify at least a db name");
1323 0 : return NULL;
1324 : }
1325 :
1326 : /**
1327 : * Case 2: There's database name, but no table name: activate a flag
1328 : * for browsing the database, fetching all the schemas that contain
1329 : * raster tables
1330 : **/
1331 0 : nPos = CSLFindName(papszParams, "table");
1332 0 : if (nPos == -1) {
1333 0 : bBrowseDatabase = true;
1334 :
1335 : /* Get schema name, if exist */
1336 0 : nPos = CSLFindName(papszParams, "schema");
1337 0 : if (nPos != -1) {
1338 0 : pszSchema = CPLStrdup(CPLParseNameValue(papszParams[nPos], NULL));
1339 : /* Delete this pair from params array */
1340 0 : papszParams = CSLRemoveStrings(papszParams, nPos, 1, NULL);
1341 : }
1342 :
1343 : /**
1344 : * Remove the rest of the parameters, if exist (they mustn't be present
1345 : * if we want a valid PQ connection string)
1346 : **/
1347 0 : nPos = CSLFindName(papszParams, "column");
1348 0 : if (nPos != -1) {
1349 : /* Delete this pair from params array */
1350 0 : papszParams = CSLRemoveStrings(papszParams, nPos, 1, NULL);
1351 : }
1352 :
1353 0 : nPos = CSLFindName(papszParams, "where");
1354 0 : if (nPos != -1) {
1355 : /* Delete this pair from params array */
1356 0 : papszParams = CSLRemoveStrings(papszParams, nPos, 1, NULL);
1357 : }
1358 : } else {
1359 0 : pszTable = CPLStrdup(CPLParseNameValue(papszParams[nPos], NULL));
1360 : /* Delete this pair from params array */
1361 0 : papszParams = CSLRemoveStrings(papszParams, nPos, 1, NULL);
1362 :
1363 : /**
1364 : * Case 3: There's database and table name, but no column
1365 : * name: Use a default column name and use the table to create the
1366 : * dataset
1367 : **/
1368 0 : nPos = CSLFindName(papszParams, "column");
1369 0 : if (nPos == -1) {
1370 0 : pszColumn = CPLStrdup(DEFAULT_COLUMN);
1371 : } /**
1372 : * Case 4: There's database, table and column name: Use the table to
1373 : * create a dataset
1374 : **/
1375 : else {
1376 0 : pszColumn = CPLStrdup(CPLParseNameValue(papszParams[nPos], NULL));
1377 : /* Delete this pair from params array */
1378 0 : papszParams = CSLRemoveStrings(papszParams, nPos, 1, NULL);
1379 : }
1380 :
1381 : /* Get the rest of the parameters, if exist */
1382 0 : nPos = CSLFindName(papszParams, "schema");
1383 0 : if (nPos == -1) {
1384 0 : pszSchema = CPLStrdup(DEFAULT_SCHEMA);
1385 : } else {
1386 0 : pszSchema = CPLStrdup(CPLParseNameValue(papszParams[nPos], NULL));
1387 : /* Delete this pair from params array */
1388 0 : papszParams = CSLRemoveStrings(papszParams, nPos, 1, NULL);
1389 : }
1390 :
1391 0 : nPos = CSLFindName(papszParams, "where");
1392 0 : if (nPos != -1) {
1393 0 : pszWhere = CPLStrdup(CPLParseNameValue(papszParams[nPos], NULL));
1394 : /* Delete this pair from params array */
1395 0 : papszParams = CSLRemoveStrings(papszParams, nPos, 1, NULL);
1396 : }
1397 :
1398 : }
1399 :
1400 : /* Parse pszWhere, if needed */
1401 0 : if (pszWhere) {
1402 0 : pszTmp = ReplaceQuotes(pszWhere, strlen(pszWhere));
1403 0 : CPLFree(pszWhere);
1404 0 : pszWhere = pszTmp;
1405 : }
1406 :
1407 :
1408 : /***************************************
1409 : * Construct a valid connection string
1410 : ***************************************/
1411 : pszValidConnectionString = (char*) CPLCalloc(strlen(poOpenInfo->pszFilename),
1412 0 : sizeof (char));
1413 0 : for (i = 0; i < CSLCount(papszParams); i++) {
1414 0 : pszValidConnectionString = strncat(pszValidConnectionString, papszParams[i], strlen(papszParams[i]));
1415 0 : pszValidConnectionString = strncat(pszValidConnectionString, " ", strlen(" "));
1416 :
1417 : //CPLStrlcat(pszValidConnectionString, papszParams[i], strlen(papszParams[i]));
1418 : //CPLStrlcat(pszValidConnectionString, " ", strlen(" "));
1419 :
1420 : }
1421 :
1422 :
1423 :
1424 : /********************************************************************
1425 : * Open a new database connection
1426 : * TODO: Use enviroment vars (PGHOST, PGPORT, PGUSER) instead of
1427 : * default values.
1428 : * TODO: USE DRIVER INSTEAD OF OPENNING A CONNECTION HERE. MEMORY
1429 : * PROBLEMS DETECTED (SEE DRIVER DESTRUCTOR FOR FURTHER INFORMATION)
1430 : ********************************************************************/
1431 : /*
1432 : poConn = poDriver->GetConnection(pszValidConnectionString,
1433 : CSLFetchNameValueDef(papszParams, "host", DEFAULT_HOST),
1434 : CSLFetchNameValueDef(papszParams, "port", DEFAULT_PORT),
1435 : CSLFetchNameValueDef(papszParams, "user", DEFAULT_USER),
1436 : CSLFetchNameValueDef(papszParams, "password", DEFAULT_PASSWORD));
1437 :
1438 : if (poConn == NULL) {
1439 : CPLError(CE_Failure, CPLE_AppDefined,
1440 : "Couldn't establish a database connection");
1441 : CSLDestroy(papszParams);
1442 : CPLFree(pszConnectionString);
1443 : CPLFree(pszValidConnectionString);
1444 : if (pszSchema)
1445 : CPLFree(pszSchema);
1446 : if (pszTable)
1447 : CPLFree(pszTable);
1448 : if (pszColumn)
1449 : CPLFree(pszColumn);
1450 : if (pszWhere)
1451 : CPLFree(pszWhere);
1452 :
1453 : delete poDriver;
1454 :
1455 : return NULL;
1456 : }
1457 : */
1458 :
1459 :
1460 : /* Frees no longer needed memory */
1461 0 : CSLDestroy(papszParams);
1462 0 : CPLFree(pszConnectionString);
1463 :
1464 : /**
1465 : * Get connection
1466 : * TODO: Try to get connection from poDriver
1467 : **/
1468 0 : poConn = PQconnectdb(pszValidConnectionString);
1469 0 : if (poConn == NULL) {
1470 : CPLError(CE_Failure, CPLE_AppDefined,
1471 0 : "Couldn't establish a database connection");
1472 0 : if (pszSchema)
1473 0 : CPLFree(pszSchema);
1474 0 : if (pszTable)
1475 0 : CPLFree(pszTable);
1476 0 : if (pszColumn)
1477 0 : CPLFree(pszColumn);
1478 0 : if (pszWhere)
1479 0 : CPLFree(pszWhere);
1480 :
1481 0 : return NULL;
1482 :
1483 : }
1484 :
1485 : /* Check geometry type existence */
1486 0 : poResult = PQexec(poConn, "SELECT oid FROM pg_type WHERE typname = 'geometry'");
1487 0 : if (
1488 : poResult == NULL ||
1489 : PQresultStatus(poResult) != PGRES_TUPLES_OK ||
1490 : PQntuples(poResult) <= 0
1491 : ) {
1492 : CPLError(CE_Failure, CPLE_AppDefined,
1493 : "Error checking geometry type existence. Is PostGIS correctly \
1494 0 : installed?: %s", PQerrorMessage(poConn));
1495 0 : if (poResult != NULL)
1496 0 : PQclear(poResult);
1497 0 : if (pszSchema)
1498 0 : CPLFree(pszSchema);
1499 0 : if (pszTable)
1500 0 : CPLFree(pszTable);
1501 0 : if (pszColumn)
1502 0 : CPLFree(pszColumn);
1503 0 : if (pszWhere)
1504 0 : CPLFree(pszWhere);
1505 :
1506 0 : PQfinish(poConn);
1507 :
1508 : //delete poDriver;
1509 :
1510 0 : return NULL;
1511 : }
1512 :
1513 0 : PQclear(poResult);
1514 :
1515 :
1516 : /* Check spatial tables existence */
1517 : poResult = PQexec(poConn, "select pg_namespace.nspname as schemaname, \
1518 : pg_class.relname as tablename from pg_class, \
1519 : pg_namespace where pg_class.relnamespace = pg_namespace.oid \
1520 : and (pg_class.relname='raster_columns' or \
1521 : pg_class.relname='raster_overviews' or \
1522 : pg_class.relname='geometry_columns' or \
1523 0 : pg_class.relname='spatial_ref_sys')");
1524 0 : if (
1525 : poResult == NULL ||
1526 : PQresultStatus(poResult) != PGRES_TUPLES_OK ||
1527 : PQntuples(poResult) <= 0
1528 : ) {
1529 : CPLError(CE_Failure, CPLE_AppDefined,
1530 : "Error checking needed tables existence: %s",
1531 0 : PQerrorMessage(poConn));
1532 0 : if (poResult != NULL)
1533 0 : PQclear(poResult);
1534 0 : if (pszSchema)
1535 0 : CPLFree(pszSchema);
1536 0 : if (pszTable)
1537 0 : CPLFree(pszTable);
1538 0 : if (pszColumn)
1539 0 : CPLFree(pszColumn);
1540 0 : if (pszWhere)
1541 0 : CPLFree(pszWhere);
1542 :
1543 0 : PQfinish(poConn);
1544 :
1545 : //delete poDriver;
1546 :
1547 0 : return NULL;
1548 :
1549 : }
1550 :
1551 0 : PQclear(poResult);
1552 :
1553 : /*************************************************************************
1554 : * No table will be read. Only shows information about the existent raster
1555 : * tables
1556 : *************************************************************************/
1557 0 : if (bBrowseDatabase) {
1558 : /**
1559 : * Creates empty dataset object, only for subdatasets
1560 : **/
1561 0 : poDS = new PostGISRasterDataset();
1562 0 : poDS->poConn = poConn;
1563 0 : poDS->eAccess = GA_ReadOnly;
1564 : //poDS->poDriver = poDriver;
1565 0 : poDS->nMode = (pszSchema) ? BROWSE_SCHEMA : BROWSE_DATABASE;
1566 0 : poDS->nRasterXSize = 0;
1567 0 : poDS->nRasterYSize = 0;
1568 0 : poDS->adfGeoTransform[0] = 0.0;
1569 0 : poDS->adfGeoTransform[1] = 0.0;
1570 0 : poDS->adfGeoTransform[2] = 0.0;
1571 0 : poDS->adfGeoTransform[3] = 0.0;
1572 0 : poDS->adfGeoTransform[4] = 0.0;
1573 0 : poDS->adfGeoTransform[5] = 0.0;
1574 :
1575 : /**
1576 : * Look for raster tables at database and
1577 : * store them as subdatasets
1578 : **/
1579 0 : if (!poDS->BrowseDatabase(pszSchema, pszValidConnectionString)) {
1580 0 : CPLFree(pszValidConnectionString);
1581 0 : delete poDS;
1582 0 : return NULL;
1583 : }
1584 : }
1585 : /***********************************************************************
1586 : * A table will be read: Fetch raster properties from db. Pay attention
1587 : * to the block size: if the raster is blocked at database, the block
1588 : * size can be fetched from each block size, if regular blocking table
1589 : **********************************************************************/
1590 : else {
1591 0 : poDS = new PostGISRasterDataset();
1592 0 : poDS->poConn = poConn;
1593 0 : poDS->eAccess = poOpenInfo->eAccess;
1594 0 : poDS->nMode = nMode;
1595 : //poDS->poDriver = poDriver;
1596 :
1597 0 : poDS->pszSchema = pszSchema;
1598 0 : poDS->pszTable = pszTable;
1599 0 : poDS->pszColumn = pszColumn;
1600 0 : poDS->pszWhere = pszWhere;
1601 :
1602 : /**
1603 : * Fetch basic raster metadata from db
1604 : **/
1605 :
1606 0 : if (!poDS->SetRasterProperties(pszValidConnectionString)) {
1607 0 : CPLFree(pszValidConnectionString);
1608 0 : delete poDS;
1609 0 : return NULL;
1610 : }
1611 :
1612 : /* Set raster bands */
1613 0 : if (!poDS->SetRasterBands()) {
1614 0 : CPLFree(pszValidConnectionString);
1615 0 : delete poDS;
1616 0 : return NULL;
1617 : }
1618 :
1619 : }
1620 :
1621 0 : CPLFree(pszValidConnectionString);
1622 0 : return poDS;
1623 :
1624 : }
1625 :
1626 : /*****************************************
1627 : * \brief Get Metadata from raster
1628 : * TODO: Add more options (the result of
1629 : * calling ST_Metadata, for example)
1630 : *****************************************/
1631 0 : char** PostGISRasterDataset::GetMetadata(const char *pszDomain) {
1632 0 : if (pszDomain != NULL && EQUALN(pszDomain, "SUBDATASETS", 11))
1633 0 : return papszSubdatasets;
1634 : else
1635 0 : return GDALDataset::GetMetadata(pszDomain);
1636 : }
1637 :
1638 : /*****************************************************
1639 : * \brief Fetch the projection definition string
1640 : * for this dataset in OpenGIS WKT format. It should
1641 : * be suitable for use with the OGRSpatialReference
1642 : * class.
1643 : *****************************************************/
1644 0 : const char* PostGISRasterDataset::GetProjectionRef() {
1645 0 : CPLString osCommand;
1646 : PGresult* poResult;
1647 :
1648 0 : if (nSrid == -1)
1649 0 : return "";
1650 :
1651 0 : if (pszProjection)
1652 0 : return pszProjection;
1653 :
1654 : /********************************************************
1655 : * Reading proj from database
1656 : ********************************************************/
1657 : osCommand.Printf("SELECT srtext FROM spatial_ref_sys where SRID=%d",
1658 0 : nSrid);
1659 0 : poResult = PQexec(this->poConn, osCommand.c_str());
1660 0 : if (poResult && PQresultStatus(poResult) == PGRES_TUPLES_OK
1661 : && PQntuples(poResult) > 0) {
1662 : /*
1663 : * TODO: Memory leak detected with valgrind caused by allocation here.
1664 : * Even when the return string is freed
1665 : */
1666 0 : pszProjection = CPLStrdup(PQgetvalue(poResult, 0, 0));
1667 : }
1668 :
1669 0 : if (poResult)
1670 0 : PQclear(poResult);
1671 :
1672 0 : return pszProjection;
1673 : }
1674 :
1675 : /**********************************************************
1676 : * \brief Set projection definition. The input string must
1677 : * be in OGC WKT or PROJ.4 format
1678 : **********************************************************/
1679 0 : CPLErr PostGISRasterDataset::SetProjection(const char * pszProjectionRef) {
1680 0 : VALIDATE_POINTER1(pszProjectionRef, "SetProjection", CE_Failure);
1681 :
1682 0 : CPLString osCommand;
1683 : PGresult * poResult;
1684 0 : int nFetchedSrid = -1;
1685 :
1686 :
1687 : /*****************************************************************
1688 : * Check if the dataset allows updating
1689 : *****************************************************************/
1690 0 : if (GetAccess() != GA_Update) {
1691 : CPLError(CE_Failure, CPLE_NoWriteAccess,
1692 0 : "This driver doesn't allow write access");
1693 0 : return CE_Failure;
1694 : }
1695 :
1696 : /*****************************************************************
1697 : * Look for projection with this text
1698 : *****************************************************************/
1699 :
1700 : // First, WKT text
1701 : osCommand.Printf("SELECT srid FROM spatial_ref_sys where srtext='%s'",
1702 0 : pszProjectionRef);
1703 0 : poResult = PQexec(poConn, osCommand.c_str());
1704 :
1705 0 : if (poResult && PQresultStatus(poResult) == PGRES_TUPLES_OK
1706 : && PQntuples(poResult) > 0) {
1707 :
1708 0 : nFetchedSrid = atoi(PQgetvalue(poResult, 0, 0));
1709 :
1710 : // update class attribute
1711 0 : nSrid = nFetchedSrid;
1712 :
1713 :
1714 : // update raster_columns table
1715 : osCommand.Printf("UPDATE raster_columns SET srid=%d WHERE \
1716 : r_table_name = '%s' AND r_column = '%s'",
1717 0 : nSrid, pszTable, pszColumn);
1718 0 : poResult = PQexec(poConn, osCommand.c_str());
1719 0 : if (poResult == NULL || PQresultStatus(poResult) != PGRES_COMMAND_OK) {
1720 : CPLError(CE_Failure, CPLE_AppDefined,
1721 : "Couldn't update raster_columns table: %s",
1722 0 : PQerrorMessage(poConn));
1723 0 : return CE_Failure;
1724 : }
1725 :
1726 : // TODO: Update ALL blocks with the new srid...
1727 :
1728 0 : return CE_None;
1729 : }
1730 : // If not, proj4 text
1731 : else {
1732 : osCommand.Printf(
1733 : "SELECT srid FROM spatial_ref_sys where proj4text='%s'",
1734 0 : pszProjectionRef);
1735 0 : poResult = PQexec(poConn, osCommand.c_str());
1736 :
1737 0 : if (poResult && PQresultStatus(poResult) == PGRES_TUPLES_OK
1738 : && PQntuples(poResult) > 0) {
1739 :
1740 0 : nFetchedSrid = atoi(PQgetvalue(poResult, 0, 0));
1741 :
1742 : // update class attribute
1743 0 : nSrid = nFetchedSrid;
1744 :
1745 : // update raster_columns table
1746 : osCommand.Printf("UPDATE raster_columns SET srid=%d WHERE \
1747 : r_table_name = '%s' AND r_column = '%s'",
1748 0 : nSrid, pszTable, pszColumn);
1749 :
1750 0 : poResult = PQexec(poConn, osCommand.c_str());
1751 0 : if (poResult == NULL ||
1752 : PQresultStatus(poResult) != PGRES_COMMAND_OK) {
1753 : CPLError(CE_Failure, CPLE_AppDefined,
1754 : "Couldn't update raster_columns table: %s",
1755 0 : PQerrorMessage(poConn));
1756 0 : return CE_Failure;
1757 : }
1758 :
1759 : // TODO: Update ALL blocks with the new srid...
1760 :
1761 0 : return CE_None;
1762 : }
1763 : else {
1764 : CPLError(CE_Failure, CPLE_WrongFormat,
1765 0 : "Couldn't find WKT neither proj4 definition");
1766 0 : return CE_Failure;
1767 : }
1768 0 : }
1769 : }
1770 :
1771 : /********************************************************
1772 : * \brief Set the affine transformation coefficients
1773 : ********************************************************/
1774 0 : CPLErr PostGISRasterDataset::SetGeoTransform(double* padfTransform) {
1775 0 : if (!padfTransform)
1776 0 : return CE_Failure;
1777 :
1778 0 : adfGeoTransform[0] = padfTransform[0];
1779 0 : adfGeoTransform[1] = padfTransform[1];
1780 0 : adfGeoTransform[2] = padfTransform[2];
1781 0 : adfGeoTransform[3] = padfTransform[3];
1782 0 : adfGeoTransform[4] = padfTransform[4];
1783 0 : adfGeoTransform[5] = padfTransform[5];
1784 :
1785 0 : return CE_None;
1786 : }
1787 :
1788 : /********************************************************
1789 : * \brief Get the affine transformation coefficients
1790 : ********************************************************/
1791 0 : CPLErr PostGISRasterDataset::GetGeoTransform(double * padfTransform) {
1792 : // copy necessary values in supplied buffer
1793 0 : padfTransform[0] = adfGeoTransform[0];
1794 0 : padfTransform[1] = adfGeoTransform[1];
1795 0 : padfTransform[2] = adfGeoTransform[2];
1796 0 : padfTransform[3] = adfGeoTransform[3];
1797 0 : padfTransform[4] = adfGeoTransform[4];
1798 0 : padfTransform[5] = adfGeoTransform[5];
1799 :
1800 0 : return CE_None;
1801 : }
1802 :
1803 :
1804 : /************************************************************************/
1805 : /* GDALRegister_PostGISRaster() */
1806 :
1807 : /************************************************************************/
1808 558 : void GDALRegister_PostGISRaster() {
1809 : GDALDriver *poDriver;
1810 :
1811 558 : if (GDALGetDriverByName("PostGISRaster") == NULL) {
1812 537 : poDriver = new GDALDriver();
1813 :
1814 537 : poDriver->SetDescription("PostGISRaster");
1815 : poDriver->SetMetadataItem(GDAL_DMD_LONGNAME,
1816 537 : "PostGIS Raster driver");
1817 :
1818 537 : poDriver->pfnOpen = PostGISRasterDataset::Open;
1819 :
1820 537 : GetGDALDriverManager()->RegisterDriver(poDriver);
1821 : }
1822 558 : }
1823 :
|
