1 : /******************************************************************************
2 : * $Id: ecrgtocdataset.cpp 23033 2011-09-03 18:46:11Z rouault $
3 : *
4 : * Project: ECRG TOC read Translator
5 : * Purpose: Implementation of ECRGTOCDataset and ECRGTOCSubDataset.
6 : * Author: Even Rouault, even.rouault at mines-paris.org
7 : *
8 : ******************************************************************************
9 : * Copyright (c) 2011, Even Rouault
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 : // g++ -g -Wall -fPIC frmts/nitf/ecrgtocdataset.cpp -shared -o gdal_ECRGTOC.so -Iport -Igcore -Iogr -Ifrmts/vrt -L. -lgdal
31 :
32 : #include "gdal_proxy.h"
33 : #include "ogr_srs_api.h"
34 : #include "vrtdataset.h"
35 : #include "cpl_minixml.h"
36 : #include <vector>
37 :
38 : CPL_CVSID("$Id: ecrgtocdataset.cpp 23033 2011-09-03 18:46:11Z rouault $");
39 :
40 : /** Overview of used classes :
41 : - ECRGTOCDataset : lists the different subdatasets, listed in the .xml,
42 : as subdatasets
43 : - ECRGTOCSubDataset : one of these subdatasets, implemented as a VRT, of
44 : the relevant NITF tiles
45 : - ECRGTOCProxyRasterDataSet : a "proxy" dataset that maps to a NITF tile
46 : */
47 :
48 : typedef struct
49 : {
50 : const char* pszName;
51 : const char* pszPath;
52 : int nScale;
53 : int nZone;
54 0 : } FrameDesc;
55 :
56 : /************************************************************************/
57 : /* ==================================================================== */
58 : /* ECRGTOCDataset */
59 : /* ==================================================================== */
60 : /************************************************************************/
61 :
62 : class ECRGTOCDataset : public GDALPamDataset
63 : {
64 : char **papszSubDatasets;
65 : double adfGeoTransform[6];
66 :
67 : char **papszFileList;
68 :
69 : public:
70 11 : ECRGTOCDataset()
71 11 : {
72 11 : papszSubDatasets = NULL;
73 11 : papszFileList = NULL;
74 11 : }
75 :
76 11 : ~ECRGTOCDataset()
77 11 : {
78 11 : CSLDestroy( papszSubDatasets );
79 11 : CSLDestroy(papszFileList);
80 11 : }
81 :
82 : virtual char **GetMetadata( const char * pszDomain = "" );
83 :
84 1 : virtual char **GetFileList() { return CSLDuplicate(papszFileList); }
85 :
86 : void AddSubDataset(const char* pszFilename,
87 : const char* pszProductTitle,
88 : const char* pszDiscId);
89 :
90 1 : virtual CPLErr GetGeoTransform( double * padfGeoTransform)
91 : {
92 1 : memcpy(padfGeoTransform, adfGeoTransform, 6 * sizeof(double));
93 1 : return CE_None;
94 : }
95 :
96 1 : virtual const char *GetProjectionRef(void)
97 : {
98 1 : return SRS_WKT_WGS84;
99 : }
100 :
101 : static GDALDataset* Build( const char* pszTOCFilename,
102 : CPLXMLNode* psXML,
103 : CPLString osProduct,
104 : CPLString osDiscId,
105 : const char* pszFilename);
106 :
107 : static int Identify( GDALOpenInfo * poOpenInfo );
108 : static GDALDataset* Open( GDALOpenInfo * poOpenInfo );
109 : };
110 :
111 : /************************************************************************/
112 : /* ==================================================================== */
113 : /* ECRGTOCSubDataset */
114 : /* ==================================================================== */
115 : /************************************************************************/
116 :
117 : class ECRGTOCSubDataset : public VRTDataset
118 : {
119 : char** papszFileList;
120 :
121 : public:
122 5 : ECRGTOCSubDataset(int nXSize, int nYSize) : VRTDataset(nXSize, nYSize)
123 : {
124 : /* Don't try to write a VRT file */
125 5 : SetWritable(FALSE);
126 :
127 : /* The driver is set to VRT in VRTDataset constructor. */
128 : /* We have to set it to the expected value ! */
129 5 : poDriver = (GDALDriver *) GDALGetDriverByName( "ECRGTOC" );
130 :
131 5 : papszFileList = NULL;
132 5 : }
133 :
134 5 : ~ECRGTOCSubDataset()
135 5 : {
136 5 : CSLDestroy(papszFileList);
137 5 : }
138 :
139 3 : virtual char **GetFileList() { return CSLDuplicate(papszFileList); }
140 :
141 : static GDALDataset* Build( const char* pszProductTitle,
142 : const char* pszDiscId,
143 : int nScale,
144 : int nCountSubDataset,
145 : const char* pszTOCFilename,
146 : const std::vector<FrameDesc>& aosFrameDesc,
147 : double dfGlobalMinX,
148 : double dfGlobalMinY,
149 : double dfGlobalMaxX,
150 : double dfGlobalMaxY,
151 : double dfGlobalPixelXSize,
152 : double dfGlobalPixelYSize);
153 : };
154 :
155 : /************************************************************************/
156 : /* AddSubDataset() */
157 : /************************************************************************/
158 :
159 6 : void ECRGTOCDataset::AddSubDataset( const char* pszFilename,
160 : const char* pszProductTitle,
161 : const char* pszDiscId )
162 :
163 : {
164 : char szName[80];
165 6 : int nCount = CSLCount(papszSubDatasets ) / 2;
166 :
167 6 : sprintf( szName, "SUBDATASET_%d_NAME", nCount+1 );
168 : papszSubDatasets =
169 : CSLSetNameValue( papszSubDatasets, szName,
170 : CPLSPrintf( "ECRG_TOC_ENTRY:%s:%s:%s",
171 6 : pszProductTitle, pszDiscId, pszFilename ) );
172 :
173 6 : sprintf( szName, "SUBDATASET_%d_DESC", nCount+1 );
174 : papszSubDatasets =
175 : CSLSetNameValue( papszSubDatasets, szName,
176 6 : CPLSPrintf( "%s:%s", pszProductTitle, pszDiscId));
177 6 : }
178 :
179 : /************************************************************************/
180 : /* GetMetadata() */
181 : /************************************************************************/
182 :
183 5 : char **ECRGTOCDataset::GetMetadata( const char *pszDomain )
184 :
185 : {
186 5 : if( pszDomain != NULL && EQUAL(pszDomain,"SUBDATASETS") )
187 5 : return papszSubDatasets;
188 :
189 0 : return GDALPamDataset::GetMetadata( pszDomain );
190 : }
191 :
192 : /************************************************************************/
193 : /* GetScaleFromString() */
194 : /************************************************************************/
195 :
196 11 : static int GetScaleFromString(const char* pszScale)
197 : {
198 11 : const char* pszPtr = strstr(pszScale, "1:");
199 11 : if (pszPtr)
200 11 : pszPtr = pszPtr + 2;
201 : else
202 0 : pszPtr = pszScale;
203 :
204 11 : int nScale = 0;
205 : char ch;
206 66 : while((ch = *pszPtr) != '\0')
207 : {
208 88 : if (ch >= '0' && ch <= '9')
209 33 : nScale = nScale * 10 + ch - '0';
210 22 : else if (ch == ' ')
211 : ;
212 11 : else if (ch == 'k' || ch == 'K')
213 11 : return nScale * 1000;
214 0 : else if (ch == 'm' || ch == 'M')
215 0 : return nScale * 1000000;
216 : else
217 0 : return 0;
218 44 : pszPtr ++;
219 : }
220 0 : return nScale;
221 : }
222 :
223 : /************************************************************************/
224 : /* GetFromBase34() */
225 : /************************************************************************/
226 :
227 42 : static GIntBig GetFromBase34(const char* pszVal, int nMaxSize)
228 : {
229 : int i;
230 42 : GIntBig nFrameNumber = 0;
231 462 : for(i=0;i<nMaxSize;i++)
232 : {
233 420 : char ch = pszVal[i];
234 420 : if (ch == '\0')
235 0 : break;
236 : int chVal;
237 420 : if (ch >= 'A' && ch <= 'Z')
238 0 : ch += 'a' - 'A';
239 : /* i and o letters are excluded, */
240 792 : if (ch >= '0' && ch <= '9')
241 372 : chVal = ch - '0';
242 57 : else if (ch >= 'a' && ch <= 'h')
243 9 : chVal = ch - 'a' + 10;
244 42 : else if (ch >= 'j' && ch < 'n')
245 3 : chVal = ch - 'a' + 10 - 1;
246 72 : else if (ch > 'p' && ch <= 'z')
247 36 : chVal = ch - 'a' + 10 - 2;
248 : else
249 : {
250 0 : CPLDebug("ECRG", "Invalid base34 value : %s", pszVal);
251 0 : break;
252 : }
253 420 : nFrameNumber = nFrameNumber * 34 + chVal;
254 : }
255 :
256 42 : return nFrameNumber;
257 : }
258 :
259 : /************************************************************************/
260 : /* GetExtent() */
261 : /************************************************************************/
262 :
263 : /* MIL-PRF-32283 - Table II. ECRG zone limits. */
264 : /* starting with a fake zone 0 for conveniency */
265 : static const int anZoneUpperLat[] = { 0, 32, 48, 56, 64, 68, 72, 76, 80 };
266 :
267 : /* APPENDIX 70, TABLE III of MIL-A-89007 */
268 : static const int anACst_ADRG[] =
269 : { 369664, 302592, 245760, 199168, 163328, 137216, 110080, 82432 };
270 : static const int nBCst_ADRG = 400384;
271 :
272 : #define CEIL_ROUND(a, b) (int)(ceil((double)(a)/(b))*(b))
273 : #define NEAR_ROUND(a, b) (int)(floor((double)(a)/(b) + 0.5)*(b))
274 :
275 : #define ECRG_PIXELS 2304
276 :
277 : static
278 42 : int GetExtent(const char* pszFrameName, int nScale, int nZone,
279 : double& dfMinX, double& dfMaxX, double& dfMinY, double& dfMaxY,
280 : double& dfPixelXSize, double& dfPixelYSize)
281 : {
282 42 : int nAbsZone = abs(nZone);
283 :
284 : /************************************************************************/
285 : /* Compute east-west constant */
286 : /************************************************************************/
287 : /* MIL-PRF-89038 - 60.1.2 - East-west pixel constant. */
288 42 : int nEW_ADRG = CEIL_ROUND(anACst_ADRG[nAbsZone-1] * (1e6 / nScale), 512);
289 42 : int nEW_CADRG = NEAR_ROUND(nEW_ADRG / (150. / 100.), 256);
290 : /* MIL-PRF-32283 - D.2.1.2 - East-west pixel constant. */
291 42 : int nEW = nEW_CADRG / 256 * 384;
292 :
293 : /************************************************************************/
294 : /* Compute number of longitudinal frames */
295 : /************************************************************************/
296 : /* MIL-PRF-32283 - D.2.1.7 - Longitudinal frames and subframes */
297 42 : int nCols = (int)ceil((double)nEW / ECRG_PIXELS);
298 :
299 : /************************************************************************/
300 : /* Compute north-south constant */
301 : /************************************************************************/
302 : /* MIL-PRF-89038 - 60.1.1 - North-south. pixel constant */
303 42 : int nNS_ADRG = CEIL_ROUND(nBCst_ADRG * (1e6 / nScale), 512) / 4;
304 42 : int nNS_CADRG = NEAR_ROUND(nNS_ADRG / (150. / 100.), 256);
305 : /* MIL-PRF-32283 - D.2.1.1 - North-south. pixel constant and Frame Width/Height */
306 42 : int nNS = nNS_CADRG / 256 * 384;
307 :
308 : /************************************************************************/
309 : /* Compute number of latitudinal frames and latitude of top of zone */
310 : /************************************************************************/
311 42 : dfPixelYSize = 90.0 / nNS;
312 :
313 42 : double dfFrameLatHeight = dfPixelYSize * ECRG_PIXELS;
314 :
315 : /* MIL-PRF-32283 - D.2.1.5 - Equatorward and poleward zone extents. */
316 42 : int nUpperZoneFrames = (int)ceil(anZoneUpperLat[nAbsZone] / dfFrameLatHeight);
317 42 : int nBottomZoneFrames = (int)floor(anZoneUpperLat[nAbsZone-1] / dfFrameLatHeight);
318 42 : int nRows = nUpperZoneFrames - nBottomZoneFrames;
319 :
320 : /* Not sure to really understand D.2.1.5.a. Testing needed */
321 42 : if (nZone < 0)
322 : {
323 0 : nUpperZoneFrames = -nBottomZoneFrames;
324 0 : nBottomZoneFrames = nUpperZoneFrames - nRows;
325 : }
326 :
327 42 : double dfUpperZoneTopLat = dfFrameLatHeight * nUpperZoneFrames;
328 :
329 : /************************************************************************/
330 : /* Compute coordinates of the frame in the zone */
331 : /************************************************************************/
332 :
333 : /* Converts the first 10 characters into a number from base 34 */
334 42 : GIntBig nFrameNumber = GetFromBase34(pszFrameName, 10);
335 :
336 : /* MIL-PRF-32283 - A.2.6.1 */
337 42 : GIntBig nY = nFrameNumber / nCols;
338 42 : GIntBig nX = nFrameNumber % nCols;
339 :
340 : /************************************************************************/
341 : /* Compute extent of the frame */
342 : /************************************************************************/
343 :
344 : /* The nY is counted from the bottom of the zone... Pfff */
345 42 : dfMaxY = dfUpperZoneTopLat - (nRows - 1 - nY) * dfFrameLatHeight;
346 42 : dfMinY = dfMaxY - dfFrameLatHeight;
347 :
348 42 : dfPixelXSize = 360.0 / nEW;
349 :
350 42 : double dfFrameLongWidth = dfPixelXSize * ECRG_PIXELS;
351 42 : dfMinX = -180.0 + nX * dfFrameLongWidth;
352 42 : dfMaxX = dfMinX + dfFrameLongWidth;
353 :
354 : //CPLDebug("ECRG", "Frame %s : minx=%.16g, maxy=%.16g, maxx=%.16g, miny=%.16g",
355 : // pszFrameName, dfMinX, dfMaxY, dfMaxX, dfMinY);
356 :
357 42 : return TRUE;
358 : }
359 :
360 : /************************************************************************/
361 : /* ==================================================================== */
362 : /* ECRGTOCProxyRasterDataSet */
363 : /* ==================================================================== */
364 : /************************************************************************/
365 :
366 : class ECRGTOCProxyRasterDataSet : public GDALProxyPoolDataset
367 14 : {
368 : /* The following parameters are only for sanity checking */
369 : int checkDone;
370 : int checkOK;
371 : double dfMinX;
372 : double dfMaxY;
373 : double dfPixelXSize;
374 : double dfPixelYSize;
375 : ECRGTOCSubDataset* poSubDataset;
376 :
377 : public:
378 : ECRGTOCProxyRasterDataSet(ECRGTOCSubDataset* poSubDataset,
379 : const char* fileName,
380 : int nXSize, int nYSize,
381 : double dfMinX, double dfMaxY,
382 : double dfPixelXSize, double dfPixelYSize);
383 :
384 4716 : GDALDataset* RefUnderlyingDataset()
385 : {
386 4716 : GDALDataset* poSourceDS = GDALProxyPoolDataset::RefUnderlyingDataset();
387 4716 : if (poSourceDS)
388 : {
389 4716 : if (!checkDone)
390 4 : SanityCheckOK(poSourceDS);
391 4716 : if (!checkOK)
392 : {
393 0 : GDALProxyPoolDataset::UnrefUnderlyingDataset(poSourceDS);
394 0 : poSourceDS = NULL;
395 : }
396 : }
397 4716 : return poSourceDS;
398 : }
399 :
400 4716 : void UnrefUnderlyingDataset(GDALDataset* poUnderlyingDataset)
401 : {
402 4716 : GDALProxyPoolDataset::UnrefUnderlyingDataset(poUnderlyingDataset);
403 4716 : }
404 :
405 : int SanityCheckOK(GDALDataset* poSourceDS);
406 : };
407 :
408 : /************************************************************************/
409 : /* ECRGTOCProxyRasterDataSet() */
410 : /************************************************************************/
411 :
412 14 : ECRGTOCProxyRasterDataSet::ECRGTOCProxyRasterDataSet
413 : (ECRGTOCSubDataset* poSubDataset,
414 : const char* fileName,
415 : int nXSize, int nYSize,
416 : double dfMinX, double dfMaxY,
417 : double dfPixelXSize, double dfPixelYSize) :
418 : /* Mark as shared since the VRT will take several references if we are in RGBA mode (4 bands for this dataset) */
419 14 : GDALProxyPoolDataset(fileName, nXSize, nYSize, GA_ReadOnly, TRUE, SRS_WKT_WGS84)
420 : {
421 : int i;
422 14 : this->poSubDataset = poSubDataset;
423 14 : this->dfMinX = dfMinX;
424 14 : this->dfMaxY = dfMaxY;
425 14 : this->dfPixelXSize = dfPixelXSize;
426 14 : this->dfPixelYSize = dfPixelYSize;
427 :
428 14 : checkDone = FALSE;
429 14 : checkOK = FALSE;
430 :
431 56 : for(i=0;i<3;i++)
432 : {
433 42 : SetBand(i + 1, new GDALProxyPoolRasterBand(this, i+1, GDT_Byte, nXSize, 1));
434 : }
435 14 : }
436 :
437 : /************************************************************************/
438 : /* SanityCheckOK() */
439 : /************************************************************************/
440 :
441 : #define WARN_CHECK_DS(x) do { if (!(x)) { CPLError(CE_Warning, CPLE_AppDefined,\
442 : "For %s, assert '" #x "' failed", GetDescription()); checkOK = FALSE; } } while(0)
443 :
444 4 : int ECRGTOCProxyRasterDataSet::SanityCheckOK(GDALDataset* poSourceDS)
445 : {
446 : /*int nSrcBlockXSize, nSrcBlockYSize;
447 : int nBlockXSize, nBlockYSize;*/
448 : double adfGeoTransform[6];
449 4 : if (checkDone)
450 0 : return checkOK;
451 :
452 4 : checkOK = TRUE;
453 4 : checkDone = TRUE;
454 :
455 4 : poSourceDS->GetGeoTransform(adfGeoTransform);
456 4 : WARN_CHECK_DS(fabs(adfGeoTransform[0] - dfMinX) < 1e-10);
457 4 : WARN_CHECK_DS(fabs(adfGeoTransform[3] - dfMaxY) < 1e-10);
458 4 : WARN_CHECK_DS(fabs(adfGeoTransform[1] - dfPixelXSize) < 1e-10);
459 4 : WARN_CHECK_DS(fabs(adfGeoTransform[5] - (-dfPixelYSize)) < 1e-10);
460 4 : WARN_CHECK_DS(adfGeoTransform[2] == 0 &&
461 : adfGeoTransform[4] == 0); /* No rotation */
462 4 : WARN_CHECK_DS(poSourceDS->GetRasterCount() == 3);
463 4 : WARN_CHECK_DS(poSourceDS->GetRasterXSize() == nRasterXSize);
464 4 : WARN_CHECK_DS(poSourceDS->GetRasterYSize() == nRasterYSize);
465 4 : WARN_CHECK_DS(EQUAL(poSourceDS->GetProjectionRef(), SRS_WKT_WGS84));
466 : /*poSourceDS->GetRasterBand(1)->GetBlockSize(&nSrcBlockXSize, &nSrcBlockYSize);
467 : GetRasterBand(1)->GetBlockSize(&nBlockXSize, &nBlockYSize);
468 : WARN_CHECK_DS(nSrcBlockXSize == nBlockXSize);
469 : WARN_CHECK_DS(nSrcBlockYSize == nBlockYSize);*/
470 4 : WARN_CHECK_DS(poSourceDS->GetRasterBand(1)->GetRasterDataType() == GDT_Byte);
471 :
472 4 : return checkOK;
473 : }
474 :
475 : /************************************************************************/
476 : /* BuildFullName() */
477 : /************************************************************************/
478 :
479 42 : static const char* BuildFullName(const char* pszTOCFilename,
480 : const char* pszFramePath,
481 : const char* pszFrameName)
482 : {
483 : char* pszPath;
484 96 : if (pszFramePath[0] == '.' &&
485 36 : (pszFramePath[1] == '/' ||pszFramePath[1] == '\\'))
486 18 : pszPath = CPLStrdup(pszFramePath + 2);
487 : else
488 24 : pszPath = CPLStrdup(pszFramePath);
489 366 : for(int i=0;pszPath[i] != '\0';i++)
490 : {
491 324 : if (pszPath[i] == '\\')
492 60 : pszPath[i] = '/';
493 : }
494 42 : const char* pszName = CPLFormFilename(pszPath, pszFrameName, NULL);
495 42 : CPLFree(pszPath);
496 42 : pszPath = NULL;
497 42 : const char* pszTOCPath = CPLGetDirname(pszTOCFilename);
498 42 : const char* pszFirstSlashInName = strchr(pszName, '/');
499 42 : if (pszFirstSlashInName != NULL)
500 : {
501 42 : int nFirstDirLen = pszFirstSlashInName - pszName;
502 108 : if ((int)strlen(pszTOCPath) >= nFirstDirLen + 1 &&
503 42 : (pszTOCPath[strlen(pszTOCPath) - (nFirstDirLen + 1)] == '/' ||
504 24 : pszTOCPath[strlen(pszTOCPath) - (nFirstDirLen + 1)] == '\\') &&
505 : strncmp(pszTOCPath + strlen(pszTOCPath) - nFirstDirLen, pszName, nFirstDirLen) == 0)
506 : {
507 18 : pszTOCPath = CPLGetDirname(pszTOCPath);
508 : }
509 : }
510 42 : return CPLProjectRelativeFilename(pszTOCPath, pszName);
511 : }
512 :
513 : /************************************************************************/
514 : /* Build() */
515 : /************************************************************************/
516 :
517 : /* Builds a ECRGTOCSubDataset from the set of files of the toc entry */
518 5 : GDALDataset* ECRGTOCSubDataset::Build( const char* pszProductTitle,
519 : const char* pszDiscId,
520 : int nScale,
521 : int nCountSubDataset,
522 : const char* pszTOCFilename,
523 : const std::vector<FrameDesc>& aosFrameDesc,
524 : double dfGlobalMinX,
525 : double dfGlobalMinY,
526 : double dfGlobalMaxX,
527 : double dfGlobalMaxY,
528 : double dfGlobalPixelXSize,
529 : double dfGlobalPixelYSize)
530 : {
531 : int i, j;
532 : GDALDriver *poDriver;
533 : ECRGTOCSubDataset *poVirtualDS;
534 : int nSizeX, nSizeY;
535 : double adfGeoTransform[6];
536 :
537 5 : poDriver = GetGDALDriverManager()->GetDriverByName("VRT");
538 5 : if( poDriver == NULL )
539 0 : return NULL;
540 :
541 5 : nSizeX = (int)((dfGlobalMaxX - dfGlobalMinX) / dfGlobalPixelXSize + 0.5);
542 5 : nSizeY = (int)((dfGlobalMaxY - dfGlobalMinY) / dfGlobalPixelYSize + 0.5);
543 :
544 : /* ------------------------------------ */
545 : /* Create the VRT with the overall size */
546 : /* ------------------------------------ */
547 5 : poVirtualDS = new ECRGTOCSubDataset( nSizeX, nSizeY );
548 :
549 5 : poVirtualDS->SetProjection(SRS_WKT_WGS84);
550 :
551 5 : adfGeoTransform[0] = dfGlobalMinX;
552 5 : adfGeoTransform[1] = dfGlobalPixelXSize;
553 5 : adfGeoTransform[2] = 0;
554 5 : adfGeoTransform[3] = dfGlobalMaxY;
555 5 : adfGeoTransform[4] = 0;
556 5 : adfGeoTransform[5] = -dfGlobalPixelYSize;
557 5 : poVirtualDS->SetGeoTransform(adfGeoTransform);
558 :
559 20 : for (i=0;i<3;i++)
560 : {
561 15 : poVirtualDS->AddBand(GDT_Byte, NULL);
562 15 : GDALRasterBand *poBand = poVirtualDS->GetRasterBand( i + 1 );
563 15 : poBand->SetColorInterpretation((GDALColorInterp)(GCI_RedBand+i));
564 : }
565 :
566 5 : poVirtualDS->SetDescription(pszTOCFilename);
567 :
568 5 : poVirtualDS->SetMetadataItem("PRODUCT_TITLE", pszProductTitle);
569 5 : poVirtualDS->SetMetadataItem("DISC_ID", pszDiscId);
570 5 : if (nScale != -1)
571 5 : poVirtualDS->SetMetadataItem("SCALE", CPLString().Printf("%d", nScale));
572 :
573 : /* -------------------------------------------------------------------- */
574 : /* Check for overviews. */
575 : /* -------------------------------------------------------------------- */
576 :
577 : poVirtualDS->oOvManager.Initialize( poVirtualDS,
578 5 : CPLString().Printf("%s.%d", pszTOCFilename, nCountSubDataset));
579 :
580 5 : poVirtualDS->papszFileList = poVirtualDS->GDALDataset::GetFileList();
581 :
582 19 : for(i=0;i<(int)aosFrameDesc.size(); i++)
583 : {
584 : const char* pszName = BuildFullName(pszTOCFilename,
585 : aosFrameDesc[i].pszPath,
586 14 : aosFrameDesc[i].pszName);
587 :
588 14 : double dfMinX = 0, dfMaxX = 0, dfMinY = 0, dfMaxY = 0,
589 14 : dfPixelXSize = 0, dfPixelYSize = 0;
590 : GetExtent(aosFrameDesc[i].pszName,
591 : aosFrameDesc[i].nScale, aosFrameDesc[i].nZone,
592 14 : dfMinX, dfMaxX, dfMinY, dfMaxY, dfPixelXSize, dfPixelYSize);
593 :
594 14 : int nFrameXSize = (int)((dfMaxX - dfMinX) / dfPixelXSize + 0.5);
595 14 : int nFrameYSize = (int)((dfMaxY - dfMinY) / dfPixelYSize + 0.5);
596 :
597 14 : poVirtualDS->papszFileList = CSLAddString(poVirtualDS->papszFileList, pszName);
598 :
599 : /* We create proxy datasets and raster bands */
600 : /* Using real datasets and raster bands is possible in theory */
601 : /* However for large datasets, a TOC entry can include several hundreds of files */
602 : /* and we finally reach the limit of maximum file descriptors open at the same time ! */
603 : /* So the idea is to warp the datasets into a proxy and open the underlying dataset only when it is */
604 : /* needed (IRasterIO operation). To improve a bit efficiency, we have a cache of opened */
605 : /* underlying datasets */
606 : ECRGTOCProxyRasterDataSet* poDS = new ECRGTOCProxyRasterDataSet(
607 : (ECRGTOCSubDataset*)poVirtualDS, pszName, nFrameXSize, nFrameYSize,
608 14 : dfMinX, dfMaxY, dfPixelXSize, dfPixelYSize);
609 :
610 56 : for(j=0;j<3;j++)
611 : {
612 : VRTSourcedRasterBand *poBand = (VRTSourcedRasterBand*)
613 42 : poVirtualDS->GetRasterBand( j + 1 );
614 : /* Place the raster band at the right position in the VRT */
615 : poBand->AddSimpleSource(poDS->GetRasterBand(j + 1),
616 : 0, 0, nFrameXSize, nFrameYSize,
617 : (int)((dfMinX - dfGlobalMinX) / dfGlobalPixelXSize + 0.5),
618 : (int)((dfGlobalMaxY - dfMaxY) / dfGlobalPixelYSize + 0.5),
619 : (int)((dfMaxX - dfMinX) / dfGlobalPixelXSize + 0.5),
620 42 : (int)((dfMaxY - dfMinY) / dfGlobalPixelYSize + 0.5));
621 : }
622 :
623 : /* The ECRGTOCProxyRasterDataSet will be destroyed when its last raster band will be */
624 : /* destroyed */
625 14 : poDS->Dereference();
626 : }
627 :
628 5 : poVirtualDS->SetMetadataItem("INTERLEAVE", "PIXEL", "IMAGE_STRUCTURE");
629 :
630 5 : return poVirtualDS;
631 : }
632 :
633 : /************************************************************************/
634 : /* Build() */
635 : /************************************************************************/
636 :
637 11 : GDALDataset* ECRGTOCDataset::Build(const char* pszTOCFilename,
638 : CPLXMLNode* psXML,
639 : CPLString osProduct,
640 : CPLString osDiscId,
641 : const char* pszOpenInfoFilename)
642 : {
643 11 : CPLXMLNode* psTOC = CPLGetXMLNode(psXML, "=Table_of_Contents");
644 11 : if (psTOC == NULL)
645 : {
646 : CPLError(CE_Warning, CPLE_AppDefined,
647 0 : "Cannot find Table_of_Contents element");
648 0 : return NULL;
649 : }
650 :
651 11 : double dfGlobalMinX = 0, dfGlobalMinY = 0, dfGlobalMaxX = 0, dfGlobalMaxY= 0;
652 11 : double dfGlobalPixelXSize = 0, dfGlobalPixelYSize = 0;
653 11 : int bGlobalExtentValid = FALSE;
654 :
655 11 : ECRGTOCDataset* poDS = new ECRGTOCDataset();
656 11 : int nSubDatasets = 0;
657 :
658 17 : int bLookForSubDataset = osProduct.size() != 0 && osDiscId.size() != 0;
659 :
660 11 : int nCountSubDataset = 0;
661 :
662 11 : poDS->SetDescription( pszOpenInfoFilename );
663 11 : poDS->papszFileList = poDS->GDALDataset::GetFileList();
664 :
665 34 : for(CPLXMLNode* psIter1 = psTOC->psChild;
666 : psIter1 != NULL;
667 : psIter1 = psIter1->psNext)
668 : {
669 28 : if (!(psIter1->eType == CXT_Element && psIter1->pszValue != NULL &&
670 : strcmp(psIter1->pszValue, "product") == 0))
671 17 : continue;
672 :
673 : const char* pszProductTitle =
674 11 : CPLGetXMLValue(psIter1, "product_title", NULL);
675 11 : if (pszProductTitle == NULL)
676 : {
677 : CPLError(CE_Warning, CPLE_AppDefined,
678 0 : "Cannot find product_title attribute");
679 0 : continue;
680 : }
681 :
682 11 : if (bLookForSubDataset && strcmp(pszProductTitle, osProduct.c_str()) != 0)
683 0 : continue;
684 :
685 11 : for(CPLXMLNode* psIter2 = psIter1->psChild;
686 : psIter2 != NULL;
687 : psIter2 = psIter2->psNext)
688 : {
689 23 : if (!(psIter2->eType == CXT_Element && psIter2->pszValue != NULL &&
690 : strcmp(psIter2->pszValue, "disc") == 0))
691 11 : continue;
692 :
693 12 : const char* pszDiscId = CPLGetXMLValue(psIter2, "id", NULL);
694 12 : if (pszDiscId == NULL)
695 : {
696 : CPLError(CE_Warning, CPLE_AppDefined,
697 0 : "Cannot find id attribute");
698 0 : continue;
699 : }
700 :
701 12 : if (bLookForSubDataset && strcmp(pszDiscId, osDiscId.c_str()) != 0)
702 1 : continue;
703 :
704 11 : nCountSubDataset ++;
705 :
706 11 : CPLXMLNode* psFrameList = CPLGetXMLNode(psIter2, "frame_list");
707 11 : if (psFrameList == NULL)
708 : {
709 : CPLError(CE_Warning, CPLE_AppDefined,
710 0 : "Cannot find frame_list element");
711 0 : continue;
712 : }
713 :
714 11 : int nValidFrames = 0;
715 :
716 11 : std::vector<FrameDesc> aosFrameDesc;
717 :
718 11 : int nSubDatasetScale = -1;
719 :
720 33 : for(CPLXMLNode* psIter3 = psFrameList->psChild;
721 : psIter3 != NULL;
722 : psIter3 = psIter3->psNext)
723 : {
724 22 : if (!(psIter3->eType == CXT_Element &&
725 : psIter3->pszValue != NULL &&
726 : strcmp(psIter3->pszValue, "scale") == 0))
727 11 : continue;
728 :
729 11 : const char* pszSize = CPLGetXMLValue(psIter3, "size", NULL);
730 11 : if (pszSize == NULL)
731 : {
732 : CPLError(CE_Warning, CPLE_AppDefined,
733 0 : "Cannot find size attribute");
734 0 : continue;
735 : }
736 :
737 11 : int nScale = GetScaleFromString(pszSize);
738 11 : if (nScale <= 0)
739 : {
740 : CPLError(CE_Warning, CPLE_AppDefined,
741 0 : "Invalid scale %s", pszSize);
742 0 : continue;
743 : }
744 :
745 11 : if (nValidFrames == 0)
746 11 : nSubDatasetScale = nScale;
747 : else
748 0 : nSubDatasetScale = -1;
749 :
750 50 : for(CPLXMLNode* psIter4 = psIter3->psChild;
751 : psIter4 != NULL;
752 : psIter4 = psIter4->psNext)
753 : {
754 39 : if (!(psIter4->eType == CXT_Element &&
755 : psIter4->pszValue != NULL &&
756 : strcmp(psIter4->pszValue, "frame") == 0))
757 11 : continue;
758 :
759 : const char* pszFrameName =
760 28 : CPLGetXMLValue(psIter4, "name", NULL);
761 28 : if (pszFrameName == NULL)
762 : {
763 : CPLError(CE_Warning, CPLE_AppDefined,
764 0 : "Cannot find name element");
765 0 : continue;
766 : }
767 :
768 28 : if (strlen(pszFrameName) != 18)
769 : {
770 : CPLError(CE_Warning, CPLE_AppDefined,
771 : "Invalid value for name element : %s",
772 0 : pszFrameName);
773 0 : continue;
774 : }
775 :
776 : const char* pszFramePath =
777 28 : CPLGetXMLValue(psIter4, "frame_path", NULL);
778 28 : if (pszFramePath == NULL)
779 : {
780 : CPLError(CE_Warning, CPLE_AppDefined,
781 0 : "Cannot find frame_path element");
782 0 : continue;
783 : }
784 :
785 : const char* pszFrameZone =
786 28 : CPLGetXMLValue(psIter4, "frame_zone", NULL);
787 28 : if (pszFrameZone == NULL)
788 : {
789 : CPLError(CE_Warning, CPLE_AppDefined,
790 0 : "Cannot find frame_zone element");
791 0 : continue;
792 : }
793 28 : if (strlen(pszFrameZone) != 1)
794 : {
795 : CPLError(CE_Warning, CPLE_AppDefined,
796 : "Invalid value for frame_zone element : %s",
797 0 : pszFrameZone);
798 0 : continue;
799 : }
800 28 : char chZone = pszFrameZone[0];
801 28 : int nZone = 0;
802 56 : if (chZone >= '1' && chZone <= '9')
803 28 : nZone = chZone - '0';
804 0 : else if (chZone >= 'a' && chZone <= 'h')
805 0 : nZone = -(chZone - 'a' + 1);
806 0 : else if (chZone >= 'A' && chZone <= 'H')
807 0 : nZone = -(chZone - 'A' + 1);
808 0 : else if (chZone == 'j' || chZone == 'J')
809 0 : nZone = -9;
810 : else
811 : {
812 : CPLError(CE_Warning, CPLE_AppDefined,
813 : "Invalid value for frame_zone element : %s",
814 0 : pszFrameZone);
815 0 : continue;
816 : }
817 28 : if (nZone == 9 || nZone == -9)
818 : {
819 : CPLError(CE_Warning, CPLE_AppDefined,
820 0 : "Polar zones unhandled by current implementation");
821 0 : continue;
822 : }
823 :
824 28 : double dfMinX = 0, dfMaxX = 0,
825 28 : dfMinY = 0, dfMaxY = 0,
826 28 : dfPixelXSize = 0, dfPixelYSize = 0;
827 28 : if (!GetExtent(pszFrameName, nScale, nZone,
828 : dfMinX, dfMaxX, dfMinY, dfMaxY,
829 : dfPixelXSize, dfPixelYSize))
830 : {
831 0 : continue;
832 : }
833 :
834 28 : nValidFrames ++;
835 :
836 : const char* pszFullName = BuildFullName(pszTOCFilename,
837 : pszFramePath,
838 28 : pszFrameName);
839 28 : poDS->papszFileList = CSLAddString(poDS->papszFileList, pszFullName);
840 :
841 28 : if (!bGlobalExtentValid)
842 : {
843 10 : dfGlobalMinX = dfMinX;
844 10 : dfGlobalMinY = dfMinY;
845 10 : dfGlobalMaxX = dfMaxX;
846 10 : dfGlobalMaxY = dfMaxY;
847 10 : dfGlobalPixelXSize = dfPixelXSize;
848 10 : dfGlobalPixelYSize = dfPixelYSize;
849 10 : bGlobalExtentValid = TRUE;
850 : }
851 : else
852 : {
853 18 : if (dfMinX < dfGlobalMinX) dfGlobalMinX = dfMinX;
854 18 : if (dfMinY < dfGlobalMinY) dfGlobalMinY = dfMinY;
855 18 : if (dfMaxX > dfGlobalMaxX) dfGlobalMaxX = dfMaxX;
856 18 : if (dfMaxY > dfGlobalMaxY) dfGlobalMaxY = dfMaxY;
857 18 : if (dfPixelXSize < dfGlobalPixelXSize)
858 0 : dfGlobalPixelXSize = dfPixelXSize;
859 18 : if (dfPixelYSize < dfGlobalPixelYSize)
860 0 : dfGlobalPixelYSize = dfPixelYSize;
861 : }
862 :
863 28 : if (bLookForSubDataset)
864 : {
865 : FrameDesc frameDesc;
866 14 : frameDesc.pszName = pszFrameName;
867 14 : frameDesc.pszPath = pszFramePath;
868 14 : frameDesc.nScale = nScale;
869 14 : frameDesc.nZone = nZone;
870 14 : aosFrameDesc.push_back(frameDesc);
871 : }
872 : }
873 : }
874 :
875 11 : if (bLookForSubDataset)
876 : {
877 5 : delete poDS;
878 5 : if (nValidFrames == 0)
879 0 : return NULL;
880 : return ECRGTOCSubDataset::Build(pszProductTitle,
881 : pszDiscId,
882 : nSubDatasetScale,
883 : nCountSubDataset,
884 : pszTOCFilename,
885 : aosFrameDesc,
886 : dfGlobalMinX,
887 : dfGlobalMinY,
888 : dfGlobalMaxX,
889 : dfGlobalMaxY,
890 : dfGlobalPixelXSize,
891 5 : dfGlobalPixelYSize);
892 : }
893 :
894 6 : if (nValidFrames)
895 : {
896 : poDS->AddSubDataset(pszOpenInfoFilename,
897 6 : pszProductTitle, pszDiscId);
898 6 : nSubDatasets ++;
899 : }
900 : }
901 : }
902 :
903 6 : if (!bGlobalExtentValid)
904 : {
905 1 : delete poDS;
906 1 : return NULL;
907 : }
908 :
909 5 : if (nSubDatasets == 1)
910 : {
911 : const char* pszSubDatasetName = CSLFetchNameValue(
912 4 : poDS->GetMetadata("SUBDATASETS"), "SUBDATASET_1_NAME");
913 4 : GDALOpenInfo oOpenInfo(pszSubDatasetName, GA_ReadOnly);
914 4 : delete poDS;
915 4 : GDALDataset* poRetDS = Open(&oOpenInfo);
916 4 : if (poRetDS)
917 4 : poRetDS->SetDescription(pszOpenInfoFilename);
918 4 : return poRetDS;
919 : }
920 :
921 1 : poDS->adfGeoTransform[0] = dfGlobalMinX;
922 1 : poDS->adfGeoTransform[1] = dfGlobalPixelXSize;
923 1 : poDS->adfGeoTransform[2] = 0.0;
924 1 : poDS->adfGeoTransform[3] = dfGlobalMaxY;
925 1 : poDS->adfGeoTransform[4] = 0.0;
926 1 : poDS->adfGeoTransform[5] = - dfGlobalPixelYSize;
927 :
928 1 : poDS->nRasterXSize = (int)(0.5 + (dfGlobalMaxX - dfGlobalMinX) / dfGlobalPixelXSize);
929 1 : poDS->nRasterYSize = (int)(0.5 + (dfGlobalMaxY - dfGlobalMinY) / dfGlobalPixelYSize);
930 :
931 : /* -------------------------------------------------------------------- */
932 : /* Initialize any PAM information. */
933 : /* -------------------------------------------------------------------- */
934 1 : poDS->TryLoadXML();
935 :
936 1 : return poDS;
937 : }
938 :
939 : /************************************************************************/
940 : /* Identify() */
941 : /************************************************************************/
942 :
943 13290 : int ECRGTOCDataset::Identify( GDALOpenInfo * poOpenInfo )
944 :
945 : {
946 13290 : const char *pszFilename = poOpenInfo->pszFilename;
947 13290 : const char *pabyHeader = (const char *) poOpenInfo->pabyHeader;
948 :
949 : /* -------------------------------------------------------------------- */
950 : /* Is this a sub-dataset selector? If so, it is obviously ECRGTOC. */
951 : /* -------------------------------------------------------------------- */
952 13290 : if( EQUALN(pszFilename, "ECRG_TOC_ENTRY:",strlen("ECRG_TOC_ENTRY:")))
953 6 : return TRUE;
954 :
955 : /* -------------------------------------------------------------------- */
956 : /* First we check to see if the file has the expected header */
957 : /* bytes. */
958 : /* -------------------------------------------------------------------- */
959 13284 : if( pabyHeader == NULL )
960 10751 : return FALSE;
961 :
962 2533 : if ( strstr(pabyHeader, "<Table_of_Contents>") != NULL &&
963 : strstr(pabyHeader, "<file_header ") != NULL)
964 4 : return TRUE;
965 :
966 2529 : if ( strstr(pabyHeader, "<!DOCTYPE Table_of_Contents [") != NULL)
967 1 : return TRUE;
968 :
969 2528 : return FALSE;
970 : }
971 :
972 : /************************************************************************/
973 : /* Open() */
974 : /************************************************************************/
975 :
976 4010 : GDALDataset *ECRGTOCDataset::Open( GDALOpenInfo * poOpenInfo )
977 :
978 : {
979 4010 : const char *pszFilename = poOpenInfo->pszFilename;
980 4010 : CPLString osProduct, osDiscId;
981 :
982 4010 : if( !Identify( poOpenInfo ) )
983 3999 : return NULL;
984 :
985 11 : if( EQUALN(pszFilename, "ECRG_TOC_ENTRY:",strlen("ECRG_TOC_ENTRY:")))
986 : {
987 6 : pszFilename += strlen("ECRG_TOC_ENTRY:");
988 6 : osProduct = pszFilename;
989 6 : size_t iPos = osProduct.find(":");
990 6 : if (iPos == std::string::npos)
991 0 : return NULL;
992 6 : osProduct.resize(iPos);
993 :
994 6 : pszFilename += iPos + 1;
995 6 : osDiscId = pszFilename;
996 6 : iPos = osDiscId.find(":");
997 6 : if (iPos == std::string::npos)
998 0 : return NULL;
999 6 : osDiscId.resize(iPos);
1000 :
1001 6 : pszFilename += iPos + 1;
1002 : }
1003 :
1004 : /* -------------------------------------------------------------------- */
1005 : /* Parse the XML file */
1006 : /* -------------------------------------------------------------------- */
1007 11 : CPLXMLNode* psXML = CPLParseXMLFile(pszFilename);
1008 11 : if (psXML == NULL)
1009 : {
1010 0 : return NULL;
1011 : }
1012 :
1013 : GDALDataset* poDS = Build( pszFilename, psXML, osProduct, osDiscId,
1014 11 : poOpenInfo->pszFilename);
1015 11 : CPLDestroyXMLNode(psXML);
1016 :
1017 11 : if (poDS && poOpenInfo->eAccess == GA_Update)
1018 : {
1019 : CPLError(CE_Failure, CPLE_NotSupported,
1020 0 : "ECRGTOC driver does not support update mode");
1021 0 : delete poDS;
1022 0 : return NULL;
1023 : }
1024 :
1025 11 : return poDS;
1026 : }
1027 :
1028 : /************************************************************************/
1029 : /* GDALRegister_ECRGTOC() */
1030 : /************************************************************************/
1031 :
1032 558 : void GDALRegister_ECRGTOC()
1033 :
1034 : {
1035 : GDALDriver *poDriver;
1036 :
1037 558 : if( GDALGetDriverByName( "ECRGTOC" ) == NULL )
1038 : {
1039 537 : poDriver = new GDALDriver();
1040 :
1041 537 : poDriver->SetDescription( "ECRGTOC" );
1042 : poDriver->SetMetadataItem( GDAL_DMD_LONGNAME,
1043 537 : "ECRG TOC format" );
1044 :
1045 537 : poDriver->pfnIdentify = ECRGTOCDataset::Identify;
1046 537 : poDriver->pfnOpen = ECRGTOCDataset::Open;
1047 :
1048 : poDriver->SetMetadataItem( GDAL_DMD_HELPTOPIC,
1049 537 : "frmt_various.html#ECRGTOC" );
1050 537 : poDriver->SetMetadataItem( GDAL_DMD_EXTENSION, "xml" );
1051 537 : poDriver->SetMetadataItem( GDAL_DCAP_VIRTUALIO, "YES" );
1052 :
1053 537 : GetGDALDriverManager()->RegisterDriver( poDriver );
1054 : }
1055 558 : }
|