Line data Source code
1 : /******************************************************************************
2 : *
3 : * Project: GView
4 : * Purpose: Implementation of Atlantis HKV labelled blob support
5 : * Author: Frank Warmerdam, warmerdam@pobox.com
6 : *
7 : ******************************************************************************
8 : * Copyright (c) 2000, Frank Warmerdam
9 : * Copyright (c) 2008-2012, Even Rouault <even dot rouault at spatialys.com>
10 : *
11 : * SPDX-License-Identifier: MIT
12 : ****************************************************************************/
13 :
14 : #include <ctype.h>
15 :
16 : #include "atlsci_spheroid.h"
17 : #include "cpl_string.h"
18 : #include "gdal_frmts.h"
19 : #include "ogr_spatialref.h"
20 : #include "rawdataset.h"
21 :
22 : #include <cmath>
23 :
24 : #include <algorithm>
25 :
26 : /************************************************************************/
27 : /* ==================================================================== */
28 : /* HKVRasterBand */
29 : /* ==================================================================== */
30 : /************************************************************************/
31 :
32 : class HKVDataset;
33 :
34 : class HKVRasterBand final : public RawRasterBand
35 : {
36 : friend class HKVDataset;
37 :
38 : public:
39 : HKVRasterBand(HKVDataset *poDS, int nBand, VSILFILE *fpRaw,
40 : unsigned int nImgOffset, int nPixelOffset, int nLineOffset,
41 : GDALDataType eDataType, int bNativeOrder);
42 :
43 188 : ~HKVRasterBand() override
44 94 : {
45 188 : }
46 :
47 : CPLErr SetNoDataValue(double) override;
48 : };
49 :
50 : /************************************************************************/
51 : /* HKV Spheroids */
52 : /************************************************************************/
53 :
54 : class HKVSpheroidList : public SpheroidList
55 : {
56 : public:
57 : HKVSpheroidList();
58 :
59 44 : ~HKVSpheroidList()
60 44 : {
61 44 : }
62 : };
63 :
64 44 : HKVSpheroidList ::HKVSpheroidList()
65 : {
66 44 : num_spheroids = 58;
67 44 : epsilonR = 0.1;
68 44 : epsilonI = 0.000001;
69 :
70 44 : spheroids[0].SetValuesByEqRadiusAndInvFlattening("airy-1830", 6377563.396,
71 : 299.3249646);
72 44 : spheroids[1].SetValuesByEqRadiusAndInvFlattening("modified-airy",
73 : 6377340.189, 299.3249646);
74 44 : spheroids[2].SetValuesByEqRadiusAndInvFlattening("australian-national",
75 : 6378160, 298.25);
76 44 : spheroids[3].SetValuesByEqRadiusAndInvFlattening("bessel-1841-namibia",
77 : 6377483.865, 299.1528128);
78 44 : spheroids[4].SetValuesByEqRadiusAndInvFlattening("bessel-1841", 6377397.155,
79 : 299.1528128);
80 44 : spheroids[5].SetValuesByEqRadiusAndInvFlattening("clarke-1858", 6378294.0,
81 : 294.297);
82 44 : spheroids[6].SetValuesByEqRadiusAndInvFlattening("clarke-1866", 6378206.4,
83 : 294.9786982);
84 44 : spheroids[7].SetValuesByEqRadiusAndInvFlattening("clarke-1880", 6378249.145,
85 : 293.465);
86 44 : spheroids[8].SetValuesByEqRadiusAndInvFlattening("everest-india-1830",
87 : 6377276.345, 300.8017);
88 44 : spheroids[9].SetValuesByEqRadiusAndInvFlattening("everest-sabah-sarawak",
89 : 6377298.556, 300.8017);
90 44 : spheroids[10].SetValuesByEqRadiusAndInvFlattening("everest-india-1956",
91 : 6377301.243, 300.8017);
92 44 : spheroids[11].SetValuesByEqRadiusAndInvFlattening("everest-malaysia-1969",
93 : 6377295.664, 300.8017);
94 44 : spheroids[12].SetValuesByEqRadiusAndInvFlattening("everest-malay-sing",
95 : 6377304.063, 300.8017);
96 44 : spheroids[13].SetValuesByEqRadiusAndInvFlattening("everest-pakistan",
97 : 6377309.613, 300.8017);
98 44 : spheroids[14].SetValuesByEqRadiusAndInvFlattening("modified-fisher-1960",
99 : 6378155, 298.3);
100 44 : spheroids[15].SetValuesByEqRadiusAndInvFlattening("helmert-1906", 6378200,
101 : 298.3);
102 44 : spheroids[16].SetValuesByEqRadiusAndInvFlattening("hough-1960", 6378270,
103 : 297);
104 44 : spheroids[17].SetValuesByEqRadiusAndInvFlattening("hughes", 6378273.0,
105 : 298.279);
106 44 : spheroids[18].SetValuesByEqRadiusAndInvFlattening("indonesian-1974",
107 : 6378160, 298.247);
108 44 : spheroids[19].SetValuesByEqRadiusAndInvFlattening("international-1924",
109 : 6378388, 297);
110 44 : spheroids[20].SetValuesByEqRadiusAndInvFlattening("iugc-67", 6378160.0,
111 : 298.254);
112 44 : spheroids[21].SetValuesByEqRadiusAndInvFlattening("iugc-75", 6378140.0,
113 : 298.25298);
114 44 : spheroids[22].SetValuesByEqRadiusAndInvFlattening("krassovsky-1940",
115 : 6378245, 298.3);
116 44 : spheroids[23].SetValuesByEqRadiusAndInvFlattening("kaula", 6378165.0,
117 : 292.308);
118 44 : spheroids[24].SetValuesByEqRadiusAndInvFlattening("grs-80", 6378137,
119 : 298.257222101);
120 44 : spheroids[25].SetValuesByEqRadiusAndInvFlattening("south-american-1969",
121 : 6378160, 298.25);
122 44 : spheroids[26].SetValuesByEqRadiusAndInvFlattening("wgs-72", 6378135,
123 : 298.26);
124 44 : spheroids[27].SetValuesByEqRadiusAndInvFlattening("wgs-84", 6378137,
125 : 298.257223563);
126 44 : spheroids[28].SetValuesByEqRadiusAndInvFlattening("ev-wgs-84", 6378137.0,
127 : 298.252841);
128 44 : spheroids[29].SetValuesByEqRadiusAndInvFlattening("ev-bessel", 6377397.0,
129 : 299.1976073);
130 :
131 44 : spheroids[30].SetValuesByEqRadiusAndInvFlattening("airy_1830", 6377563.396,
132 : 299.3249646);
133 44 : spheroids[31].SetValuesByEqRadiusAndInvFlattening("modified_airy",
134 : 6377340.189, 299.3249646);
135 44 : spheroids[32].SetValuesByEqRadiusAndInvFlattening("australian_national",
136 : 6378160, 298.25);
137 44 : spheroids[33].SetValuesByEqRadiusAndInvFlattening("bessel_1841_namibia",
138 : 6377483.865, 299.1528128);
139 44 : spheroids[34].SetValuesByEqRadiusAndInvFlattening("bessel_1841",
140 : 6377397.155, 299.1528128);
141 44 : spheroids[35].SetValuesByEqRadiusAndInvFlattening("clarke_1858", 6378294.0,
142 : 294.297);
143 44 : spheroids[36].SetValuesByEqRadiusAndInvFlattening("clarke_1866", 6378206.4,
144 : 294.9786982);
145 44 : spheroids[37].SetValuesByEqRadiusAndInvFlattening("clarke_1880",
146 : 6378249.145, 293.465);
147 44 : spheroids[38].SetValuesByEqRadiusAndInvFlattening("everest_india_1830",
148 : 6377276.345, 300.8017);
149 44 : spheroids[39].SetValuesByEqRadiusAndInvFlattening("everest_sabah_sarawak",
150 : 6377298.556, 300.8017);
151 44 : spheroids[40].SetValuesByEqRadiusAndInvFlattening("everest_india_1956",
152 : 6377301.243, 300.8017);
153 44 : spheroids[41].SetValuesByEqRadiusAndInvFlattening("everest_malaysia_1969",
154 : 6377295.664, 300.8017);
155 44 : spheroids[42].SetValuesByEqRadiusAndInvFlattening("everest_malay_sing",
156 : 6377304.063, 300.8017);
157 44 : spheroids[43].SetValuesByEqRadiusAndInvFlattening("everest_pakistan",
158 : 6377309.613, 300.8017);
159 44 : spheroids[44].SetValuesByEqRadiusAndInvFlattening("modified_fisher_1960",
160 : 6378155, 298.3);
161 44 : spheroids[45].SetValuesByEqRadiusAndInvFlattening("helmert_1906", 6378200,
162 : 298.3);
163 44 : spheroids[46].SetValuesByEqRadiusAndInvFlattening("hough_1960", 6378270,
164 : 297);
165 44 : spheroids[47].SetValuesByEqRadiusAndInvFlattening("indonesian_1974",
166 : 6378160, 298.247);
167 44 : spheroids[48].SetValuesByEqRadiusAndInvFlattening("international_1924",
168 : 6378388, 297);
169 44 : spheroids[49].SetValuesByEqRadiusAndInvFlattening("iugc_67", 6378160.0,
170 : 298.254);
171 44 : spheroids[50].SetValuesByEqRadiusAndInvFlattening("iugc_75", 6378140.0,
172 : 298.25298);
173 44 : spheroids[51].SetValuesByEqRadiusAndInvFlattening("krassovsky_1940",
174 : 6378245, 298.3);
175 44 : spheroids[52].SetValuesByEqRadiusAndInvFlattening("grs_80", 6378137,
176 : 298.257222101);
177 44 : spheroids[53].SetValuesByEqRadiusAndInvFlattening("south_american_1969",
178 : 6378160, 298.25);
179 44 : spheroids[54].SetValuesByEqRadiusAndInvFlattening("wgs_72", 6378135,
180 : 298.26);
181 44 : spheroids[55].SetValuesByEqRadiusAndInvFlattening("wgs_84", 6378137,
182 : 298.257223563);
183 44 : spheroids[56].SetValuesByEqRadiusAndInvFlattening("ev_wgs_84", 6378137.0,
184 : 298.252841);
185 44 : spheroids[57].SetValuesByEqRadiusAndInvFlattening("ev_bessel", 6377397.0,
186 : 299.1976073);
187 44 : }
188 :
189 : CPLErr SaveHKVAttribFile(const char *pszFilenameIn, int nXSize, int nYSize,
190 : int nBands, GDALDataType eType, int bNoDataSet,
191 : double dfNoDataValue);
192 :
193 : /************************************************************************/
194 : /* ==================================================================== */
195 : /* HKVDataset */
196 : /* ==================================================================== */
197 : /************************************************************************/
198 :
199 : class HKVDataset final : public RawDataset
200 : {
201 : friend class HKVRasterBand;
202 :
203 : char *pszPath;
204 : VSILFILE *fpBlob;
205 :
206 : int nGCPCount;
207 : GDAL_GCP *pasGCPList;
208 :
209 : void ProcessGeoref(const char *);
210 : void ProcessGeorefGCP(char **, const char *, double, double);
211 :
212 44 : void SetVersion(float version_number)
213 : {
214 : // Update stored info.
215 44 : MFF2version = version_number;
216 44 : }
217 :
218 : float MFF2version;
219 :
220 : GDALDataType eRasterType;
221 :
222 : void SetNoDataValue(double);
223 :
224 : OGRSpatialReference m_oSRS{};
225 : OGRSpatialReference m_oGCPSRS{};
226 : double adfGeoTransform[6];
227 :
228 : char **papszAttrib;
229 :
230 : bool bGeorefChanged;
231 : char **papszGeoref;
232 :
233 : // NOTE: The MFF2 format goes against GDAL's API in that nodata values are
234 : // set per-dataset rather than per-band. To compromise, for writing out,
235 : // the dataset's nodata value will be set to the last value set on any of
236 : // the raster bands.
237 :
238 : bool bNoDataSet;
239 : bool bNoDataChanged;
240 : double dfNoDataValue;
241 :
242 : CPL_DISALLOW_COPY_ASSIGN(HKVDataset)
243 :
244 : CPLErr Close() override;
245 :
246 : public:
247 : HKVDataset();
248 : ~HKVDataset() override;
249 :
250 2 : int GetGCPCount() override /* const */
251 : {
252 2 : return nGCPCount;
253 : }
254 :
255 0 : const OGRSpatialReference *GetGCPSpatialRef() const override
256 : {
257 0 : return m_oGCPSRS.IsEmpty() ? nullptr : &m_oGCPSRS;
258 : }
259 :
260 : const GDAL_GCP *GetGCPs() override;
261 :
262 15 : const OGRSpatialReference *GetSpatialRef() const override
263 : {
264 15 : return m_oSRS.IsEmpty() ? nullptr : &m_oSRS;
265 : }
266 :
267 : CPLErr GetGeoTransform(double *) override;
268 :
269 : CPLErr SetGeoTransform(double *) override;
270 : CPLErr SetSpatialRef(const OGRSpatialReference *poSRS) override;
271 :
272 : static GDALDataset *Open(GDALOpenInfo *);
273 : static GDALDataset *Create(const char *pszFilename, int nXSize, int nYSize,
274 : int nBands, GDALDataType eType,
275 : char **papszParamList);
276 : static GDALDataset *CreateCopy(const char *pszFilename,
277 : GDALDataset *poSrcDS, int bStrict,
278 : char **papszOptions,
279 : GDALProgressFunc pfnProgress,
280 : void *pProgressData);
281 :
282 : static CPLErr Delete(const char *pszName);
283 : };
284 :
285 : /************************************************************************/
286 : /* ==================================================================== */
287 : /* HKVRasterBand */
288 : /* ==================================================================== */
289 : /************************************************************************/
290 :
291 : /************************************************************************/
292 : /* HKVRasterBand() */
293 : /************************************************************************/
294 :
295 94 : HKVRasterBand::HKVRasterBand(HKVDataset *poDSIn, int nBandIn, VSILFILE *fpRawIn,
296 : unsigned int nImgOffsetIn, int nPixelOffsetIn,
297 : int nLineOffsetIn, GDALDataType eDataTypeIn,
298 94 : int bNativeOrderIn)
299 : : RawRasterBand(GDALDataset::FromHandle(poDSIn), nBandIn, fpRawIn,
300 : nImgOffsetIn, nPixelOffsetIn, nLineOffsetIn, eDataTypeIn,
301 94 : bNativeOrderIn, RawRasterBand::OwnFP::NO)
302 :
303 : {
304 94 : poDS = poDSIn;
305 94 : nBand = nBandIn;
306 :
307 94 : nBlockXSize = poDS->GetRasterXSize();
308 94 : nBlockYSize = 1;
309 94 : }
310 :
311 : /************************************************************************/
312 : /* SetNoDataValue() */
313 : /************************************************************************/
314 :
315 0 : CPLErr HKVRasterBand::SetNoDataValue(double dfNewValue)
316 :
317 : {
318 0 : HKVDataset *poHKVDS = reinterpret_cast<HKVDataset *>(poDS);
319 0 : RawRasterBand::SetNoDataValue(dfNewValue);
320 0 : poHKVDS->SetNoDataValue(dfNewValue);
321 :
322 0 : return CE_None;
323 : }
324 :
325 : /************************************************************************/
326 : /* ==================================================================== */
327 : /* HKVDataset */
328 : /* ==================================================================== */
329 : /************************************************************************/
330 :
331 : /************************************************************************/
332 : /* HKVDataset() */
333 : /************************************************************************/
334 :
335 44 : HKVDataset::HKVDataset()
336 : : pszPath(nullptr), fpBlob(nullptr), nGCPCount(0), pasGCPList(nullptr),
337 : // Initialize datasets to new version; change if necessary.
338 : MFF2version(1.1f), eRasterType(GDT_Unknown), papszAttrib(nullptr),
339 : bGeorefChanged(false), papszGeoref(nullptr), bNoDataSet(false),
340 44 : bNoDataChanged(false), dfNoDataValue(0.0)
341 : {
342 44 : m_oSRS.SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
343 44 : m_oGCPSRS.SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
344 44 : adfGeoTransform[0] = 0.0;
345 44 : adfGeoTransform[1] = 1.0;
346 44 : adfGeoTransform[2] = 0.0;
347 44 : adfGeoTransform[3] = 0.0;
348 44 : adfGeoTransform[4] = 0.0;
349 44 : adfGeoTransform[5] = 1.0;
350 44 : }
351 :
352 : /************************************************************************/
353 : /* ~HKVDataset() */
354 : /************************************************************************/
355 :
356 88 : HKVDataset::~HKVDataset()
357 :
358 : {
359 44 : HKVDataset::Close();
360 88 : }
361 :
362 : /************************************************************************/
363 : /* Close() */
364 : /************************************************************************/
365 :
366 88 : CPLErr HKVDataset::Close()
367 : {
368 88 : CPLErr eErr = CE_None;
369 88 : if (nOpenFlags != OPEN_FLAGS_CLOSED)
370 : {
371 44 : if (HKVDataset::FlushCache(true) != CE_None)
372 0 : eErr = CE_Failure;
373 :
374 44 : if (bGeorefChanged)
375 : {
376 : const char *pszFilename =
377 26 : CPLFormFilename(pszPath, "georef", nullptr);
378 26 : CSLSave(papszGeoref, pszFilename);
379 : }
380 :
381 44 : if (bNoDataChanged)
382 : {
383 0 : SaveHKVAttribFile(pszPath, nRasterXSize, nRasterYSize, nBands,
384 0 : eRasterType, bNoDataSet, dfNoDataValue);
385 : }
386 :
387 44 : if (fpBlob)
388 : {
389 44 : if (VSIFCloseL(fpBlob) != 0)
390 : {
391 0 : eErr = CE_Failure;
392 0 : CPLError(CE_Failure, CPLE_FileIO, "I/O error");
393 : }
394 : }
395 :
396 44 : if (nGCPCount > 0)
397 : {
398 14 : GDALDeinitGCPs(nGCPCount, pasGCPList);
399 14 : CPLFree(pasGCPList);
400 : }
401 :
402 44 : CPLFree(pszPath);
403 44 : CSLDestroy(papszGeoref);
404 44 : CSLDestroy(papszAttrib);
405 :
406 44 : if (GDALPamDataset::Close() != CE_None)
407 0 : eErr = CE_Failure;
408 : }
409 88 : return eErr;
410 : }
411 :
412 : /************************************************************************/
413 : /* SetNoDataValue() */
414 : /************************************************************************/
415 :
416 0 : void HKVDataset::SetNoDataValue(double dfNewValue)
417 :
418 : {
419 0 : bNoDataSet = true;
420 0 : bNoDataChanged = true;
421 0 : dfNoDataValue = dfNewValue;
422 0 : }
423 :
424 : /************************************************************************/
425 : /* SaveHKVAttribFile() */
426 : /************************************************************************/
427 :
428 26 : CPLErr SaveHKVAttribFile(const char *pszFilenameIn, int nXSize, int nYSize,
429 : int nBands, GDALDataType eType, int bNoDataSet,
430 : double dfNoDataValue)
431 :
432 : {
433 26 : const char *pszFilename = CPLFormFilename(pszFilenameIn, "attrib", nullptr);
434 :
435 26 : FILE *fp = VSIFOpen(pszFilename, "wt");
436 26 : if (fp == nullptr)
437 : {
438 0 : CPLError(CE_Failure, CPLE_OpenFailed, "Couldn't create %s.",
439 : pszFilename);
440 0 : return CE_Failure;
441 : }
442 :
443 26 : fprintf(fp, "channel.enumeration = %d\n", nBands);
444 26 : fprintf(fp, "channel.interleave = { *pixel tile sequential }\n");
445 26 : fprintf(fp, "extent.cols = %d\n", nXSize);
446 26 : fprintf(fp, "extent.rows = %d\n", nYSize);
447 :
448 26 : switch (eType)
449 : {
450 11 : case GDT_Byte:
451 11 : fprintf(fp, "pixel.encoding = "
452 : "{ *unsigned twos-complement ieee-754 }\n");
453 11 : break;
454 :
455 3 : case GDT_UInt16:
456 3 : fprintf(fp, "pixel.encoding = "
457 : "{ *unsigned twos-complement ieee-754 }\n");
458 3 : break;
459 :
460 6 : case GDT_CInt16:
461 : case GDT_Int16:
462 6 : fprintf(fp, "pixel.encoding = "
463 : "{ unsigned *twos-complement ieee-754 }\n");
464 6 : break;
465 :
466 6 : case GDT_CFloat32:
467 : case GDT_Float32:
468 6 : fprintf(fp, "pixel.encoding = "
469 : "{ unsigned twos-complement *ieee-754 }\n");
470 6 : break;
471 :
472 0 : default:
473 0 : CPLAssert(false);
474 : }
475 :
476 26 : fprintf(fp, "pixel.size = %d\n", GDALGetDataTypeSizeBits(eType));
477 26 : if (GDALDataTypeIsComplex(eType))
478 6 : fprintf(fp, "pixel.field = { real *complex }\n");
479 : else
480 20 : fprintf(fp, "pixel.field = { *real complex }\n");
481 :
482 : #ifdef CPL_MSB
483 : fprintf(fp, "pixel.order = { lsbf *msbf }\n");
484 : #else
485 26 : fprintf(fp, "pixel.order = { *lsbf msbf }\n");
486 : #endif
487 :
488 26 : if (bNoDataSet)
489 0 : fprintf(fp, "pixel.no_data = %s\n", CPLSPrintf("%f", dfNoDataValue));
490 :
491 : // Version information- only create the new style.
492 26 : fprintf(fp, "version = 1.1");
493 :
494 26 : if (VSIFClose(fp) != 0)
495 0 : return CE_Failure;
496 26 : return CE_None;
497 : }
498 :
499 : /************************************************************************/
500 : /* GetGeoTransform() */
501 : /************************************************************************/
502 :
503 29 : CPLErr HKVDataset::GetGeoTransform(double *padfTransform)
504 :
505 : {
506 29 : memcpy(padfTransform, adfGeoTransform, sizeof(double) * 6);
507 29 : return CE_None;
508 : }
509 :
510 : /************************************************************************/
511 : /* SetGeoTransform() */
512 : /************************************************************************/
513 :
514 26 : CPLErr HKVDataset::SetGeoTransform(double *padfTransform)
515 :
516 : {
517 : // NOTE: Geotransform coordinates must match the current projection
518 : // of the dataset being changed (not the geotransform source).
519 : // i.e. be in lat/longs for LL projected; UTM for UTM projected.
520 : // SET PROJECTION BEFORE SETTING GEOTRANSFORM TO AVOID SYNCHRONIZATION
521 : // PROBLEMS.
522 :
523 : // Update the geotransform itself.
524 26 : memcpy(adfGeoTransform, padfTransform, sizeof(double) * 6);
525 :
526 : // Clear previous gcps.
527 26 : if (nGCPCount > 0)
528 : {
529 0 : GDALDeinitGCPs(nGCPCount, pasGCPList);
530 0 : CPLFree(pasGCPList);
531 : }
532 26 : nGCPCount = 0;
533 26 : pasGCPList = nullptr;
534 :
535 : // Return if the identity transform is set.
536 26 : if (adfGeoTransform[0] == 0.0 && adfGeoTransform[1] == 1.0 &&
537 0 : adfGeoTransform[2] == 0.0 && adfGeoTransform[3] == 0.0 &&
538 0 : adfGeoTransform[4] == 0.0 && adfGeoTransform[5] == 1.0)
539 0 : return CE_None;
540 :
541 : // Update georef text info for saving later, and update GCPs to match
542 : // geotransform.
543 :
544 26 : OGRCoordinateTransformation *poTransform = nullptr;
545 26 : bool bSuccess = true;
546 :
547 : // Projection parameter checking will have been done in SetProjection.
548 37 : if ((CSLFetchNameValue(papszGeoref, "projection.name") != nullptr) &&
549 11 : (EQUAL(CSLFetchNameValue(papszGeoref, "projection.name"), "UTM")))
550 : {
551 1 : auto poLLSRS = m_oSRS.CloneGeogCS();
552 1 : if (poLLSRS)
553 : {
554 1 : poLLSRS->SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
555 1 : poTransform = OGRCreateCoordinateTransformation(&m_oSRS, poLLSRS);
556 1 : delete poLLSRS;
557 1 : if (poTransform == nullptr)
558 : {
559 0 : bSuccess = false;
560 0 : CPLErrorReset();
561 : }
562 : }
563 : else
564 : {
565 0 : bSuccess = false;
566 : }
567 : }
568 25 : else if (((CSLFetchNameValue(papszGeoref, "projection.name") != nullptr) &&
569 10 : (!EQUAL(CSLFetchNameValue(papszGeoref, "projection.name"),
570 50 : "LL"))) ||
571 25 : (CSLFetchNameValue(papszGeoref, "projection.name") == nullptr))
572 : {
573 15 : return CE_Failure;
574 : }
575 :
576 11 : nGCPCount = 0;
577 11 : pasGCPList = static_cast<GDAL_GCP *>(CPLCalloc(sizeof(GDAL_GCP), 5));
578 :
579 : /* -------------------------------------------------------------------- */
580 : /* top left */
581 : /* -------------------------------------------------------------------- */
582 11 : GDALInitGCPs(1, pasGCPList + nGCPCount);
583 11 : CPLFree(pasGCPList[nGCPCount].pszId);
584 11 : pasGCPList[nGCPCount].pszId = CPLStrdup("top_left");
585 :
586 11 : double temp_lat = 0.0;
587 11 : double temp_long = 0.0;
588 11 : if (MFF2version > 1.0)
589 : {
590 11 : temp_lat = padfTransform[3];
591 11 : temp_long = padfTransform[0];
592 11 : pasGCPList[nGCPCount].dfGCPPixel = 0.0;
593 11 : pasGCPList[nGCPCount].dfGCPLine = 0.0;
594 : }
595 : else
596 : {
597 0 : temp_lat =
598 0 : padfTransform[3] + 0.5 * padfTransform[4] + 0.5 * padfTransform[5];
599 0 : temp_long =
600 0 : padfTransform[0] + 0.5 * padfTransform[1] + 0.5 * padfTransform[2];
601 0 : pasGCPList[nGCPCount].dfGCPPixel = 0.5;
602 0 : pasGCPList[nGCPCount].dfGCPLine = 0.5;
603 : }
604 11 : pasGCPList[nGCPCount].dfGCPX = temp_long;
605 11 : pasGCPList[nGCPCount].dfGCPY = temp_lat;
606 11 : pasGCPList[nGCPCount].dfGCPZ = 0.0;
607 11 : nGCPCount++;
608 :
609 11 : if (poTransform != nullptr)
610 : {
611 1 : if (!bSuccess || !poTransform->Transform(1, &temp_long, &temp_lat))
612 0 : bSuccess = false;
613 : }
614 :
615 11 : if (bSuccess)
616 : {
617 11 : char szValue[128] = {'\0'};
618 11 : CPLsnprintf(szValue, sizeof(szValue), "%.10f", temp_lat);
619 11 : papszGeoref =
620 11 : CSLSetNameValue(papszGeoref, "top_left.latitude", szValue);
621 :
622 11 : CPLsnprintf(szValue, sizeof(szValue), "%.10f", temp_long);
623 11 : papszGeoref =
624 11 : CSLSetNameValue(papszGeoref, "top_left.longitude", szValue);
625 : }
626 :
627 : /* -------------------------------------------------------------------- */
628 : /* top_right */
629 : /* -------------------------------------------------------------------- */
630 11 : GDALInitGCPs(1, pasGCPList + nGCPCount);
631 11 : CPLFree(pasGCPList[nGCPCount].pszId);
632 11 : pasGCPList[nGCPCount].pszId = CPLStrdup("top_right");
633 :
634 11 : if (MFF2version > 1.0)
635 : {
636 11 : temp_lat = padfTransform[3] + GetRasterXSize() * padfTransform[4];
637 11 : temp_long = padfTransform[0] + GetRasterXSize() * padfTransform[1];
638 11 : pasGCPList[nGCPCount].dfGCPPixel = GetRasterXSize();
639 11 : pasGCPList[nGCPCount].dfGCPLine = 0.0;
640 : }
641 : else
642 : {
643 0 : temp_lat = padfTransform[3] +
644 0 : (GetRasterXSize() - 0.5) * padfTransform[4] +
645 0 : 0.5 * padfTransform[5];
646 0 : temp_long = padfTransform[0] +
647 0 : (GetRasterXSize() - 0.5) * padfTransform[1] +
648 0 : 0.5 * padfTransform[2];
649 0 : pasGCPList[nGCPCount].dfGCPPixel = GetRasterXSize() - 0.5;
650 0 : pasGCPList[nGCPCount].dfGCPLine = 0.5;
651 : }
652 11 : pasGCPList[nGCPCount].dfGCPX = temp_long;
653 11 : pasGCPList[nGCPCount].dfGCPY = temp_lat;
654 11 : pasGCPList[nGCPCount].dfGCPZ = 0.0;
655 11 : nGCPCount++;
656 :
657 11 : if (poTransform != nullptr)
658 : {
659 1 : if (!bSuccess || !poTransform->Transform(1, &temp_long, &temp_lat))
660 0 : bSuccess = false;
661 : }
662 :
663 11 : if (bSuccess)
664 : {
665 11 : char szValue[128] = {'\0'};
666 11 : CPLsnprintf(szValue, sizeof(szValue), "%.10f", temp_lat);
667 11 : papszGeoref =
668 11 : CSLSetNameValue(papszGeoref, "top_right.latitude", szValue);
669 :
670 11 : CPLsnprintf(szValue, sizeof(szValue), "%.10f", temp_long);
671 11 : papszGeoref =
672 11 : CSLSetNameValue(papszGeoref, "top_right.longitude", szValue);
673 : }
674 :
675 : /* -------------------------------------------------------------------- */
676 : /* bottom_left */
677 : /* -------------------------------------------------------------------- */
678 11 : GDALInitGCPs(1, pasGCPList + nGCPCount);
679 11 : CPLFree(pasGCPList[nGCPCount].pszId);
680 11 : pasGCPList[nGCPCount].pszId = CPLStrdup("bottom_left");
681 :
682 11 : if (MFF2version > 1.0)
683 : {
684 11 : temp_lat = padfTransform[3] + GetRasterYSize() * padfTransform[5];
685 11 : temp_long = padfTransform[0] + GetRasterYSize() * padfTransform[2];
686 11 : pasGCPList[nGCPCount].dfGCPPixel = 0.0;
687 11 : pasGCPList[nGCPCount].dfGCPLine = GetRasterYSize();
688 : }
689 : else
690 : {
691 0 : temp_lat = padfTransform[3] + 0.5 * padfTransform[4] +
692 0 : (GetRasterYSize() - 0.5) * padfTransform[5];
693 0 : temp_long = padfTransform[0] + 0.5 * padfTransform[1] +
694 0 : (GetRasterYSize() - 0.5) * padfTransform[2];
695 0 : pasGCPList[nGCPCount].dfGCPPixel = 0.5;
696 0 : pasGCPList[nGCPCount].dfGCPLine = GetRasterYSize() - 0.5;
697 : }
698 11 : pasGCPList[nGCPCount].dfGCPX = temp_long;
699 11 : pasGCPList[nGCPCount].dfGCPY = temp_lat;
700 11 : pasGCPList[nGCPCount].dfGCPZ = 0.0;
701 11 : nGCPCount++;
702 :
703 11 : if (poTransform != nullptr)
704 : {
705 1 : if (!bSuccess || !poTransform->Transform(1, &temp_long, &temp_lat))
706 0 : bSuccess = false;
707 : }
708 :
709 11 : if (bSuccess)
710 : {
711 11 : char szValue[128] = {'\0'};
712 11 : CPLsnprintf(szValue, sizeof(szValue), "%.10f", temp_lat);
713 11 : papszGeoref =
714 11 : CSLSetNameValue(papszGeoref, "bottom_left.latitude", szValue);
715 :
716 11 : CPLsnprintf(szValue, sizeof(szValue), "%.10f", temp_long);
717 11 : papszGeoref =
718 11 : CSLSetNameValue(papszGeoref, "bottom_left.longitude", szValue);
719 : }
720 :
721 : /* -------------------------------------------------------------------- */
722 : /* bottom_right */
723 : /* -------------------------------------------------------------------- */
724 11 : GDALInitGCPs(1, pasGCPList + nGCPCount);
725 11 : CPLFree(pasGCPList[nGCPCount].pszId);
726 11 : pasGCPList[nGCPCount].pszId = CPLStrdup("bottom_right");
727 :
728 11 : if (MFF2version > 1.0)
729 : {
730 11 : temp_lat = padfTransform[3] + GetRasterXSize() * padfTransform[4] +
731 11 : GetRasterYSize() * padfTransform[5];
732 11 : temp_long = padfTransform[0] + GetRasterXSize() * padfTransform[1] +
733 11 : GetRasterYSize() * padfTransform[2];
734 11 : pasGCPList[nGCPCount].dfGCPPixel = GetRasterXSize();
735 11 : pasGCPList[nGCPCount].dfGCPLine = GetRasterYSize();
736 : }
737 : else
738 : {
739 0 : temp_lat = padfTransform[3] +
740 0 : (GetRasterXSize() - 0.5) * padfTransform[4] +
741 0 : (GetRasterYSize() - 0.5) * padfTransform[5];
742 0 : temp_long = padfTransform[0] +
743 0 : (GetRasterXSize() - 0.5) * padfTransform[1] +
744 0 : (GetRasterYSize() - 0.5) * padfTransform[2];
745 0 : pasGCPList[nGCPCount].dfGCPPixel = GetRasterXSize() - 0.5;
746 0 : pasGCPList[nGCPCount].dfGCPLine = GetRasterYSize() - 0.5;
747 : }
748 11 : pasGCPList[nGCPCount].dfGCPX = temp_long;
749 11 : pasGCPList[nGCPCount].dfGCPY = temp_lat;
750 11 : pasGCPList[nGCPCount].dfGCPZ = 0.0;
751 11 : nGCPCount++;
752 :
753 11 : if (poTransform != nullptr)
754 : {
755 1 : if (!bSuccess || !poTransform->Transform(1, &temp_long, &temp_lat))
756 0 : bSuccess = false;
757 : }
758 :
759 11 : if (bSuccess)
760 : {
761 11 : char szValue[128] = {'\0'};
762 11 : CPLsnprintf(szValue, sizeof(szValue), "%.10f", temp_lat);
763 11 : papszGeoref =
764 11 : CSLSetNameValue(papszGeoref, "bottom_right.latitude", szValue);
765 :
766 11 : CPLsnprintf(szValue, sizeof(szValue), "%.10f", temp_long);
767 11 : papszGeoref =
768 11 : CSLSetNameValue(papszGeoref, "bottom_right.longitude", szValue);
769 : }
770 :
771 : /* -------------------------------------------------------------------- */
772 : /* Center */
773 : /* -------------------------------------------------------------------- */
774 11 : GDALInitGCPs(1, pasGCPList + nGCPCount);
775 11 : CPLFree(pasGCPList[nGCPCount].pszId);
776 11 : pasGCPList[nGCPCount].pszId = CPLStrdup("centre");
777 :
778 11 : temp_lat = padfTransform[3] + GetRasterXSize() * padfTransform[4] * 0.5 +
779 11 : GetRasterYSize() * padfTransform[5] * 0.5;
780 11 : temp_long = padfTransform[0] + GetRasterXSize() * padfTransform[1] * 0.5 +
781 11 : GetRasterYSize() * padfTransform[2] * 0.5;
782 11 : pasGCPList[nGCPCount].dfGCPPixel = GetRasterXSize() / 2.0;
783 11 : pasGCPList[nGCPCount].dfGCPLine = GetRasterYSize() / 2.0;
784 :
785 11 : pasGCPList[nGCPCount].dfGCPX = temp_long;
786 11 : pasGCPList[nGCPCount].dfGCPY = temp_lat;
787 11 : pasGCPList[nGCPCount].dfGCPZ = 0.0;
788 11 : nGCPCount++;
789 :
790 11 : if (poTransform != nullptr)
791 : {
792 1 : if (!bSuccess || !poTransform->Transform(1, &temp_long, &temp_lat))
793 0 : bSuccess = false;
794 : }
795 :
796 11 : if (bSuccess)
797 : {
798 11 : char szValue[128] = {'\0'};
799 11 : CPLsnprintf(szValue, sizeof(szValue), "%.10f", temp_lat);
800 11 : papszGeoref = CSLSetNameValue(papszGeoref, "centre.latitude", szValue);
801 :
802 11 : CPLsnprintf(szValue, sizeof(szValue), "%.10f", temp_long);
803 11 : papszGeoref = CSLSetNameValue(papszGeoref, "centre.longitude", szValue);
804 : }
805 :
806 11 : if (!bSuccess)
807 : {
808 0 : CPLError(CE_Warning, CPLE_AppDefined,
809 : "Error setting header info in SetGeoTransform. "
810 : "Changes may not be saved properly.");
811 : }
812 :
813 11 : if (poTransform != nullptr)
814 1 : delete poTransform;
815 :
816 11 : bGeorefChanged = true;
817 :
818 11 : return CE_None;
819 : }
820 :
821 : /************************************************************************/
822 : /* SetSpatialRef() */
823 : /* */
824 : /* We provide very limited support for setting the projection. */
825 : /************************************************************************/
826 :
827 26 : CPLErr HKVDataset::SetSpatialRef(const OGRSpatialReference *poSRS)
828 :
829 : {
830 : // Update a georef file.
831 26 : if (poSRS == nullptr)
832 : {
833 0 : m_oSRS.Clear();
834 0 : return CE_None;
835 : }
836 26 : m_oSRS = *poSRS;
837 :
838 27 : if ((m_oSRS.GetAttrValue("PROJECTION") != nullptr) &&
839 1 : (EQUAL(m_oSRS.GetAttrValue("PROJECTION"), SRS_PT_TRANSVERSE_MERCATOR)))
840 : {
841 1 : papszGeoref = CSLSetNameValue(papszGeoref, "projection.name", "utm");
842 1 : OGRErr ogrerrorOl = OGRERR_NONE;
843 1 : papszGeoref = CSLSetNameValue(
844 : papszGeoref, "projection.origin_longitude",
845 : CPLSPrintf("%f", m_oSRS.GetProjParm(SRS_PP_CENTRAL_MERIDIAN, 0.0,
846 : &ogrerrorOl)));
847 : }
848 50 : else if (m_oSRS.GetAttrValue("PROJECTION") == nullptr &&
849 25 : m_oSRS.IsGeographic())
850 : {
851 25 : papszGeoref = CSLSetNameValue(papszGeoref, "projection.name", "LL");
852 : }
853 : else
854 : {
855 0 : CPLError(CE_Warning, CPLE_AppDefined, "Unrecognized projection.");
856 0 : return CE_Failure;
857 : }
858 :
859 26 : OGRErr ogrerrorEq = OGRERR_NONE;
860 26 : const double eq_radius = m_oSRS.GetSemiMajor(&ogrerrorEq);
861 :
862 26 : OGRErr ogrerrorInvf = OGRERR_NONE;
863 26 : const double inv_flattening = m_oSRS.GetInvFlattening(&ogrerrorInvf);
864 :
865 26 : if ((ogrerrorEq == OGRERR_NONE) && (ogrerrorInvf == OGRERR_NONE))
866 : {
867 26 : HKVSpheroidList *hkvEllipsoids = new HKVSpheroidList;
868 : char *spheroid_name =
869 26 : hkvEllipsoids->GetSpheroidNameByEqRadiusAndInvFlattening(
870 : eq_radius, inv_flattening);
871 26 : if (spheroid_name != nullptr)
872 : {
873 26 : papszGeoref =
874 26 : CSLSetNameValue(papszGeoref, "spheroid.name", spheroid_name);
875 : }
876 26 : CPLFree(spheroid_name);
877 26 : delete hkvEllipsoids;
878 : }
879 : else
880 : {
881 : // Default to previous behavior if spheroid not found by radius and
882 : // inverse flattening.
883 0 : char *pszProjection = nullptr;
884 0 : m_oSRS.exportToWkt(&pszProjection);
885 0 : if (pszProjection && strstr(pszProjection, "Bessel") != nullptr)
886 : {
887 0 : papszGeoref =
888 0 : CSLSetNameValue(papszGeoref, "spheroid.name", "ev-bessel");
889 : }
890 : else
891 : {
892 0 : papszGeoref =
893 0 : CSLSetNameValue(papszGeoref, "spheroid.name", "ev-wgs-84");
894 : }
895 0 : CPLFree(pszProjection);
896 : }
897 26 : bGeorefChanged = true;
898 26 : return CE_None;
899 : }
900 :
901 : /************************************************************************/
902 : /* GetGCP() */
903 : /************************************************************************/
904 :
905 0 : const GDAL_GCP *HKVDataset::GetGCPs()
906 :
907 : {
908 0 : return pasGCPList;
909 : }
910 :
911 : /************************************************************************/
912 : /* ProcessGeorefGCP() */
913 : /************************************************************************/
914 :
915 90 : void HKVDataset::ProcessGeorefGCP(char **papszGeorefIn, const char *pszBase,
916 : double dfRasterX, double dfRasterY)
917 :
918 : {
919 : /* -------------------------------------------------------------------- */
920 : /* Fetch the GCP from the string list. */
921 : /* -------------------------------------------------------------------- */
922 90 : char szFieldName[128] = {'\0'};
923 90 : snprintf(szFieldName, sizeof(szFieldName), "%s.latitude", pszBase);
924 90 : double dfLat = 0.0;
925 90 : if (CSLFetchNameValue(papszGeorefIn, szFieldName) == nullptr)
926 75 : return;
927 : else
928 15 : dfLat = CPLAtof(CSLFetchNameValue(papszGeorefIn, szFieldName));
929 :
930 15 : snprintf(szFieldName, sizeof(szFieldName), "%s.longitude", pszBase);
931 15 : double dfLong = 0.0;
932 15 : if (CSLFetchNameValue(papszGeorefIn, szFieldName) == nullptr)
933 0 : return;
934 : else
935 15 : dfLong = CPLAtof(CSLFetchNameValue(papszGeorefIn, szFieldName));
936 :
937 : /* -------------------------------------------------------------------- */
938 : /* Add the gcp to the internal list. */
939 : /* -------------------------------------------------------------------- */
940 15 : GDALInitGCPs(1, pasGCPList + nGCPCount);
941 :
942 15 : CPLFree(pasGCPList[nGCPCount].pszId);
943 :
944 15 : pasGCPList[nGCPCount].pszId = CPLStrdup(pszBase);
945 :
946 15 : pasGCPList[nGCPCount].dfGCPX = dfLong;
947 15 : pasGCPList[nGCPCount].dfGCPY = dfLat;
948 15 : pasGCPList[nGCPCount].dfGCPZ = 0.0;
949 :
950 15 : pasGCPList[nGCPCount].dfGCPPixel = dfRasterX;
951 15 : pasGCPList[nGCPCount].dfGCPLine = dfRasterY;
952 :
953 15 : nGCPCount++;
954 : }
955 :
956 : /************************************************************************/
957 : /* ProcessGeoref() */
958 : /************************************************************************/
959 :
960 18 : void HKVDataset::ProcessGeoref(const char *pszFilename)
961 :
962 : {
963 : /* -------------------------------------------------------------------- */
964 : /* Load the georef file, and boil white space away from around */
965 : /* the equal sign. */
966 : /* -------------------------------------------------------------------- */
967 18 : CSLDestroy(papszGeoref);
968 18 : papszGeoref = CSLLoad(pszFilename);
969 18 : if (papszGeoref == nullptr)
970 0 : return;
971 :
972 18 : HKVSpheroidList *hkvEllipsoids = new HKVSpheroidList;
973 :
974 87 : for (int i = 0; papszGeoref[i] != nullptr; i++)
975 : {
976 69 : int iDst = 0;
977 69 : char *pszLine = papszGeoref[i];
978 :
979 1899 : for (int iSrc = 0; pszLine[iSrc] != '\0'; iSrc++)
980 : {
981 1830 : if (pszLine[iSrc] != ' ')
982 : {
983 1830 : pszLine[iDst++] = pszLine[iSrc];
984 : }
985 : }
986 69 : pszLine[iDst] = '\0';
987 : }
988 :
989 : /* -------------------------------------------------------------------- */
990 : /* Try to get GCPs, in lat/longs . */
991 : /* -------------------------------------------------------------------- */
992 18 : nGCPCount = 0;
993 18 : pasGCPList = reinterpret_cast<GDAL_GCP *>(CPLCalloc(sizeof(GDAL_GCP), 5));
994 :
995 18 : if (MFF2version > 1.0)
996 : {
997 18 : ProcessGeorefGCP(papszGeoref, "top_left", 0, 0);
998 18 : ProcessGeorefGCP(papszGeoref, "top_right", GetRasterXSize(), 0);
999 18 : ProcessGeorefGCP(papszGeoref, "bottom_left", 0, GetRasterYSize());
1000 18 : ProcessGeorefGCP(papszGeoref, "bottom_right", GetRasterXSize(),
1001 18 : GetRasterYSize());
1002 18 : ProcessGeorefGCP(papszGeoref, "centre", GetRasterXSize() / 2.0,
1003 18 : GetRasterYSize() / 2.0);
1004 : }
1005 : else
1006 : {
1007 0 : ProcessGeorefGCP(papszGeoref, "top_left", 0.5, 0.5);
1008 0 : ProcessGeorefGCP(papszGeoref, "top_right", GetRasterXSize() - 0.5, 0.5);
1009 0 : ProcessGeorefGCP(papszGeoref, "bottom_left", 0.5,
1010 0 : GetRasterYSize() - 0.5);
1011 0 : ProcessGeorefGCP(papszGeoref, "bottom_right", GetRasterXSize() - 0.5,
1012 0 : GetRasterYSize() - 0.5);
1013 0 : ProcessGeorefGCP(papszGeoref, "centre", GetRasterXSize() / 2.0,
1014 0 : GetRasterYSize() / 2.0);
1015 : }
1016 :
1017 18 : if (nGCPCount == 0)
1018 : {
1019 15 : CPLFree(pasGCPList);
1020 15 : pasGCPList = nullptr;
1021 : }
1022 :
1023 : /* -------------------------------------------------------------------- */
1024 : /* Do we have a recognised projection? */
1025 : /* -------------------------------------------------------------------- */
1026 18 : const char *pszProjName = CSLFetchNameValue(papszGeoref, "projection.name");
1027 : const char *pszOriginLong =
1028 18 : CSLFetchNameValue(papszGeoref, "projection.origin_longitude");
1029 : const char *pszSpheroidName =
1030 18 : CSLFetchNameValue(papszGeoref, "spheroid.name");
1031 :
1032 36 : if (pszSpheroidName != nullptr &&
1033 18 : hkvEllipsoids->SpheroidInList(pszSpheroidName))
1034 : {
1035 : #if 0
1036 : // TODO(schwehr): Enable in trunk after 2.1 branch and fix.
1037 : // Breaks tests on some platforms.
1038 : CPLError( CE_Failure, CPLE_AppDefined,
1039 : "Unrecognized ellipsoid. Not handled. "
1040 : "Spheroid name not in spheroid list: '%s'",
1041 : pszSpheroidName );
1042 : #endif
1043 : // Why were eq_radius and inv_flattening never used?
1044 : // eq_radius = hkvEllipsoids->GetSpheroidEqRadius(pszSpheroidName);
1045 : // inv_flattening =
1046 : // hkvEllipsoids->GetSpheroidInverseFlattening(pszSpheroidName);
1047 : }
1048 0 : else if (pszProjName != nullptr)
1049 : {
1050 0 : CPLError(CE_Warning, CPLE_AppDefined,
1051 : "Unrecognized ellipsoid. Not handled.");
1052 : // TODO(schwehr): This error is was never what was happening.
1053 : // CPLError( CE_Warning, CPLE_AppDefined,
1054 : // "Unrecognized ellipsoid. Using wgs-84 parameters.");
1055 : // eq_radius=hkvEllipsoids->GetSpheroidEqRadius("wgs-84"); */
1056 : // inv_flattening=hkvEllipsoids->GetSpheroidInverseFlattening("wgs-84");
1057 : }
1058 :
1059 18 : if (pszProjName != nullptr && EQUAL(pszProjName, "utm") && nGCPCount == 5)
1060 : {
1061 : // int nZone = (int)((CPLAtof(pszOriginLong)+184.5) / 6.0);
1062 3 : int nZone = 31; // TODO(schwehr): Where does 31 come from?
1063 :
1064 3 : if (pszOriginLong == nullptr)
1065 : {
1066 : // If origin not specified, assume 0.0.
1067 0 : CPLError(
1068 : CE_Warning, CPLE_AppDefined,
1069 : "No projection origin longitude specified. Assuming 0.0.");
1070 : }
1071 : else
1072 : {
1073 3 : nZone = 31 + static_cast<int>(floor(CPLAtof(pszOriginLong) / 6.0));
1074 : }
1075 :
1076 6 : OGRSpatialReference oUTM;
1077 :
1078 3 : if (pasGCPList[4].dfGCPY < 0)
1079 0 : oUTM.SetUTM(nZone, 0);
1080 : else
1081 3 : oUTM.SetUTM(nZone, 1);
1082 :
1083 6 : OGRSpatialReference oLL;
1084 3 : oLL.SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
1085 3 : if (pszOriginLong != nullptr)
1086 : {
1087 3 : oUTM.SetProjParm(SRS_PP_CENTRAL_MERIDIAN, CPLAtof(pszOriginLong));
1088 3 : oLL.SetProjParm(SRS_PP_LONGITUDE_OF_ORIGIN, CPLAtof(pszOriginLong));
1089 : }
1090 :
1091 3 : if ((pszSpheroidName == nullptr) ||
1092 3 : (EQUAL(pszSpheroidName, "wgs-84")) ||
1093 3 : (EQUAL(pszSpheroidName, "wgs_84")))
1094 : {
1095 0 : oUTM.SetWellKnownGeogCS("WGS84");
1096 0 : oLL.SetWellKnownGeogCS("WGS84");
1097 : }
1098 : else
1099 : {
1100 3 : if (hkvEllipsoids->SpheroidInList(pszSpheroidName))
1101 : {
1102 3 : oUTM.SetGeogCS(
1103 : "unknown", "unknown", pszSpheroidName,
1104 : hkvEllipsoids->GetSpheroidEqRadius(pszSpheroidName),
1105 : hkvEllipsoids->GetSpheroidInverseFlattening(
1106 : pszSpheroidName));
1107 3 : oLL.SetGeogCS(
1108 : "unknown", "unknown", pszSpheroidName,
1109 : hkvEllipsoids->GetSpheroidEqRadius(pszSpheroidName),
1110 : hkvEllipsoids->GetSpheroidInverseFlattening(
1111 : pszSpheroidName));
1112 : }
1113 : else
1114 : {
1115 0 : CPLError(CE_Warning, CPLE_AppDefined,
1116 : "Unrecognized ellipsoid. Using wgs-84 parameters.");
1117 0 : oUTM.SetWellKnownGeogCS("WGS84");
1118 0 : oLL.SetWellKnownGeogCS("WGS84");
1119 : }
1120 : }
1121 :
1122 : OGRCoordinateTransformation *poTransform =
1123 3 : OGRCreateCoordinateTransformation(&oLL, &oUTM);
1124 :
1125 3 : bool bSuccess = true;
1126 3 : if (poTransform == nullptr)
1127 : {
1128 0 : CPLErrorReset();
1129 0 : bSuccess = false;
1130 : }
1131 :
1132 3 : double dfUtmX[5] = {0.0};
1133 3 : double dfUtmY[5] = {0.0};
1134 :
1135 3 : if (poTransform != nullptr)
1136 : {
1137 18 : for (int gcp_index = 0; gcp_index < 5; gcp_index++)
1138 : {
1139 15 : dfUtmX[gcp_index] = pasGCPList[gcp_index].dfGCPX;
1140 15 : dfUtmY[gcp_index] = pasGCPList[gcp_index].dfGCPY;
1141 :
1142 15 : if (bSuccess && !poTransform->Transform(1, &(dfUtmX[gcp_index]),
1143 : &(dfUtmY[gcp_index])))
1144 0 : bSuccess = false;
1145 : }
1146 : }
1147 :
1148 3 : if (bSuccess)
1149 : {
1150 : // Update GCPS to proper projection.
1151 18 : for (int gcp_index = 0; gcp_index < 5; gcp_index++)
1152 : {
1153 15 : pasGCPList[gcp_index].dfGCPX = dfUtmX[gcp_index];
1154 15 : pasGCPList[gcp_index].dfGCPY = dfUtmY[gcp_index];
1155 : }
1156 :
1157 3 : m_oGCPSRS = oUTM;
1158 :
1159 3 : bool transform_ok = CPL_TO_BOOL(
1160 3 : GDALGCPsToGeoTransform(5, pasGCPList, adfGeoTransform, 0));
1161 :
1162 3 : if (!transform_ok)
1163 : {
1164 : // Transform may not be sufficient in all cases (slant range
1165 : // projection).
1166 0 : adfGeoTransform[0] = 0.0;
1167 0 : adfGeoTransform[1] = 1.0;
1168 0 : adfGeoTransform[2] = 0.0;
1169 0 : adfGeoTransform[3] = 0.0;
1170 0 : adfGeoTransform[4] = 0.0;
1171 0 : adfGeoTransform[5] = 1.0;
1172 0 : m_oGCPSRS.Clear();
1173 : }
1174 : else
1175 : {
1176 3 : m_oSRS = std::move(oUTM);
1177 : }
1178 : }
1179 :
1180 3 : if (poTransform != nullptr)
1181 6 : delete poTransform;
1182 : }
1183 15 : else if (pszProjName != nullptr && nGCPCount == 5)
1184 : {
1185 0 : OGRSpatialReference oLL;
1186 0 : oLL.SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
1187 :
1188 0 : if (pszOriginLong != nullptr)
1189 : {
1190 0 : oLL.SetProjParm(SRS_PP_LONGITUDE_OF_ORIGIN, CPLAtof(pszOriginLong));
1191 : }
1192 :
1193 0 : if (pszSpheroidName == nullptr ||
1194 0 : EQUAL(pszSpheroidName, "wgs-84") || // Dash.
1195 0 : EQUAL(pszSpheroidName, "wgs_84")) // Underscore.
1196 : {
1197 0 : oLL.SetWellKnownGeogCS("WGS84");
1198 : }
1199 : else
1200 : {
1201 0 : if (hkvEllipsoids->SpheroidInList(pszSpheroidName))
1202 : {
1203 0 : oLL.SetGeogCS(
1204 : "", "", pszSpheroidName,
1205 : hkvEllipsoids->GetSpheroidEqRadius(pszSpheroidName),
1206 : hkvEllipsoids->GetSpheroidInverseFlattening(
1207 : pszSpheroidName));
1208 : }
1209 : else
1210 : {
1211 0 : CPLError(CE_Warning, CPLE_AppDefined,
1212 : "Unrecognized ellipsoid. "
1213 : "Using wgs-84 parameters.");
1214 0 : oLL.SetWellKnownGeogCS("WGS84");
1215 : }
1216 : }
1217 :
1218 0 : const bool transform_ok = CPL_TO_BOOL(
1219 0 : GDALGCPsToGeoTransform(5, pasGCPList, adfGeoTransform, 0));
1220 :
1221 0 : m_oSRS.Clear();
1222 :
1223 0 : if (!transform_ok)
1224 : {
1225 0 : adfGeoTransform[0] = 0.0;
1226 0 : adfGeoTransform[1] = 1.0;
1227 0 : adfGeoTransform[2] = 0.0;
1228 0 : adfGeoTransform[3] = 0.0;
1229 0 : adfGeoTransform[4] = 0.0;
1230 0 : adfGeoTransform[5] = 1.0;
1231 : }
1232 : else
1233 : {
1234 0 : m_oSRS = oLL;
1235 : }
1236 :
1237 0 : m_oGCPSRS = std::move(oLL);
1238 : }
1239 :
1240 18 : delete hkvEllipsoids;
1241 : }
1242 :
1243 : /************************************************************************/
1244 : /* Open() */
1245 : /************************************************************************/
1246 :
1247 30698 : GDALDataset *HKVDataset::Open(GDALOpenInfo *poOpenInfo)
1248 :
1249 : {
1250 : /* -------------------------------------------------------------------- */
1251 : /* We assume the dataset is passed as a directory. Check for */
1252 : /* an attrib and blob file as a minimum. */
1253 : /* -------------------------------------------------------------------- */
1254 30698 : if (!poOpenInfo->bIsDirectory)
1255 30431 : return nullptr;
1256 :
1257 : const char *pszFilename =
1258 267 : CPLFormFilename(poOpenInfo->pszFilename, "image_data", nullptr);
1259 : VSIStatBuf sStat;
1260 267 : if (VSIStat(pszFilename, &sStat) != 0)
1261 223 : pszFilename = CPLFormFilename(poOpenInfo->pszFilename, "blob", nullptr);
1262 267 : if (VSIStat(pszFilename, &sStat) != 0)
1263 223 : return nullptr;
1264 :
1265 44 : pszFilename = CPLFormFilename(poOpenInfo->pszFilename, "attrib", nullptr);
1266 44 : if (VSIStat(pszFilename, &sStat) != 0)
1267 0 : return nullptr;
1268 :
1269 : /* -------------------------------------------------------------------- */
1270 : /* Load the attrib file, and boil white space away from around */
1271 : /* the equal sign. */
1272 : /* -------------------------------------------------------------------- */
1273 44 : char **papszAttrib = CSLLoad(pszFilename);
1274 44 : if (papszAttrib == nullptr)
1275 0 : return nullptr;
1276 :
1277 440 : for (int i = 0; papszAttrib[i] != nullptr; i++)
1278 : {
1279 396 : int iDst = 0;
1280 396 : char *pszLine = papszAttrib[i];
1281 :
1282 11173 : for (int iSrc = 0; pszLine[iSrc] != '\0'; iSrc++)
1283 : {
1284 10777 : if (pszLine[iSrc] != ' ')
1285 : {
1286 9369 : pszLine[iDst++] = pszLine[iSrc];
1287 : }
1288 : }
1289 396 : pszLine[iDst] = '\0';
1290 : }
1291 :
1292 : /* -------------------------------------------------------------------- */
1293 : /* Create a corresponding GDALDataset. */
1294 : /* -------------------------------------------------------------------- */
1295 88 : auto poDS = std::make_unique<HKVDataset>();
1296 :
1297 44 : poDS->pszPath = CPLStrdup(poOpenInfo->pszFilename);
1298 44 : poDS->papszAttrib = papszAttrib;
1299 :
1300 44 : poDS->eAccess = poOpenInfo->eAccess;
1301 :
1302 : /* -------------------------------------------------------------------- */
1303 : /* Set some dataset wide information. */
1304 : /* -------------------------------------------------------------------- */
1305 44 : bool bNative = false;
1306 44 : bool bComplex = false;
1307 44 : int nRawBands = 0;
1308 :
1309 88 : if (CSLFetchNameValue(papszAttrib, "extent.cols") == nullptr ||
1310 44 : CSLFetchNameValue(papszAttrib, "extent.rows") == nullptr)
1311 : {
1312 0 : return nullptr;
1313 : }
1314 :
1315 44 : poDS->nRasterXSize = atoi(CSLFetchNameValue(papszAttrib, "extent.cols"));
1316 44 : poDS->nRasterYSize = atoi(CSLFetchNameValue(papszAttrib, "extent.rows"));
1317 :
1318 44 : if (!GDALCheckDatasetDimensions(poDS->nRasterXSize, poDS->nRasterYSize))
1319 : {
1320 0 : return nullptr;
1321 : }
1322 :
1323 44 : const char *pszValue = CSLFetchNameValue(papszAttrib, "pixel.order");
1324 44 : if (pszValue == nullptr)
1325 0 : bNative = true;
1326 : else
1327 : {
1328 : #ifdef CPL_MSB
1329 : bNative = strstr(pszValue, "*msbf") != NULL;
1330 : #else
1331 44 : bNative = strstr(pszValue, "*lsbf") != nullptr;
1332 : #endif
1333 : }
1334 :
1335 44 : bool bNoDataSet = false;
1336 44 : double dfNoDataValue = 0.0;
1337 44 : pszValue = CSLFetchNameValue(papszAttrib, "pixel.no_data");
1338 44 : if (pszValue != nullptr)
1339 : {
1340 0 : bNoDataSet = true;
1341 0 : dfNoDataValue = CPLAtof(pszValue);
1342 : }
1343 :
1344 44 : pszValue = CSLFetchNameValue(papszAttrib, "channel.enumeration");
1345 44 : if (pszValue != nullptr)
1346 44 : nRawBands = atoi(pszValue);
1347 : else
1348 0 : nRawBands = 1;
1349 :
1350 44 : if (!GDALCheckBandCount(nRawBands, TRUE))
1351 : {
1352 0 : return nullptr;
1353 : }
1354 :
1355 44 : pszValue = CSLFetchNameValue(papszAttrib, "pixel.field");
1356 44 : if (pszValue != nullptr && strstr(pszValue, "*complex") != nullptr)
1357 10 : bComplex = true;
1358 : else
1359 34 : bComplex = false;
1360 :
1361 : /* Get the version number, if present (if not, assume old version. */
1362 : /* Versions differ in their interpretation of corner coordinates. */
1363 :
1364 44 : if (CSLFetchNameValue(papszAttrib, "version") != nullptr)
1365 44 : poDS->SetVersion(static_cast<float>(
1366 44 : CPLAtof(CSLFetchNameValue(papszAttrib, "version"))));
1367 : else
1368 0 : poDS->SetVersion(1.0);
1369 :
1370 : /* -------------------------------------------------------------------- */
1371 : /* Figure out the datatype */
1372 : /* -------------------------------------------------------------------- */
1373 44 : const char *pszEncoding = CSLFetchNameValue(papszAttrib, "pixel.encoding");
1374 44 : if (pszEncoding == nullptr)
1375 0 : pszEncoding = "{ *unsigned }";
1376 :
1377 44 : int nSize = 1;
1378 44 : if (CSLFetchNameValue(papszAttrib, "pixel.size") != nullptr)
1379 44 : nSize = atoi(CSLFetchNameValue(papszAttrib, "pixel.size")) / 8;
1380 : #if 0
1381 : int nPseudoBands;
1382 : if( bComplex )
1383 : nPseudoBands = 2;
1384 : else
1385 : nPseudoBands = 1;
1386 : #endif
1387 :
1388 : GDALDataType eType;
1389 44 : if (nSize == 1)
1390 19 : eType = GDT_Byte;
1391 25 : else if (nSize == 2 && strstr(pszEncoding, "*unsigned") != nullptr)
1392 5 : eType = GDT_UInt16;
1393 20 : else if (nSize == 4 && bComplex)
1394 5 : eType = GDT_CInt16;
1395 15 : else if (nSize == 2)
1396 5 : eType = GDT_Int16;
1397 10 : else if (nSize == 4 && strstr(pszEncoding, "*unsigned") != nullptr)
1398 0 : eType = GDT_UInt32;
1399 10 : else if (nSize == 8 && strstr(pszEncoding, "*two") != nullptr && bComplex)
1400 0 : eType = GDT_CInt32;
1401 10 : else if (nSize == 4 && strstr(pszEncoding, "*two") != nullptr)
1402 0 : eType = GDT_Int32;
1403 10 : else if (nSize == 8 && bComplex)
1404 5 : eType = GDT_CFloat32;
1405 5 : else if (nSize == 4)
1406 5 : eType = GDT_Float32;
1407 0 : else if (nSize == 16 && bComplex)
1408 0 : eType = GDT_CFloat64;
1409 0 : else if (nSize == 8)
1410 0 : eType = GDT_Float64;
1411 : else
1412 : {
1413 0 : CPLError(CE_Failure, CPLE_AppDefined,
1414 : "Unsupported pixel data type in %s.\n"
1415 : "pixel.size=%d pixel.encoding=%s",
1416 0 : poDS->pszPath, nSize, pszEncoding);
1417 0 : return nullptr;
1418 : }
1419 :
1420 : /* -------------------------------------------------------------------- */
1421 : /* Open the blob file. */
1422 : /* -------------------------------------------------------------------- */
1423 44 : pszFilename = CPLFormFilename(poDS->pszPath, "image_data", nullptr);
1424 44 : if (VSIStat(pszFilename, &sStat) != 0)
1425 0 : pszFilename = CPLFormFilename(poDS->pszPath, "blob", nullptr);
1426 44 : if (poOpenInfo->eAccess == GA_ReadOnly)
1427 : {
1428 18 : poDS->fpBlob = VSIFOpenL(pszFilename, "rb");
1429 18 : if (poDS->fpBlob == nullptr)
1430 : {
1431 0 : CPLError(CE_Failure, CPLE_OpenFailed,
1432 : "Unable to open file %s for read access.", pszFilename);
1433 0 : return nullptr;
1434 : }
1435 : }
1436 : else
1437 : {
1438 26 : poDS->fpBlob = VSIFOpenL(pszFilename, "rb+");
1439 26 : if (poDS->fpBlob == nullptr)
1440 : {
1441 0 : CPLError(CE_Failure, CPLE_OpenFailed,
1442 : "Unable to open file %s for update access.", pszFilename);
1443 0 : return nullptr;
1444 : }
1445 : }
1446 :
1447 : /* -------------------------------------------------------------------- */
1448 : /* Build the overview filename, as blob file = "_ovr". */
1449 : /* -------------------------------------------------------------------- */
1450 88 : std::string osOvrFilename(pszFilename);
1451 44 : osOvrFilename += "_ovr";
1452 :
1453 : /* -------------------------------------------------------------------- */
1454 : /* Define the bands. */
1455 : /* -------------------------------------------------------------------- */
1456 44 : const int nPixelOffset = nRawBands * nSize;
1457 44 : const int nLineOffset = nPixelOffset * poDS->GetRasterXSize();
1458 44 : int nOffset = 0;
1459 :
1460 138 : for (int iRawBand = 0; iRawBand < nRawBands; iRawBand++)
1461 : {
1462 : auto poBand = std::make_unique<HKVRasterBand>(
1463 94 : poDS.get(), poDS->GetRasterCount() + 1, poDS->fpBlob, nOffset,
1464 94 : nPixelOffset, nLineOffset, eType, bNative);
1465 94 : if (!poBand->IsValid())
1466 0 : return nullptr;
1467 :
1468 94 : if (bNoDataSet)
1469 0 : poBand->SetNoDataValue(dfNoDataValue);
1470 94 : poDS->SetBand(poDS->GetRasterCount() + 1, std::move(poBand));
1471 94 : nOffset += GDALGetDataTypeSizeBytes(eType);
1472 : }
1473 :
1474 44 : poDS->eRasterType = eType;
1475 :
1476 : /* -------------------------------------------------------------------- */
1477 : /* Process the georef file if there is one. */
1478 : /* -------------------------------------------------------------------- */
1479 44 : pszFilename = CPLFormFilename(poDS->pszPath, "georef", nullptr);
1480 44 : if (VSIStat(pszFilename, &sStat) == 0)
1481 18 : poDS->ProcessGeoref(pszFilename);
1482 :
1483 : /* -------------------------------------------------------------------- */
1484 : /* Initialize any PAM information. */
1485 : /* -------------------------------------------------------------------- */
1486 44 : poDS->SetDescription(osOvrFilename.c_str());
1487 44 : poDS->TryLoadXML();
1488 :
1489 : /* -------------------------------------------------------------------- */
1490 : /* Handle overviews. */
1491 : /* -------------------------------------------------------------------- */
1492 44 : poDS->oOvManager.Initialize(poDS.get(), osOvrFilename.c_str(), nullptr,
1493 : TRUE);
1494 :
1495 44 : return poDS.release();
1496 : }
1497 :
1498 : /************************************************************************/
1499 : /* Create() */
1500 : /************************************************************************/
1501 :
1502 51 : GDALDataset *HKVDataset::Create(const char *pszFilenameIn, int nXSize,
1503 : int nYSize, int nBandsIn, GDALDataType eType,
1504 : char ** /* papszParamList */)
1505 :
1506 : {
1507 : /* -------------------------------------------------------------------- */
1508 : /* Verify input options. */
1509 : /* -------------------------------------------------------------------- */
1510 51 : if (nBandsIn <= 0)
1511 : {
1512 1 : CPLError(CE_Failure, CPLE_NotSupported,
1513 : "HKV driver does not support %d bands.", nBandsIn);
1514 1 : return nullptr;
1515 : }
1516 :
1517 50 : if (eType != GDT_Byte && eType != GDT_UInt16 && eType != GDT_Int16 &&
1518 27 : eType != GDT_CInt16 && eType != GDT_Float32 && eType != GDT_CFloat32)
1519 : {
1520 21 : CPLError(CE_Failure, CPLE_AppDefined,
1521 : "Attempt to create HKV file with currently unsupported\n"
1522 : "data type (%s).",
1523 : GDALGetDataTypeName(eType));
1524 :
1525 21 : return nullptr;
1526 : }
1527 :
1528 : /* -------------------------------------------------------------------- */
1529 : /* Establish the name of the directory we will be creating the */
1530 : /* new HKV directory in. Verify that this is a directory. */
1531 : /* -------------------------------------------------------------------- */
1532 29 : char *pszBaseDir = nullptr;
1533 :
1534 29 : if (strlen(CPLGetPath(pszFilenameIn)) == 0)
1535 0 : pszBaseDir = CPLStrdup(".");
1536 : else
1537 29 : pszBaseDir = CPLStrdup(CPLGetPath(pszFilenameIn));
1538 :
1539 : VSIStatBuf sStat;
1540 29 : if (CPLStat(pszBaseDir, &sStat) != 0 || !VSI_ISDIR(sStat.st_mode))
1541 : {
1542 3 : CPLError(CE_Failure, CPLE_AppDefined,
1543 : "Attempt to create HKV dataset under %s,\n"
1544 : "but this is not a valid directory.",
1545 : pszBaseDir);
1546 3 : CPLFree(pszBaseDir);
1547 3 : return nullptr;
1548 : }
1549 :
1550 26 : CPLFree(pszBaseDir);
1551 26 : pszBaseDir = nullptr;
1552 :
1553 26 : if (VSIMkdir(pszFilenameIn, 0755) != 0)
1554 : {
1555 0 : CPLError(CE_Failure, CPLE_AppDefined, "Unable to create directory %s.",
1556 : pszFilenameIn);
1557 0 : return nullptr;
1558 : }
1559 :
1560 : /* -------------------------------------------------------------------- */
1561 : /* Create the header file. */
1562 : /* -------------------------------------------------------------------- */
1563 26 : CPLErr CEHeaderCreated = SaveHKVAttribFile(pszFilenameIn, nXSize, nYSize,
1564 : nBandsIn, eType, FALSE, 0.0);
1565 :
1566 26 : if (CEHeaderCreated != CE_None)
1567 0 : return nullptr;
1568 :
1569 : /* -------------------------------------------------------------------- */
1570 : /* Create the blob file. */
1571 : /* -------------------------------------------------------------------- */
1572 :
1573 : const char *pszFilename =
1574 26 : CPLFormFilename(pszFilenameIn, "image_data", nullptr);
1575 26 : FILE *fp = VSIFOpen(pszFilename, "wb");
1576 26 : if (fp == nullptr)
1577 : {
1578 0 : CPLError(CE_Failure, CPLE_OpenFailed, "Couldn't create %s.\n",
1579 : pszFilename);
1580 0 : return nullptr;
1581 : }
1582 :
1583 26 : bool bOK = VSIFWrite(reinterpret_cast<void *>(const_cast<char *>("")), 1, 1,
1584 26 : fp) == 1;
1585 26 : if (VSIFClose(fp) != 0)
1586 0 : bOK &= false;
1587 :
1588 26 : if (!bOK)
1589 0 : return nullptr;
1590 : /* -------------------------------------------------------------------- */
1591 : /* Open the dataset normally. */
1592 : /* -------------------------------------------------------------------- */
1593 26 : return GDALDataset::FromHandle(GDALOpen(pszFilenameIn, GA_Update));
1594 : }
1595 :
1596 : /************************************************************************/
1597 : /* Delete() */
1598 : /* */
1599 : /* An HKV Blob dataset consists of a bunch of files in a */
1600 : /* directory. Try to delete all the files, then the */
1601 : /* directory. */
1602 : /************************************************************************/
1603 :
1604 1 : CPLErr HKVDataset::Delete(const char *pszName)
1605 :
1606 : {
1607 : VSIStatBuf sStat;
1608 1 : if (CPLStat(pszName, &sStat) != 0 || !VSI_ISDIR(sStat.st_mode))
1609 : {
1610 0 : CPLError(CE_Failure, CPLE_AppDefined,
1611 : "%s does not appear to be an HKV Dataset, as it is not "
1612 : "a path to a directory.",
1613 : pszName);
1614 0 : return CE_Failure;
1615 : }
1616 :
1617 1 : char **papszFiles = VSIReadDir(pszName);
1618 6 : for (int i = 0; i < CSLCount(papszFiles); i++)
1619 : {
1620 5 : if (EQUAL(papszFiles[i], ".") || EQUAL(papszFiles[i], ".."))
1621 2 : continue;
1622 :
1623 : const char *pszTarget =
1624 3 : CPLFormFilename(pszName, papszFiles[i], nullptr);
1625 3 : if (VSIUnlink(pszTarget) != 0)
1626 : {
1627 0 : CPLError(CE_Failure, CPLE_AppDefined,
1628 : "Unable to delete file %s,"
1629 : "HKVDataset Delete(%s) failed.",
1630 : pszTarget, pszName);
1631 0 : CSLDestroy(papszFiles);
1632 0 : return CE_Failure;
1633 : }
1634 : }
1635 :
1636 1 : CSLDestroy(papszFiles);
1637 :
1638 1 : if (VSIRmdir(pszName) != 0)
1639 : {
1640 0 : CPLError(CE_Failure, CPLE_AppDefined,
1641 : "Unable to delete directory %s,"
1642 : "HKVDataset Delete() failed.",
1643 : pszName);
1644 0 : return CE_Failure;
1645 : }
1646 :
1647 1 : return CE_None;
1648 : }
1649 :
1650 : /************************************************************************/
1651 : /* CreateCopy() */
1652 : /************************************************************************/
1653 :
1654 20 : GDALDataset *HKVDataset::CreateCopy(const char *pszFilename,
1655 : GDALDataset *poSrcDS,
1656 : CPL_UNUSED int bStrict, char **papszOptions,
1657 : GDALProgressFunc pfnProgress,
1658 : void *pProgressData)
1659 : {
1660 20 : int nBands = poSrcDS->GetRasterCount();
1661 20 : if (nBands == 0)
1662 : {
1663 1 : CPLError(CE_Failure, CPLE_NotSupported,
1664 : "HKV driver does not support source dataset with zero band.");
1665 1 : return nullptr;
1666 : }
1667 :
1668 19 : GDALDataType eType = poSrcDS->GetRasterBand(1)->GetRasterDataType();
1669 :
1670 19 : if (!pfnProgress(0.0, nullptr, pProgressData))
1671 0 : return nullptr;
1672 :
1673 : /* check that other bands match type- sets type */
1674 : /* to unknown if they differ. */
1675 29 : for (int iBand = 1; iBand < poSrcDS->GetRasterCount(); iBand++)
1676 : {
1677 10 : GDALRasterBand *poBand = poSrcDS->GetRasterBand(iBand + 1);
1678 10 : eType = GDALDataTypeUnion(eType, poBand->GetRasterDataType());
1679 : }
1680 :
1681 19 : HKVDataset *poDS = reinterpret_cast<HKVDataset *>(Create(
1682 : pszFilename, poSrcDS->GetRasterXSize(), poSrcDS->GetRasterYSize(),
1683 : poSrcDS->GetRasterCount(), eType, papszOptions));
1684 :
1685 : /* Check that Create worked- return Null if it didn't */
1686 19 : if (poDS == nullptr)
1687 8 : return nullptr;
1688 :
1689 : /* -------------------------------------------------------------------- */
1690 : /* Copy the image data. */
1691 : /* -------------------------------------------------------------------- */
1692 11 : const int nXSize = poDS->GetRasterXSize();
1693 11 : const int nYSize = poDS->GetRasterYSize();
1694 :
1695 : int nBlockXSize, nBlockYSize;
1696 11 : poDS->GetRasterBand(1)->GetBlockSize(&nBlockXSize, &nBlockYSize);
1697 :
1698 11 : const int nBlockTotal = ((nXSize + nBlockXSize - 1) / nBlockXSize) *
1699 11 : ((nYSize + nBlockYSize - 1) / nBlockYSize) *
1700 11 : poSrcDS->GetRasterCount();
1701 :
1702 11 : int nBlocksDone = 0;
1703 32 : for (int iBand = 0; iBand < poSrcDS->GetRasterCount(); iBand++)
1704 : {
1705 21 : GDALRasterBand *poSrcBand = poSrcDS->GetRasterBand(iBand + 1);
1706 21 : GDALRasterBand *poDstBand = poDS->GetRasterBand(iBand + 1);
1707 :
1708 : /* Get nodata value, if relevant */
1709 21 : int pbSuccess = FALSE;
1710 21 : double dfSrcNoDataValue = poSrcBand->GetNoDataValue(&pbSuccess);
1711 21 : if (pbSuccess)
1712 0 : poDS->SetNoDataValue(dfSrcNoDataValue);
1713 :
1714 63 : void *pData = CPLMalloc(nBlockXSize * nBlockYSize *
1715 21 : GDALGetDataTypeSize(eType) / 8);
1716 :
1717 241 : for (int iYOffset = 0; iYOffset < nYSize; iYOffset += nBlockYSize)
1718 : {
1719 440 : for (int iXOffset = 0; iXOffset < nXSize; iXOffset += nBlockXSize)
1720 : {
1721 220 : if (!pfnProgress((nBlocksDone++) /
1722 220 : static_cast<float>(nBlockTotal),
1723 : nullptr, pProgressData))
1724 : {
1725 0 : CPLError(CE_Failure, CPLE_UserInterrupt, "User terminated");
1726 0 : delete poDS;
1727 0 : CPLFree(pData);
1728 :
1729 : GDALDriver *poHKVDriver = reinterpret_cast<GDALDriver *>(
1730 0 : GDALGetDriverByName("MFF2"));
1731 0 : poHKVDriver->Delete(pszFilename);
1732 0 : return nullptr;
1733 : }
1734 :
1735 220 : const int nTBXSize = std::min(nBlockXSize, nXSize - iXOffset);
1736 220 : const int nTBYSize = std::min(nBlockYSize, nYSize - iYOffset);
1737 :
1738 220 : CPLErr eErr = poSrcBand->RasterIO(
1739 : GF_Read, iXOffset, iYOffset, nTBXSize, nTBYSize, pData,
1740 : nTBXSize, nTBYSize, eType, 0, 0, nullptr);
1741 220 : if (eErr != CE_None)
1742 : {
1743 0 : delete poDS;
1744 0 : CPLFree(pData);
1745 0 : return nullptr;
1746 : }
1747 :
1748 220 : eErr = poDstBand->RasterIO(GF_Write, iXOffset, iYOffset,
1749 : nTBXSize, nTBYSize, pData, nTBXSize,
1750 : nTBYSize, eType, 0, 0, nullptr);
1751 :
1752 220 : if (eErr != CE_None)
1753 : {
1754 0 : delete poDS;
1755 0 : CPLFree(pData);
1756 0 : return nullptr;
1757 : }
1758 : }
1759 : }
1760 :
1761 21 : CPLFree(pData);
1762 : }
1763 :
1764 : /* -------------------------------------------------------------------- */
1765 : /* Copy georeferencing information, if enough is available. */
1766 : /* Only copy geotransform-style info (won't work for slant range). */
1767 : /* -------------------------------------------------------------------- */
1768 :
1769 : double *tempGeoTransform =
1770 11 : static_cast<double *>(CPLMalloc(6 * sizeof(double)));
1771 :
1772 22 : if ((poSrcDS->GetGeoTransform(tempGeoTransform) == CE_None) &&
1773 11 : (tempGeoTransform[0] != 0.0 || tempGeoTransform[1] != 1.0 ||
1774 0 : tempGeoTransform[2] != 0.0 || tempGeoTransform[3] != 0.0 ||
1775 0 : tempGeoTransform[4] != 0.0 || std::abs(tempGeoTransform[5]) != 1.0))
1776 : {
1777 11 : auto poSrcSRS = poSrcDS->GetSpatialRef();
1778 11 : if (poSrcSRS)
1779 : {
1780 11 : poDS->SetSpatialRef(poSrcSRS);
1781 11 : poDS->m_oGCPSRS = *poSrcSRS;
1782 : }
1783 11 : poDS->SetGeoTransform(tempGeoTransform);
1784 :
1785 11 : CPLFree(tempGeoTransform);
1786 :
1787 : // georef file will be saved automatically when dataset is deleted
1788 : // because SetProjection sets a flag to indicate it is necessary.
1789 : }
1790 : else
1791 : {
1792 0 : CPLFree(tempGeoTransform);
1793 : }
1794 :
1795 : // Make sure image data gets flushed.
1796 32 : for (int iBand = 0; iBand < poDS->GetRasterCount(); iBand++)
1797 : {
1798 : RawRasterBand *poDstBand =
1799 21 : reinterpret_cast<RawRasterBand *>(poDS->GetRasterBand(iBand + 1));
1800 21 : poDstBand->FlushCache(false);
1801 : }
1802 :
1803 11 : if (!pfnProgress(1.0, nullptr, pProgressData))
1804 : {
1805 0 : CPLError(CE_Failure, CPLE_UserInterrupt, "User terminated");
1806 0 : delete poDS;
1807 :
1808 : GDALDriver *poHKVDriver =
1809 0 : reinterpret_cast<GDALDriver *>(GDALGetDriverByName("MFF2"));
1810 0 : poHKVDriver->Delete(pszFilename);
1811 0 : return nullptr;
1812 : }
1813 :
1814 11 : poDS->CloneInfo(poSrcDS, GCIF_PAM_DEFAULT);
1815 :
1816 11 : return poDS;
1817 : }
1818 :
1819 : /************************************************************************/
1820 : /* GDALRegister_HKV() */
1821 : /************************************************************************/
1822 :
1823 1595 : void GDALRegister_HKV()
1824 :
1825 : {
1826 1595 : if (GDALGetDriverByName("MFF2") != nullptr)
1827 302 : return;
1828 :
1829 1293 : GDALDriver *poDriver = new GDALDriver();
1830 :
1831 1293 : poDriver->SetDescription("MFF2");
1832 1293 : poDriver->SetMetadataItem(GDAL_DCAP_RASTER, "YES");
1833 1293 : poDriver->SetMetadataItem(GDAL_DMD_LONGNAME, "Vexcel MFF2 (HKV) Raster");
1834 1293 : poDriver->SetMetadataItem(GDAL_DMD_HELPTOPIC, "drivers/raster/mff2.html");
1835 1293 : poDriver->SetMetadataItem(GDAL_DMD_CREATIONDATATYPES,
1836 : "Byte Int16 UInt16 Int32 UInt32 CInt16 "
1837 1293 : "CInt32 Float32 Float64 CFloat32 CFloat64");
1838 :
1839 1293 : poDriver->pfnOpen = HKVDataset::Open;
1840 1293 : poDriver->pfnCreate = HKVDataset::Create;
1841 1293 : poDriver->pfnDelete = HKVDataset::Delete;
1842 1293 : poDriver->pfnCreateCopy = HKVDataset::CreateCopy;
1843 :
1844 1293 : GetGDALDriverManager()->RegisterDriver(poDriver);
1845 : }
|
