Line data Source code
1 : /******************************************************************************
2 : *
3 : * Project: S-101 driver
4 : * Purpose: Implements OGRS101Reader
5 : * Author: Even Rouault <even dot rouault at spatialys.com>
6 : *
7 : ******************************************************************************
8 : * Copyright (c) 2026, Even Rouault <even dot rouault at spatialys.com>
9 : *
10 : * SPDX-License-Identifier: MIT
11 : ****************************************************************************/
12 :
13 : #include "ogr_s101.h"
14 : #include "ogrs101readerconstants.h"
15 :
16 : #include <memory>
17 :
18 : /************************************************************************/
19 : /* GetPointLayerName() */
20 : /************************************************************************/
21 :
22 : /* static */ std::string
23 501 : OGRS101Reader::GetPointLayerName(const OGRSpatialReference &oSRS)
24 : {
25 501 : return oSRS.GetAxesCount() == 2
26 : ? "Point2D"
27 1503 : : CPLSPrintf("Point3D_%s", LaunderCRSName(oSRS).c_str());
28 : }
29 :
30 : /************************************************************************/
31 : /* CreatePointFeatureDefns() */
32 : /************************************************************************/
33 :
34 : /** Create the feature definition(s) for the Point layer(s)
35 : *
36 : * There is a layer per CRS used by points.
37 : */
38 257 : bool OGRS101Reader::CreatePointFeatureDefns()
39 : {
40 257 : bool bError = false;
41 257 : m_oMapCRSIdToPointRecordIdx = CreateMapCRSIdToRecordIdxForPoints(bError);
42 257 : if (bError)
43 2 : return false;
44 420 : for (const auto &[nCRSId, anRecordIdx] : m_oMapCRSIdToPointRecordIdx)
45 : {
46 169 : const auto &oSRS = m_oMapSRS[nCRSId];
47 169 : const bool bIs2D = nCRSId == HORIZONTAL_CRS_ID;
48 : auto poFDefn = OGRFeatureDefnRefCountedPtr::makeInstance(
49 169 : GetPointLayerName(oSRS).c_str());
50 169 : poFDefn->SetGeomType(bIs2D ? wkbPoint : wkbPoint25D);
51 338 : poFDefn->GetGeomFieldDefn(0)->SetSpatialRef(
52 338 : OGRSpatialReferenceRefCountedPtr::makeClone(&oSRS).get());
53 : {
54 338 : OGRFieldDefn oFieldDefn(OGR_FIELD_NAME_RECORD_ID, OFTInteger);
55 169 : poFDefn->AddFieldDefn(&oFieldDefn);
56 : }
57 : {
58 338 : OGRFieldDefn oFieldDefn(OGR_FIELD_NAME_RECORD_VERSION, OFTInteger);
59 169 : poFDefn->AddFieldDefn(&oFieldDefn);
60 : }
61 :
62 169 : if (!InferFeatureDefn(m_oPointRecordIndex, PRID_FIELD, INAS_FIELD,
63 169 : anRecordIdx, *poFDefn, m_oMapFieldDomains))
64 : {
65 4 : return false;
66 : }
67 165 : m_oMapPointFeatureDefn[nCRSId] = std::move(poFDefn);
68 : }
69 :
70 251 : return true;
71 : }
72 :
73 : /************************************************************************/
74 : /* GetCRSIdForPointRecord() */
75 : /************************************************************************/
76 :
77 : /** Return the CRS id for a given Point record, or INVALID_CRS_ID on error */
78 : OGRS101Reader::CRSId
79 398 : OGRS101Reader::GetCRSIdForPointRecord(const DDFRecord *poRecord, int iRecord,
80 : int nRecordID) const
81 : {
82 398 : if (nRecordID < 0)
83 302 : nRecordID = poRecord->GetIntSubfield(PRID_FIELD, 0, RCID_SUBFIELD, 0);
84 :
85 14 : const auto GetErrorContext = [iRecord, nRecordID]()
86 : {
87 7 : if (iRecord >= 0)
88 7 : return CPLSPrintf("Record index=%d of PRID", iRecord);
89 : else
90 0 : return CPLSPrintf("Record ID=%d of PRID", nRecordID);
91 398 : };
92 :
93 398 : if (poRecord->FindField(C3IT_FIELD))
94 : {
95 : const CRSId nVCID =
96 40 : poRecord->GetIntSubfield(C3IT_FIELD, 0, VCID_SUBFIELD, 0);
97 40 : if (nVCID == HORIZONTAL_CRS_ID)
98 : {
99 2 : CPL_IGNORE_RET_VAL(EMIT_ERROR_OR_WARNING(
100 : CPLSPrintf("%s: VCID subfield = %d of C3IT "
101 : "field points to a non-3D CRS.",
102 : GetErrorContext(), static_cast<int>(nVCID))));
103 2 : return INVALID_CRS_ID;
104 : }
105 38 : else if (!cpl::contains(m_oMapSRS, nVCID))
106 : {
107 2 : CPL_IGNORE_RET_VAL(EMIT_ERROR_OR_WARNING(
108 : CPLSPrintf("%s: Unknown value %d for VCID subfield of C3IT "
109 : "field.",
110 : GetErrorContext(), static_cast<int>(nVCID))));
111 2 : return INVALID_CRS_ID;
112 : }
113 : else
114 : {
115 36 : return nVCID;
116 : }
117 : }
118 358 : else if (poRecord->FindField(C2IT_FIELD))
119 : {
120 355 : return HORIZONTAL_CRS_ID;
121 : }
122 : else
123 : {
124 3 : CPL_IGNORE_RET_VAL(EMIT_ERROR_OR_WARNING(
125 : CPLSPrintf("%s: No C2IT or C3IT field found.", GetErrorContext())));
126 3 : return INVALID_CRS_ID;
127 : }
128 : }
129 :
130 : /************************************************************************/
131 : /* CreateMapCRSIdToRecordIdxForPoints() */
132 : /************************************************************************/
133 :
134 : /** Browse through m_oPointRecordIndex to identify which record belongs to
135 : * each CRS and create a map from each CRS id to the record indices that use
136 : * it.
137 : */
138 : std::map<OGRS101Reader::CRSId, std::vector<int>>
139 257 : OGRS101Reader::CreateMapCRSIdToRecordIdxForPoints(bool &bError) const
140 : {
141 514 : std::map<OGRS101Reader::CRSId, std::vector<int>> map;
142 :
143 257 : const int nRecords = m_oPointRecordIndex.GetCount();
144 557 : for (int iRecord = 0; iRecord < nRecords; ++iRecord)
145 : {
146 302 : const auto poRecord = m_oPointRecordIndex.GetByIndex(iRecord);
147 : const CRSId nCRSId =
148 302 : GetCRSIdForPointRecord(poRecord, iRecord, /* nRecordID = */ -1);
149 302 : if (nCRSId == INVALID_CRS_ID)
150 : {
151 7 : if (m_bStrict)
152 : {
153 2 : bError = true;
154 2 : return {};
155 : }
156 : }
157 : else
158 : {
159 295 : map[nCRSId].push_back(iRecord);
160 : }
161 : }
162 :
163 255 : return map;
164 : }
165 :
166 : /************************************************************************/
167 : /* ReadPointGeometry() */
168 : /************************************************************************/
169 :
170 10615 : std::unique_ptr<OGRPoint> OGRS101Reader::ReadPointGeometryInternal(
171 : const DDFRecord *poRecord, int iRecord, int nRecordID, int iPnt,
172 : const OGRSpatialReference *poSRS, const bool bIs3D,
173 : const DDFField *poCoordField, const char *pszRecordFieldName) const
174 : {
175 58 : const auto GetErrorContext = [iRecord, nRecordID, pszRecordFieldName]()
176 : {
177 29 : if (iRecord >= 0)
178 22 : return CPLSPrintf("Record index=%d of %s", iRecord,
179 22 : pszRecordFieldName);
180 : else
181 7 : return CPLSPrintf("Record ID=%d of %s", nRecordID,
182 7 : pszRecordFieldName);
183 10615 : };
184 :
185 10615 : if (!poCoordField)
186 : {
187 0 : CPL_IGNORE_RET_VAL(EMIT_ERROR_OR_WARNING(CPLSPrintf(
188 : "%s: cannot find coordinate field.", GetErrorContext())));
189 0 : return nullptr;
190 : }
191 :
192 10615 : int bSuccess = false;
193 : const double dfX =
194 10615 : poRecord->GetIntSubfield(poCoordField, XCOO_SUBFIELD, iPnt, &bSuccess) /
195 10615 : static_cast<double>(m_nXScale) +
196 10615 : m_dfXShift;
197 10615 : if (!bSuccess &&
198 0 : !EMIT_ERROR_OR_WARNING(CPLSPrintf("%s: cannot read %s.",
199 : GetErrorContext(), XCOO_SUBFIELD)))
200 : {
201 0 : return nullptr;
202 : }
203 : const double dfY =
204 10615 : poRecord->GetIntSubfield(poCoordField, YCOO_SUBFIELD, iPnt, &bSuccess) /
205 10615 : static_cast<double>(m_nYScale) +
206 10615 : m_dfYShift;
207 10615 : if (!bSuccess &&
208 0 : !EMIT_ERROR_OR_WARNING(CPLSPrintf("%s: cannot read %s.",
209 : GetErrorContext(), YCOO_SUBFIELD)))
210 : {
211 0 : return nullptr;
212 : }
213 10644 : if (poSRS && !(std::fabs(dfX) <= 180 && std::fabs(dfY) <= 90) &&
214 29 : poSRS->IsGeographic())
215 : {
216 29 : if (!EMIT_ERROR_OR_WARNING(
217 : CPLSPrintf("%s: wrong coordinate value: lon=%f, lat=%f.",
218 : GetErrorContext(), dfX, dfY)))
219 : {
220 8 : return nullptr;
221 : }
222 : }
223 10607 : std::unique_ptr<OGRPoint> poPoint;
224 10607 : if (bIs3D)
225 : {
226 264 : const double dfZ = poRecord->GetIntSubfield(poCoordField, ZCOO_SUBFIELD,
227 264 : iPnt, &bSuccess) /
228 264 : static_cast<double>(m_nZScale) +
229 264 : m_dfZShift;
230 264 : if (!bSuccess &&
231 0 : !EMIT_ERROR_OR_WARNING(CPLSPrintf(
232 : "%s: cannot read %s.", GetErrorContext(), ZCOO_SUBFIELD)))
233 : {
234 0 : return nullptr;
235 : }
236 264 : poPoint = std::make_unique<OGRPoint>(dfX, dfY, dfZ);
237 : }
238 : else
239 : {
240 10343 : poPoint = std::make_unique<OGRPoint>(dfX, dfY);
241 : }
242 10607 : poPoint->assignSpatialReference(poSRS);
243 10607 : return poPoint;
244 : }
245 :
246 : std::unique_ptr<OGRPoint>
247 992 : OGRS101Reader::ReadPointGeometry(const DDFRecord *poRecord, int iRecord,
248 : int nRecordID,
249 : const OGRSpatialReference *poSRS) const
250 : {
251 992 : const DDFField *poCoordField = poRecord->FindField(C3IT_FIELD);
252 992 : const bool bIs3D = poCoordField != nullptr;
253 992 : if (!poCoordField)
254 855 : poCoordField = poRecord->FindField(C2IT_FIELD);
255 :
256 : return ReadPointGeometryInternal(poRecord, iRecord, nRecordID,
257 : /* iPnt = */ 0, poSRS, bIs3D, poCoordField,
258 992 : PRID_FIELD);
259 : }
260 :
261 : /************************************************************************/
262 : /* FillFeaturePoint() */
263 : /************************************************************************/
264 :
265 : /** Fill the content of the provided feature from the identified record
266 : * (of m_oPointRecordIndex).
267 : */
268 691 : bool OGRS101Reader::FillFeaturePoint(const DDFRecordIndex &oIndex, int iRecord,
269 : OGRFeature &oFeature) const
270 : {
271 691 : const auto poRecord = oIndex.GetByIndex(iRecord);
272 691 : CPLAssert(poRecord);
273 :
274 : const OGRSpatialReference *poSRS =
275 691 : oFeature.GetDefnRef()->GetGeomFieldDefn(0)->GetSpatialRef();
276 : auto poPoint =
277 1382 : ReadPointGeometry(poRecord, iRecord, /* nRecordID = */ -1, poSRS);
278 691 : if (poPoint)
279 686 : oFeature.SetGeometry(std::move(poPoint));
280 5 : else if (m_bStrict)
281 5 : return false;
282 :
283 1372 : return FillFeatureAttributes(oIndex, iRecord, INAS_FIELD, oFeature) &&
284 686 : FillFeatureWithNonAttrAssocSubfields(poRecord, iRecord, INAS_FIELD,
285 686 : oFeature);
286 : }
|
