Line data Source code
1 : /******************************************************************************
2 : *
3 : * Project: DRDC Configurable Airborne SAR Processor (COASP) data reader
4 : * Purpose: Support in GDAL for the DRDC COASP format data, both Metadata
5 : * and complex imagery.
6 : * Author: Philippe Vachon <philippe@cowpig.ca>
7 : * Notes: I have seen a grand total of 2 COASP scenes (3 sets of headers).
8 : * This is based on my best observations, some educated guesses and
9 : * such. So if you have a scene that doesn't work, send it to me
10 : * please and I will make it work... with violence.
11 : *
12 : ******************************************************************************
13 : * Copyright (c) 2007, Philippe Vachon
14 : * Copyright (c) 2009-2012, Even Rouault <even dot rouault at spatialys.com>
15 : *
16 : * SPDX-License-Identifier: MIT
17 : ****************************************************************************/
18 :
19 : #include "cpl_conv.h"
20 : #include "cpl_port.h"
21 : #include "cpl_string.h"
22 : #include "cpl_vsi.h"
23 : #include "gdal_frmts.h"
24 : #include "gdal_priv.h"
25 :
26 : constexpr int TYPE_GENERIC = 0;
27 : constexpr int TYPE_GEOREF = 1;
28 :
29 : enum ePolarization
30 : {
31 : hh = 0,
32 : hv,
33 : vh,
34 : vv
35 : };
36 :
37 : /*******************************************************************
38 : * Declaration of the COASPMetadata classes *
39 : *******************************************************************/
40 :
41 : class COASPMetadataItem;
42 :
43 : class COASPMetadataReader
44 : {
45 : char **papszMetadata;
46 : int nMetadataCount;
47 : int nCurrentItem;
48 :
49 : public:
50 : explicit COASPMetadataReader(char *pszFname);
51 : ~COASPMetadataReader();
52 : COASPMetadataItem *GetNextItem();
53 : COASPMetadataItem *GetItem(int nItem);
54 : int GotoMetadataItem(const char *pszName);
55 :
56 : int GetCurrentItem() const
57 : {
58 : return nCurrentItem;
59 : }
60 : };
61 :
62 : /* Your average metadata item */
63 : class COASPMetadataItem
64 : {
65 : protected:
66 : char *pszItemName;
67 : char *pszItemValue;
68 :
69 : public:
70 0 : COASPMetadataItem() : pszItemName(nullptr), pszItemValue(nullptr)
71 : {
72 0 : }
73 :
74 : COASPMetadataItem(char *pszItemName, char *pszItemValue);
75 : ~COASPMetadataItem();
76 :
77 : char *GetItemName();
78 : char *GetItemValue();
79 :
80 : static int GetType()
81 : {
82 : return TYPE_GENERIC;
83 : }
84 : };
85 :
86 : /* Same as MetadataItem class except parses GCP properly and returns
87 : * a GDAL_GCP struct
88 : */
89 : class COASPMetadataGeorefGridItem : public COASPMetadataItem
90 : {
91 : #ifdef unused
92 : int nId;
93 : int nPixels;
94 : int nLines;
95 : double ndLat;
96 : double ndLong;
97 : #endif
98 :
99 : public:
100 : COASPMetadataGeorefGridItem(int nId, int nPixels, int nLines, double ndLat,
101 : double ndLong);
102 :
103 : static const char *GetItemName()
104 : {
105 : return "georef_grid";
106 : }
107 :
108 : static GDAL_GCP *GetItemValue();
109 :
110 : static int GetType()
111 : {
112 : return TYPE_GEOREF;
113 : }
114 : };
115 :
116 : /********************************************************************
117 : * ================================================================ *
118 : * Implementation of the COASPMetadataItem Classes *
119 : * ================================================================ *
120 : ********************************************************************/
121 :
122 0 : COASPMetadataItem::COASPMetadataItem(char *pszItemName_, char *pszItemValue_)
123 0 : : pszItemName(VSIStrdup(pszItemName_)),
124 0 : pszItemValue(VSIStrdup(pszItemValue_))
125 : {
126 0 : }
127 :
128 0 : COASPMetadataItem::~COASPMetadataItem()
129 : {
130 0 : CPLFree(pszItemName);
131 0 : CPLFree(pszItemValue);
132 0 : }
133 :
134 0 : char *COASPMetadataItem::GetItemName()
135 : {
136 0 : return VSIStrdup(pszItemName);
137 : }
138 :
139 0 : char *COASPMetadataItem::GetItemValue()
140 : {
141 0 : return VSIStrdup(pszItemValue);
142 : }
143 :
144 0 : COASPMetadataGeorefGridItem::COASPMetadataGeorefGridItem(int /*nIdIn*/,
145 : int /*nPixelsIn*/,
146 : int /*nLinesIn*/,
147 : double /*ndLatIn*/,
148 0 : double /*ndLongIn*/)
149 : #ifdef unused
150 : : nId(nIdIn), nPixels(nPixelsIn), nLines(nLinesIn), ndLat(ndLatIn),
151 : ndLong(ndLongIn)
152 : #endif
153 : {
154 0 : pszItemName = VSIStrdup("georef_grid");
155 0 : }
156 :
157 0 : GDAL_GCP *COASPMetadataGeorefGridItem::GetItemValue()
158 : {
159 0 : return nullptr;
160 : }
161 :
162 : /********************************************************************
163 : * ================================================================ *
164 : * Implementation of the COASPMetadataReader Class *
165 : * ================================================================ *
166 : ********************************************************************/
167 :
168 0 : COASPMetadataReader::COASPMetadataReader(char *pszFname)
169 0 : : papszMetadata(CSLLoad(pszFname)), nMetadataCount(0), nCurrentItem(0)
170 : {
171 0 : nMetadataCount = CSLCount(papszMetadata);
172 0 : }
173 :
174 0 : COASPMetadataReader::~COASPMetadataReader()
175 : {
176 0 : CSLDestroy(papszMetadata);
177 0 : }
178 :
179 0 : COASPMetadataItem *COASPMetadataReader::GetNextItem()
180 : {
181 0 : if (nCurrentItem < 0 || nCurrentItem >= nMetadataCount)
182 0 : return nullptr;
183 :
184 0 : COASPMetadataItem *poMetadata = nullptr;
185 :
186 0 : char **papszMDTokens = CSLTokenizeString2(papszMetadata[nCurrentItem], " ",
187 : CSLT_HONOURSTRINGS);
188 0 : char *pszItemName = papszMDTokens[0];
189 0 : if (STARTS_WITH_CI(pszItemName, "georef_grid") &&
190 0 : CSLCount(papszMDTokens) >= 8)
191 : {
192 : // georef_grid ( pixels lines ) ( lat long )
193 : // 0 1 2 3 4 5 6 7 8
194 0 : int nPixels = atoi(papszMDTokens[2]);
195 0 : int nLines = atoi(papszMDTokens[3]);
196 0 : double dfLat = CPLAtof(papszMDTokens[6]);
197 0 : double dfLong = CPLAtof(papszMDTokens[7]);
198 0 : poMetadata = new COASPMetadataGeorefGridItem(nCurrentItem, nPixels,
199 0 : nLines, dfLat, dfLong);
200 : }
201 : else
202 : {
203 0 : int nCount = CSLCount(papszMDTokens);
204 0 : if (nCount >= 2)
205 : {
206 0 : char *pszItemValue = CPLStrdup(papszMDTokens[1]);
207 0 : for (int i = 2; i < nCount; i++)
208 : {
209 0 : const size_t nSize =
210 0 : strlen(pszItemValue) + 1 + strlen(papszMDTokens[i]);
211 0 : pszItemValue = (char *)CPLRealloc(pszItemValue, nSize);
212 0 : snprintf(pszItemValue + strlen(pszItemValue),
213 0 : nSize - strlen(pszItemValue), " %s", papszMDTokens[i]);
214 : }
215 :
216 0 : poMetadata = new COASPMetadataItem(pszItemName, pszItemValue);
217 :
218 0 : CPLFree(pszItemValue);
219 : }
220 : }
221 0 : CSLDestroy(papszMDTokens);
222 0 : nCurrentItem++;
223 0 : return poMetadata;
224 : }
225 :
226 : /* Goto the first metadata item with a particular name */
227 0 : int COASPMetadataReader::GotoMetadataItem(const char *pszName)
228 : {
229 0 : nCurrentItem = CSLPartialFindString(papszMetadata, pszName);
230 0 : return nCurrentItem;
231 : }
232 :
233 : /*******************************************************************
234 : * Declaration of the COASPDataset class *
235 : *******************************************************************/
236 :
237 : class COASPRasterBand;
238 :
239 : /* A couple of observations based on the data I have available to me:
240 : * a) the headers don't really change, beyond indicating data sources
241 : * and such. As such, I only read the first header specified by the
242 : * user. Note that this is agnostic: you can specify hh, vv, vh, hv and
243 : * all the data needed will be immediately available.
244 : * b) Lots of GCPs are present in the headers. This is most excellent.
245 : * c) There is no documentation on this format. All the knowledge contained
246 : * herein is from harassing various Defence Scientists at DRDC Ottawa.
247 : */
248 :
249 : class COASPDataset final : public GDALDataset
250 : {
251 : friend class COASPRasterBand;
252 : VSILFILE *fpHdr; /* File pointer for the header file */
253 : VSILFILE *fpBinHH; /* File pointer for the binary matrix */
254 : VSILFILE *fpBinHV;
255 : VSILFILE *fpBinVH;
256 : VSILFILE *fpBinVV;
257 :
258 : char *pszFileName; /* line and mission ID, mostly, i.e. l27p7 */
259 :
260 : public:
261 0 : COASPDataset()
262 0 : : fpHdr(nullptr), fpBinHH(nullptr), fpBinHV(nullptr), fpBinVH(nullptr),
263 0 : fpBinVV(nullptr), pszFileName(nullptr)
264 : {
265 0 : }
266 :
267 : ~COASPDataset();
268 :
269 : static GDALDataset *Open(GDALOpenInfo *);
270 : static int Identify(GDALOpenInfo *poOpenInfo);
271 : };
272 :
273 : /********************************************************************
274 : * ================================================================ *
275 : * Declaration and implementation of the COASPRasterBand Class *
276 : * ================================================================ *
277 : ********************************************************************/
278 :
279 : class COASPRasterBand final : public GDALRasterBand
280 : {
281 : VSILFILE *fp;
282 : // int ePol;
283 : public:
284 : COASPRasterBand(COASPDataset *poDS, GDALDataType eDataType, int ePol,
285 : VSILFILE *fp);
286 : CPLErr IReadBlock(int nBlockXOff, int nBlockYOff, void *pImage) override;
287 : };
288 :
289 0 : COASPRasterBand::COASPRasterBand(COASPDataset *poDSIn, GDALDataType eDataTypeIn,
290 0 : int /*ePolIn*/, VSILFILE *fpIn)
291 0 : : fp(fpIn) /*,
292 : ePol(ePolIn)*/
293 : {
294 0 : poDS = poDSIn;
295 0 : eDataType = eDataTypeIn;
296 0 : nBlockXSize = poDS->GetRasterXSize();
297 0 : nBlockYSize = 1;
298 0 : }
299 :
300 0 : CPLErr COASPRasterBand::IReadBlock(CPL_UNUSED int nBlockXOff, int nBlockYOff,
301 : void *pImage)
302 : {
303 0 : if (this->fp == nullptr)
304 : {
305 0 : CPLError(CE_Fatal, CPLE_AppDefined, "File pointer freed unexpectedly");
306 0 : return CE_Fatal;
307 : }
308 :
309 : /* 8 bytes per pixel: 4 bytes I, 4 bytes Q */
310 : const vsi_l_offset nByteNum =
311 0 : static_cast<vsi_l_offset>(poDS->GetRasterXSize()) * 8 * nBlockYOff;
312 :
313 0 : VSIFSeekL(this->fp, nByteNum, SEEK_SET);
314 : int nReadSize =
315 0 : (GDALGetDataTypeSize(eDataType) / 8) * poDS->GetRasterXSize();
316 0 : VSIFReadL((char *)pImage, 1, nReadSize, this->fp);
317 :
318 : #ifdef CPL_LSB
319 0 : GDALSwapWords(pImage, 4, nBlockXSize * 2, 4);
320 : #endif
321 0 : return CE_None;
322 : }
323 :
324 : /********************************************************************
325 : * ================================================================ *
326 : * Implementation of the COASPDataset Class *
327 : * ================================================================ *
328 : ********************************************************************/
329 :
330 : /************************************************************************/
331 : /* ~COASPDataset() */
332 : /************************************************************************/
333 :
334 0 : COASPDataset::~COASPDataset()
335 : {
336 0 : CPLFree(pszFileName);
337 0 : if (fpHdr)
338 0 : VSIFCloseL(fpHdr);
339 0 : if (fpBinHH)
340 0 : VSIFCloseL(fpBinHH);
341 0 : if (fpBinHV)
342 0 : VSIFCloseL(fpBinHV);
343 0 : if (fpBinVH)
344 0 : VSIFCloseL(fpBinVH);
345 0 : if (fpBinVV)
346 0 : VSIFCloseL(fpBinVV);
347 0 : }
348 :
349 : /************************************************************************/
350 : /* Identify() */
351 : /************************************************************************/
352 :
353 52819 : int COASPDataset::Identify(GDALOpenInfo *poOpenInfo)
354 : {
355 52819 : if (poOpenInfo->fpL == nullptr || poOpenInfo->nHeaderBytes < 256)
356 49154 : return 0;
357 :
358 : // With a COASP .hdr file, the first line or so is: time_first_datarec
359 :
360 3665 : if (STARTS_WITH_CI((char *)poOpenInfo->pabyHeader, "time_first_datarec"))
361 0 : return 1;
362 :
363 3665 : return 0;
364 : }
365 :
366 : /************************************************************************/
367 : /* Open() */
368 : /************************************************************************/
369 :
370 0 : GDALDataset *COASPDataset::Open(GDALOpenInfo *poOpenInfo)
371 : {
372 0 : if (!COASPDataset::Identify(poOpenInfo))
373 0 : return nullptr;
374 :
375 : /* -------------------------------------------------------------------- */
376 : /* Confirm the requested access is supported. */
377 : /* -------------------------------------------------------------------- */
378 0 : if (poOpenInfo->eAccess == GA_Update)
379 : {
380 0 : CPLError(CE_Failure, CPLE_NotSupported,
381 : "The COASP driver does not support update access to existing"
382 : " datasets.\n");
383 0 : return nullptr;
384 : }
385 :
386 : /* Create a fresh dataset for us to work with */
387 0 : COASPDataset *poDS = new COASPDataset();
388 :
389 : /* Steal the file pointer for the header */
390 0 : poDS->fpHdr = poOpenInfo->fpL;
391 0 : poOpenInfo->fpL = nullptr;
392 :
393 0 : poDS->pszFileName = VSIStrdup(poOpenInfo->pszFilename);
394 :
395 : /* determine the file name prefix */
396 0 : char *pszBaseName = VSIStrdup(CPLGetBasename(poDS->pszFileName));
397 0 : char *pszDir = VSIStrdup(CPLGetPath(poDS->pszFileName));
398 0 : const char *pszExt = "rc";
399 0 : int nNull = static_cast<int>(strlen(pszBaseName)) - 1;
400 0 : if (nNull <= 0)
401 : {
402 0 : VSIFree(pszDir);
403 0 : VSIFree(pszBaseName);
404 0 : delete poDS;
405 0 : return nullptr;
406 : }
407 0 : char *pszBase = (char *)CPLMalloc(nNull);
408 0 : strncpy(pszBase, pszBaseName, nNull);
409 0 : pszBase[nNull - 1] = '\0';
410 0 : VSIFree(pszBaseName);
411 :
412 0 : char *psChan = strstr(pszBase, "hh");
413 0 : if (psChan == nullptr)
414 : {
415 0 : psChan = strstr(pszBase, "hv");
416 : }
417 0 : if (psChan == nullptr)
418 : {
419 0 : psChan = strstr(pszBase, "vh");
420 : }
421 0 : if (psChan == nullptr)
422 : {
423 0 : psChan = strstr(pszBase, "vv");
424 : }
425 :
426 0 : if (psChan == nullptr)
427 : {
428 0 : CPLError(CE_Failure, CPLE_AppDefined,
429 : "Unable to recognize file as COASP.");
430 0 : VSIFree(pszBase);
431 0 : VSIFree(pszDir);
432 0 : delete poDS;
433 0 : return nullptr;
434 : }
435 :
436 : /* Read Metadata, set GCPs as is appropriate */
437 0 : COASPMetadataReader oReader(poDS->pszFileName);
438 :
439 : /* Get Image X and Y widths */
440 0 : oReader.GotoMetadataItem("number_lines");
441 0 : COASPMetadataItem *poItem = oReader.GetNextItem();
442 0 : if (poItem == nullptr)
443 : {
444 0 : VSIFree(pszBase);
445 0 : VSIFree(pszDir);
446 0 : delete poDS;
447 0 : return nullptr;
448 : }
449 0 : char *nValue = poItem->GetItemValue();
450 0 : poDS->nRasterYSize = atoi(nValue);
451 0 : delete poItem;
452 0 : VSIFree(nValue);
453 :
454 0 : oReader.GotoMetadataItem("number_samples");
455 0 : poItem = oReader.GetNextItem();
456 0 : if (poItem == nullptr)
457 : {
458 0 : VSIFree(pszBase);
459 0 : VSIFree(pszDir);
460 0 : delete poDS;
461 0 : return nullptr;
462 : }
463 0 : nValue = poItem->GetItemValue();
464 0 : poDS->nRasterXSize = atoi(nValue);
465 0 : delete poItem;
466 0 : VSIFree(nValue);
467 :
468 0 : if (!GDALCheckDatasetDimensions(poDS->nRasterXSize, poDS->nRasterYSize))
469 : {
470 0 : VSIFree(pszBase);
471 0 : VSIFree(pszDir);
472 0 : delete poDS;
473 0 : return nullptr;
474 : }
475 :
476 : /* Horizontal transmit, horizontal receive */
477 0 : psChan[0] = 'h';
478 0 : psChan[1] = 'h';
479 0 : const char *pszFilename = CPLFormFilename(pszDir, pszBase, pszExt);
480 :
481 0 : poDS->fpBinHH = VSIFOpenL(pszFilename, "r");
482 :
483 0 : if (poDS->fpBinHH != nullptr)
484 : {
485 : /* Set raster band */
486 0 : poDS->SetBand(
487 0 : 1, new COASPRasterBand(poDS, GDT_CFloat32, hh, poDS->fpBinHH));
488 : }
489 :
490 : /* Horizontal transmit, vertical receive */
491 0 : psChan[0] = 'h';
492 0 : psChan[1] = 'v';
493 0 : pszFilename = CPLFormFilename(pszDir, pszBase, pszExt);
494 :
495 0 : poDS->fpBinHV = VSIFOpenL(pszFilename, "r");
496 :
497 0 : if (poDS->fpBinHV != nullptr)
498 : {
499 0 : poDS->SetBand(
500 0 : 2, new COASPRasterBand(poDS, GDT_CFloat32, hv, poDS->fpBinHV));
501 : }
502 :
503 : /* Vertical transmit, horizontal receive */
504 0 : psChan[0] = 'v';
505 0 : psChan[1] = 'h';
506 0 : pszFilename = CPLFormFilename(pszDir, pszBase, pszExt);
507 :
508 0 : poDS->fpBinVH = VSIFOpenL(pszFilename, "r");
509 :
510 0 : if (poDS->fpBinVH != nullptr)
511 : {
512 0 : poDS->SetBand(
513 0 : 3, new COASPRasterBand(poDS, GDT_CFloat32, vh, poDS->fpBinVH));
514 : }
515 :
516 : /* Vertical transmit, vertical receive */
517 0 : psChan[0] = 'v';
518 0 : psChan[1] = 'v';
519 0 : pszFilename = CPLFormFilename(pszDir, pszBase, pszExt);
520 :
521 0 : poDS->fpBinVV = VSIFOpenL(pszFilename, "r");
522 :
523 0 : if (poDS->fpBinVV != nullptr)
524 : {
525 0 : poDS->SetBand(
526 0 : 4, new COASPRasterBand(poDS, GDT_CFloat32, vv, poDS->fpBinVV));
527 : }
528 :
529 : /* Oops, missing all the data? */
530 0 : if (poDS->fpBinHH == nullptr && poDS->fpBinHV == nullptr &&
531 0 : poDS->fpBinVH == nullptr && poDS->fpBinVV == nullptr)
532 : {
533 0 : CPLError(CE_Failure, CPLE_AppDefined, "Unable to find any data!");
534 0 : VSIFree(pszBase);
535 0 : VSIFree(pszDir);
536 0 : delete poDS;
537 0 : return nullptr;
538 : }
539 :
540 0 : if (poDS->GetRasterCount() == 4)
541 : {
542 0 : poDS->SetMetadataItem("MATRIX_REPRESENTATION", "SCATTERING");
543 : }
544 :
545 0 : VSIFree(pszBase);
546 0 : VSIFree(pszDir);
547 :
548 0 : return poDS;
549 : }
550 :
551 : /************************************************************************/
552 : /* GDALRegister_COASP() */
553 : /************************************************************************/
554 :
555 1595 : void GDALRegister_COASP()
556 : {
557 1595 : if (GDALGetDriverByName("COASP") != nullptr)
558 302 : return;
559 :
560 1293 : GDALDriver *poDriver = new GDALDriver();
561 :
562 1293 : poDriver->SetDescription("COASP");
563 1293 : poDriver->SetMetadataItem(GDAL_DCAP_RASTER, "YES");
564 1293 : poDriver->SetMetadataItem(GDAL_DMD_LONGNAME,
565 1293 : "DRDC COASP SAR Processor Raster");
566 1293 : poDriver->SetMetadataItem(GDAL_DMD_EXTENSION, "hdr");
567 1293 : poDriver->SetMetadataItem(GDAL_DMD_HELPTOPIC, "drivers/raster/coasp.html");
568 1293 : poDriver->pfnIdentify = COASPDataset::Identify;
569 1293 : poDriver->pfnOpen = COASPDataset::Open;
570 1293 : GetGDALDriverManager()->RegisterDriver(poDriver);
571 : }
|