GDAL

gdalgeoloc_dataset_accessor.h
1/******************************************************************************
2 *
3 * Project: GDAL
4 * Purpose: Dataset storage of geolocation array and backmap
5 * Author: Even Rouault, <even.rouault at spatialys.com>
6 *
7 ******************************************************************************
8 * Copyright (c) 2022, Planet Labs
9 *
10 * SPDX-License-Identifier: MIT
11 ****************************************************************************/
12
13#include "gdalcachedpixelaccessor.h"
14
17/************************************************************************/
18/* GDALGeoLocDatasetAccessors */
19/************************************************************************/
20
21 class GDALGeoLocDatasetAccessors
22{
23 typedef class GDALGeoLocDatasetAccessors AccessorType;
24
25 GDALGeoLocTransformInfo *m_psTransform;
26
27 CPLStringList m_aosGTiffCreationOptions{};
28
29 GDALDataset *m_poGeolocTmpDataset = nullptr;
30 GDALDataset *m_poBackmapTmpDataset = nullptr;
31 GDALDataset *m_poBackmapWeightsTmpDataset = nullptr;
32
33 GDALGeoLocDatasetAccessors(const GDALGeoLocDatasetAccessors &) = delete;
34 GDALGeoLocDatasetAccessors &
35 operator=(const GDALGeoLocDatasetAccessors &) = delete;
36
37 bool LoadGeoloc(bool bIsRegularGrid);
38
39 public:
40 static constexpr int TILE_SIZE = 256;
41 static constexpr int TILE_COUNT = 64;
42
43 GDALCachedPixelAccessor<double, TILE_SIZE, TILE_COUNT> geolocXAccessor;
44 GDALCachedPixelAccessor<double, TILE_SIZE, TILE_COUNT> geolocYAccessor;
45 GDALCachedPixelAccessor<float, TILE_SIZE, TILE_COUNT> backMapXAccessor;
46 GDALCachedPixelAccessor<float, TILE_SIZE, TILE_COUNT> backMapYAccessor;
47 GDALCachedPixelAccessor<float, TILE_SIZE, TILE_COUNT> backMapWeightAccessor;
48
49 explicit GDALGeoLocDatasetAccessors(GDALGeoLocTransformInfo *psTransform)
50 : m_psTransform(psTransform), geolocXAccessor(nullptr),
51 geolocYAccessor(nullptr), backMapXAccessor(nullptr),
52 backMapYAccessor(nullptr), backMapWeightAccessor(nullptr)
53 {
54 m_aosGTiffCreationOptions.SetNameValue("TILED", "YES");
55 m_aosGTiffCreationOptions.SetNameValue("INTERLEAVE", "BAND");
56 m_aosGTiffCreationOptions.SetNameValue("BLOCKXSIZE",
57 CPLSPrintf("%d", TILE_SIZE));
58 m_aosGTiffCreationOptions.SetNameValue("BLOCKYSIZE",
59 CPLSPrintf("%d", TILE_SIZE));
60 }
61
62 ~GDALGeoLocDatasetAccessors();
63
64 bool Load(bool bIsRegularGrid, bool bUseQuadtree);
65
66 bool AllocateBackMap();
67
68 GDALDataset *GetBackmapDataset();
69 void FlushBackmapCaches();
70
71 static void ReleaseBackmapDataset(GDALDataset *)
72 {
73 }
74
75 void FreeWghtsBackMap();
76};
77
78/************************************************************************/
79/* ~GDALGeoLocDatasetAccessors() */
80/************************************************************************/
81
82GDALGeoLocDatasetAccessors::~GDALGeoLocDatasetAccessors()
83{
84 geolocXAccessor.ResetModifiedFlag();
85 geolocYAccessor.ResetModifiedFlag();
86 backMapXAccessor.ResetModifiedFlag();
87 backMapYAccessor.ResetModifiedFlag();
88
89 FreeWghtsBackMap();
90
91 delete m_poGeolocTmpDataset;
92 delete m_poBackmapTmpDataset;
93}
94
95/************************************************************************/
96/* AllocateBackMap() */
97/************************************************************************/
98
99 bool GDALGeoLocDatasetAccessors::AllocateBackMap()
100{
101 auto poDriver = GDALDriver::FromHandle(GDALGetDriverByName("GTiff"));
102 if (poDriver == nullptr)
103 return false;
104
105 // CPLResetExtension / CPLGenerateTempFilename generate short-lived strings,
106 // so store them in a long-lived std::string
107 const std::string osBackmapTmpFilename =
108 CPLResetExtension(CPLGenerateTempFilename(nullptr), "tif");
109 m_poBackmapTmpDataset = poDriver->Create(
110 osBackmapTmpFilename.c_str(), m_psTransform->nBackMapWidth,
111 m_psTransform->nBackMapHeight, 2, GDT_Float32,
112 m_aosGTiffCreationOptions.List());
113 if (m_poBackmapTmpDataset == nullptr)
114 {
115 return false;
116 }
117 m_poBackmapTmpDataset->MarkSuppressOnClose();
118 VSIUnlink(m_poBackmapTmpDataset->GetDescription());
119 auto poBandX = m_poBackmapTmpDataset->GetRasterBand(1);
120 auto poBandY = m_poBackmapTmpDataset->GetRasterBand(2);
121
122 backMapXAccessor.SetBand(poBandX);
123 backMapYAccessor.SetBand(poBandY);
124
125 // CPLResetExtension / CPLGenerateTempFilename generate short-lived strings,
126 // so store them in a long-lived std::string
127 const std::string osBackmapWeightsTmpFilename =
128 CPLResetExtension(CPLGenerateTempFilename(nullptr), "tif");
129 m_poBackmapWeightsTmpDataset = poDriver->Create(
130 osBackmapWeightsTmpFilename.c_str(), m_psTransform->nBackMapWidth,
131 m_psTransform->nBackMapHeight, 1, GDT_Float32,
132 m_aosGTiffCreationOptions.List());
133 if (m_poBackmapWeightsTmpDataset == nullptr)
134 {
135 return false;
136 }
137 m_poBackmapWeightsTmpDataset->MarkSuppressOnClose();
138 VSIUnlink(m_poBackmapWeightsTmpDataset->GetDescription());
139 backMapWeightAccessor.SetBand(
140 m_poBackmapWeightsTmpDataset->GetRasterBand(1));
141
142 return true;
143}
144
145/************************************************************************/
146/* FreeWghtsBackMap() */
147/************************************************************************/
148
149 void GDALGeoLocDatasetAccessors::FreeWghtsBackMap()
150{
151 if (m_poBackmapWeightsTmpDataset)
152 {
153 backMapWeightAccessor.ResetModifiedFlag();
154 delete m_poBackmapWeightsTmpDataset;
155 m_poBackmapWeightsTmpDataset = nullptr;
156 }
157}
158
159/************************************************************************/
160/* GetBackmapDataset() */
161/************************************************************************/
162
163 GDALDataset *GDALGeoLocDatasetAccessors::GetBackmapDataset()
164{
165 auto poBandX = m_poBackmapTmpDataset->GetRasterBand(1);
166 auto poBandY = m_poBackmapTmpDataset->GetRasterBand(2);
167 poBandX->SetNoDataValue(INVALID_BMXY);
168 poBandY->SetNoDataValue(INVALID_BMXY);
169 return m_poBackmapTmpDataset;
170}
171
172/************************************************************************/
173/* FlushBackmapCaches() */
174/************************************************************************/
175
176 void GDALGeoLocDatasetAccessors::FlushBackmapCaches()
177{
178 backMapXAccessor.FlushCache();
179 backMapYAccessor.FlushCache();
180}
181
182/************************************************************************/
183/* Load() */
184/************************************************************************/
185
186 bool GDALGeoLocDatasetAccessors::Load(bool bIsRegularGrid, bool bUseQuadtree)
187{
188 return LoadGeoloc(bIsRegularGrid) &&
189 ((bUseQuadtree && GDALGeoLocBuildQuadTree(m_psTransform)) ||
190 (!bUseQuadtree &&
191 GDALGeoLoc<AccessorType>::GenerateBackMap(m_psTransform)));
192}
193
194/************************************************************************/
195/* LoadGeoloc() */
196/************************************************************************/
197
198 bool GDALGeoLocDatasetAccessors::LoadGeoloc(bool bIsRegularGrid)
199
200{
201 if (bIsRegularGrid)
202 {
203 const int nXSize = m_psTransform->nGeoLocXSize;
204 const int nYSize = m_psTransform->nGeoLocYSize;
205
206 auto poDriver = GDALDriver::FromHandle(GDALGetDriverByName("GTiff"));
207 if (poDriver == nullptr)
208 return false;
209
210 // CPLResetExtension / CPLGenerateTempFilename generate short-lived
211 // strings, so store them in a long-lived std::string
212 const std::string osGeolocTmpFilename =
213 CPLResetExtension(CPLGenerateTempFilename(nullptr), "tif");
214 m_poGeolocTmpDataset =
215 poDriver->Create(osGeolocTmpFilename.c_str(), nXSize, nYSize, 2,
216 GDT_Float64, m_aosGTiffCreationOptions.List());
217 if (m_poGeolocTmpDataset == nullptr)
218 {
219 return false;
220 }
221 m_poGeolocTmpDataset->MarkSuppressOnClose();
222 VSIUnlink(m_poGeolocTmpDataset->GetDescription());
223
224 auto poXBand = m_poGeolocTmpDataset->GetRasterBand(1);
225 auto poYBand = m_poGeolocTmpDataset->GetRasterBand(2);
226
227 // Case of regular grid.
228 // The XBAND contains the x coordinates for all lines.
229 // The YBAND contains the y coordinates for all columns.
230
231 double *padfTempX =
232 static_cast<double *>(VSI_MALLOC2_VERBOSE(nXSize, sizeof(double)));
233 double *padfTempY =
234 static_cast<double *>(VSI_MALLOC2_VERBOSE(nYSize, sizeof(double)));
235 if (padfTempX == nullptr || padfTempY == nullptr)
236 {
237 CPLFree(padfTempX);
238 CPLFree(padfTempY);
239 return false;
240 }
241
242 CPLErr eErr =
243 GDALRasterIO(m_psTransform->hBand_X, GF_Read, 0, 0, nXSize, 1,
244 padfTempX, nXSize, 1, GDT_Float64, 0, 0);
245
246 for (int j = 0; j < nYSize; j++)
247 {
248 if (poXBand->RasterIO(GF_Write, 0, j, nXSize, 1, padfTempX, nXSize,
249 1, GDT_Float64, 0, 0, nullptr) != CE_None)
250 {
251 eErr = CE_Failure;
252 break;
253 }
254 }
255
256 if (eErr == CE_None)
257 {
258 eErr = GDALRasterIO(m_psTransform->hBand_Y, GF_Read, 0, 0, nYSize,
259 1, padfTempY, nYSize, 1, GDT_Float64, 0, 0);
260
261 for (int i = 0; i < nXSize; i++)
262 {
263 if (poYBand->RasterIO(GF_Write, i, 0, 1, nYSize, padfTempY, 1,
264 nYSize, GDT_Float64, 0, 0,
265 nullptr) != CE_None)
266 {
267 eErr = CE_Failure;
268 break;
269 }
270 }
271 }
272
273 CPLFree(padfTempX);
274 CPLFree(padfTempY);
275
276 if (eErr != CE_None)
277 return false;
278
279 geolocXAccessor.SetBand(poXBand);
280 geolocYAccessor.SetBand(poYBand);
281 }
282 else
283 {
284 geolocXAccessor.SetBand(
285 GDALRasterBand::FromHandle(m_psTransform->hBand_X));
286 geolocYAccessor.SetBand(
287 GDALRasterBand::FromHandle(m_psTransform->hBand_Y));
288 }
289
290 GDALGeoLoc<GDALGeoLocDatasetAccessors>::LoadGeolocFinish(m_psTransform);
291 return true;
292}
293
CPLStringList
String list class designed around our use of C "char**" string lists.
Definition: cpl_string.h:436
GDALCachedPixelAccessor
Class to have reasonably fast random pixel access to a raster band, when accessing multiple pixels th...
Definition: gdalcachedpixelaccessor.h:37
GDALCachedPixelAccessor::FlushCache
bool FlushCache()
Flush content of modified tiles and drop caches.
Definition: gdalcachedpixelaccessor.h:265
GDALCachedPixelAccessor::SetBand
void SetBand(GDALRasterBand *poBand)
Assign the raster band if not known at construction time.
Definition: gdalcachedpixelaccessor.h:68
GDALCachedPixelAccessor::ResetModifiedFlag
void ResetModifiedFlag()
Reset the modified flag for cached tiles.
Definition: gdalcachedpixelaccessor.h:286
GDALDataset
A set of associated raster bands, usually from one file.
Definition: gdal_priv.h:495
GDALDataset::MarkSuppressOnClose
void MarkSuppressOnClose()
Set that the dataset must be deleted on close.
Definition: gdaldataset.cpp:1727
GDALDataset::GetRasterBand
GDALRasterBand * GetRasterBand(int)
Fetch a band object for a dataset.
Definition: gdaldataset.cpp:1012
GDALDriver::FromHandle
static GDALDriver * FromHandle(GDALDriverH hDriver)
Convert a GDALDriverH to a GDALDriver*.
Definition: gdal_priv.h:2348
GDALMajorObject::GetDescription
virtual const char * GetDescription() const
Fetch object description.
Definition: gdalmajorobject.cpp:61
GDALRasterBand::SetNoDataValue
virtual CPLErr SetNoDataValue(double dfNoData)
Set the no data value for this band.
Definition: gdalrasterband.cpp:2500
GDALRasterBand::FromHandle
static GDALRasterBand * FromHandle(GDALRasterBandH hBand)
Convert a GDALRasterBandH to a GDALRasterBand*.
Definition: gdal_priv.h:1890
CPLFree
#define CPLFree
Alias of VSIFree()
Definition: cpl_conv.h:82
CPLGenerateTempFilename
const char * CPLGenerateTempFilename(const char *pszStem)
Generate temporary file name.
Definition: cpl_path.cpp:1094
CPLResetExtension
const char * CPLResetExtension(const char *, const char *)
Replace the extension with the provided one.
Definition: cpl_path.cpp:450
CPLErr
CPLErr
Error category.
Definition: cpl_error.h:37
CPLSPrintf
const char * CPLSPrintf(const char *fmt,...)
CPLSPrintf() that works with 10 static buffer.
Definition: cpl_string.cpp:967
VSIUnlink
int VSIUnlink(const char *pszFilename)
Delete a file.
Definition: cpl_vsil.cpp:434
VSI_MALLOC2_VERBOSE
#define VSI_MALLOC2_VERBOSE(nSize1, nSize2)
VSI_MALLOC2_VERBOSE.
Definition: cpl_vsi.h:339
GDT_Float64
@ GDT_Float64
Definition: gdal.h:59
GDT_Float32
@ GDT_Float32
Definition: gdal.h:58
GF_Write
@ GF_Write
Definition: gdal.h:119
GF_Read
@ GF_Read
Definition: gdal.h:118
GDALGetDriverByName
GDALDriverH GDALGetDriverByName(const char *)
Fetch a driver based on the short name.
Definition: gdaldrivermanager.cpp:696
GDALRasterIO
CPLErr GDALRasterIO(GDALRasterBandH hRBand, GDALRWFlag eRWFlag, int nDSXOff, int nDSYOff, int nDSXSize, int nDSYSize, void *pBuffer, int nBXSize, int nBYSize, GDALDataType eBDataType, int nPixelSpace, int nLineSpace)
Read/write a region of image data for this band.
Definition: gdalrasterband.cpp:457