Line data Source code
1 : /******************************************************************************
2 : *
3 : * Project: WCS Client Driver
4 : * Purpose: Implementation of Dataset class for WCS 2.0.
5 : * Author: Ari Jolma <ari dot jolma at gmail dot com>
6 : *
7 : ******************************************************************************
8 : * Copyright (c) 2006, Frank Warmerdam
9 : * Copyright (c) 2008-2013, Even Rouault <even dot rouault at spatialys.com>
10 : * Copyright (c) 2017, Ari Jolma
11 : * Copyright (c) 2017, Finnish Environment Institute
12 : *
13 : * SPDX-License-Identifier: MIT
14 : ****************************************************************************/
15 :
16 : #include "cpl_string.h"
17 : #include "cpl_minixml.h"
18 : #include "cpl_http.h"
19 : #include "cpl_conv.h"
20 : #include "gmlutils.h"
21 : #include "gdal_frmts.h"
22 : #include "gdal_pam.h"
23 : #include "ogr_spatialref.h"
24 : #include "gmlcoverage.h"
25 :
26 : #include <algorithm>
27 :
28 : #include "wcsdataset.h"
29 : #include "wcsutils.h"
30 :
31 : using namespace WCSUtils;
32 :
33 : /************************************************************************/
34 : /* CoverageSubtype() */
35 : /* */
36 : /************************************************************************/
37 :
38 21 : static std::string CoverageSubtype(CPLXMLNode *coverage)
39 : {
40 : std::string subtype =
41 21 : CPLGetXMLValue(coverage, "ServiceParameters.CoverageSubtype", "");
42 21 : size_t pos = subtype.find("Coverage");
43 21 : if (pos != std::string::npos)
44 : {
45 21 : subtype.erase(pos);
46 : }
47 21 : return subtype;
48 : }
49 :
50 : /************************************************************************/
51 : /* GetGridNode() */
52 : /* */
53 : /************************************************************************/
54 :
55 21 : static CPLXMLNode *GetGridNode(CPLXMLNode *coverage, const std::string &subtype)
56 : {
57 21 : CPLXMLNode *grid = nullptr;
58 : // Construct the name of the node that we look under domainSet.
59 : // For now we can handle RectifiedGrid and ReferenceableGridByVectors.
60 : // Note that if this is called at GetCoverage stage, the grid should not be
61 : // NULL.
62 21 : std::string path = "domainSet";
63 21 : if (subtype == "RectifiedGrid")
64 : {
65 18 : grid = CPLGetXMLNode(coverage, (path + "." + subtype).c_str());
66 : }
67 3 : else if (subtype == "ReferenceableGrid")
68 : {
69 3 : grid = CPLGetXMLNode(coverage,
70 6 : (path + "." + subtype + "ByVectors").c_str());
71 : }
72 21 : if (!grid)
73 : {
74 0 : CPLError(CE_Failure, CPLE_AppDefined,
75 : "Can't handle coverages of type '%s'.", subtype.c_str());
76 : }
77 42 : return grid;
78 : }
79 :
80 : /************************************************************************/
81 : /* ParseParameters() */
82 : /* */
83 : /************************************************************************/
84 :
85 42 : static void ParseParameters(CPLXMLNode *service,
86 : std::vector<std::string> &dimensions,
87 : std::string &range,
88 : std::vector<std::vector<std::string>> &others)
89 : {
90 : std::vector<std::string> parameters =
91 84 : Split(CPLGetXMLValue(service, "Parameters", ""), "&");
92 84 : for (unsigned int i = 0; i < parameters.size(); ++i)
93 : {
94 42 : std::vector<std::string> kv = Split(parameters[i].c_str(), "=");
95 42 : if (kv.size() < 2)
96 : {
97 0 : continue;
98 : }
99 42 : kv[0] = CPLString(kv[0]).toupper();
100 42 : if (kv[0] == "RANGESUBSET")
101 : {
102 0 : range = kv[1];
103 : }
104 42 : else if (kv[0] == "SUBSET")
105 : {
106 0 : dimensions = Split(kv[1].c_str(), ";");
107 : }
108 : else
109 : {
110 126 : others.push_back(std::vector<std::string>{kv[0], kv[1]});
111 : }
112 : }
113 : // fallback to service values, if any
114 42 : if (range == "")
115 : {
116 42 : range = CPLGetXMLValue(service, "RangeSubset", "");
117 : }
118 42 : if (dimensions.size() == 0)
119 : {
120 42 : dimensions = Split(CPLGetXMLValue(service, "Subset", ""), ";");
121 : }
122 42 : }
123 :
124 : /************************************************************************/
125 : /* GetNativeExtent() */
126 : /* */
127 : /************************************************************************/
128 :
129 21 : std::vector<double> WCSDataset201::GetNativeExtent(int nXOff, int nYOff,
130 : int nXSize, int nYSize,
131 : CPL_UNUSED int nBufXSize,
132 : CPL_UNUSED int nBufYSize)
133 : {
134 21 : std::vector<double> extent;
135 : // WCS 2.0 extents are the outer edges of outer pixels.
136 21 : extent.push_back(m_gt[0] + (nXOff)*m_gt[1]);
137 21 : extent.push_back(m_gt[3] + (nYOff + nYSize) * m_gt[5]);
138 21 : extent.push_back(m_gt[0] + (nXOff + nXSize) * m_gt[1]);
139 21 : extent.push_back(m_gt[3] + (nYOff)*m_gt[5]);
140 21 : return extent;
141 : }
142 :
143 : /************************************************************************/
144 : /* GetCoverageRequest() */
145 : /* */
146 : /************************************************************************/
147 :
148 : std::string
149 21 : WCSDataset201::GetCoverageRequest(bool scaled, int nBufXSize, int nBufYSize,
150 : const std::vector<double> &extent,
151 : const std::string & /*osBandList*/)
152 : {
153 21 : std::string request = CPLGetXMLValue(psService, "ServiceURL", "");
154 21 : request = CPLURLAddKVP(request.c_str(), "SERVICE", "WCS");
155 21 : request += "&REQUEST=GetCoverage";
156 : request +=
157 21 : "&VERSION=" + std::string(CPLGetXMLValue(psService, "Version", ""));
158 42 : request += "&COVERAGEID=" +
159 63 : URLEncode(CPLGetXMLValue(psService, "CoverageName", ""));
160 :
161 : // note: native_crs is not really supported
162 21 : if (!native_crs)
163 : {
164 0 : std::string crs = URLEncode(CPLGetXMLValue(psService, "SRS", ""));
165 0 : request += "&OUTPUTCRS=" + crs;
166 0 : request += "&SUBSETTINGCRS=" + crs;
167 : }
168 :
169 : std::vector<std::string> domain =
170 42 : Split(CPLGetXMLValue(psService, "Domain", ""), ",");
171 21 : if (domain.size() < 2)
172 : {
173 : // eek!
174 0 : domain.push_back("E");
175 0 : domain.push_back("N");
176 : }
177 21 : const char *x = domain[0].c_str();
178 21 : const char *y = domain[1].c_str();
179 21 : if (CPLGetXMLBoolean(psService, "SubsetAxisSwap"))
180 : {
181 6 : const char *tmp = x;
182 6 : x = y;
183 6 : y = tmp;
184 : }
185 :
186 : std::vector<std::string> low =
187 42 : Split(CPLGetXMLValue(psService, "Low", ""), ",");
188 : std::vector<std::string> high =
189 42 : Split(CPLGetXMLValue(psService, "High", ""), ",");
190 42 : std::string a = CPLString().Printf("%.17g", extent[0]);
191 21 : if (low.size() > 1 && CPLAtof(low[0].c_str()) > extent[0])
192 : {
193 0 : a = low[0];
194 : }
195 42 : std::string b = CPLString().Printf("%.17g", extent[2]);
196 21 : if (high.size() > 1 && CPLAtof(high[0].c_str()) < extent[2])
197 : {
198 1 : b = high[0];
199 : }
200 : /*
201 : std::string a = CPLString().Printf(
202 : "%.17g", MAX(m_gt[0], extent[0]));
203 : std::string b = CPLString().Printf(
204 : "%.17g", MIN(m_gt[0] + nRasterXSize * m_gt[1],
205 : extent[2]));
206 : */
207 :
208 : // 09-147 KVP Protocol: subset keys must be unique
209 : // GeoServer: seems to require plain SUBSET for x and y
210 :
211 : request +=
212 21 : CPLString().Printf("&SUBSET=%s%%28%s,%s%%29", x, a.c_str(), b.c_str());
213 :
214 21 : a = CPLString().Printf("%.17g", extent[1]);
215 21 : if (low.size() > 1 && CPLAtof(low[1].c_str()) > extent[1])
216 : {
217 0 : a = low[1];
218 : }
219 21 : b = CPLString().Printf("%.17g", extent[3]);
220 21 : if (high.size() > 1 && CPLAtof(high[1].c_str()) < extent[3])
221 : {
222 0 : b = high[1];
223 : }
224 : /*
225 : a = CPLString().Printf(
226 : "%.17g", MAX(m_gt[3] + nRasterYSize * m_gt[5],
227 : extent[1])); b = CPLString().Printf(
228 : "%.17g", MIN(m_gt[3], extent[3]));
229 : */
230 :
231 : request +=
232 21 : CPLString().Printf("&SUBSET=%s%%28%s,%s%%29", y, a.c_str(), b.c_str());
233 :
234 : // Dimension and range parameters:
235 42 : std::vector<std::string> dimensions;
236 42 : std::string range;
237 42 : std::vector<std::vector<std::string>> others;
238 21 : ParseParameters(psService, dimensions, range, others);
239 :
240 : // set subsets for axis other than x/y
241 24 : for (unsigned int i = 0; i < dimensions.size(); ++i)
242 : {
243 3 : size_t pos = dimensions[i].find("(");
244 3 : std::string dim = dimensions[i].substr(0, pos);
245 3 : if (IndexOf(dim, domain) != -1)
246 : {
247 0 : continue;
248 : }
249 : std::vector<std::string> params =
250 6 : Split(FromParenthesis(dimensions[i]).c_str(), ",");
251 : request +=
252 3 : "&SUBSET" + CPLString().Printf("%i", i) + "=" + dim + "%28"; // (
253 6 : for (unsigned int j = 0; j < params.size(); ++j)
254 : {
255 : // todo: %22 (") should be used only for non-numbers
256 3 : request += "%22" + params[j] + "%22";
257 : }
258 3 : request += "%29"; // )
259 : }
260 :
261 21 : if (scaled)
262 : {
263 14 : CPLString tmp;
264 : // scaling is expressed in grid axes
265 7 : if (CPLGetXMLBoolean(psService, "UseScaleFactor"))
266 : {
267 1 : double fx = fabs((extent[2] - extent[0]) / m_gt[1] /
268 1 : ((double)nBufXSize + 0.5));
269 1 : double fy = fabs((extent[3] - extent[1]) / m_gt[5] /
270 1 : ((double)nBufYSize + 0.5));
271 1 : tmp.Printf("&SCALEFACTOR=%.15g", MIN(fx, fy));
272 : }
273 : else
274 : {
275 : std::vector<std::string> grid_axes =
276 12 : Split(CPLGetXMLValue(psService, "GridAxes", ""), ",");
277 6 : if (grid_axes.size() < 2)
278 : {
279 : // eek!
280 0 : grid_axes.push_back("E");
281 0 : grid_axes.push_back("N");
282 : }
283 : tmp.Printf("&SCALESIZE=%s%%28%i%%29,%s%%28%i%%29",
284 12 : grid_axes[0].c_str(), nBufXSize, grid_axes[1].c_str(),
285 12 : nBufYSize);
286 : }
287 7 : request += tmp;
288 : }
289 :
290 21 : if (range != "" && range != "*")
291 : {
292 0 : request += "&RANGESUBSET=" + range;
293 : }
294 :
295 : // other parameters may come from
296 : // 1) URL (others)
297 : // 2) Service file
298 21 : const char *keys[] = {WCS_URL_PARAMETERS};
299 231 : for (unsigned int i = 0; i < CPL_ARRAYSIZE(keys); i++)
300 : {
301 420 : std::string value;
302 210 : int ix = IndexOf(CPLString(keys[i]).toupper(), others);
303 210 : if (ix >= 0)
304 : {
305 0 : value = others[ix][1];
306 : }
307 : else
308 : {
309 210 : value = CPLGetXMLValue(psService, keys[i], "");
310 : }
311 210 : if (value != "")
312 : {
313 21 : request = CPLURLAddKVP(request.c_str(), keys[i], value.c_str());
314 : }
315 : }
316 : // add extra parameters
317 42 : std::string extra = CPLGetXMLValue(psService, "Parameters", "");
318 21 : if (extra != "")
319 : {
320 42 : std::vector<std::string> pairs = Split(extra.c_str(), "&");
321 42 : for (unsigned int i = 0; i < pairs.size(); ++i)
322 : {
323 21 : std::vector<std::string> pair = Split(pairs[i].c_str(), "=");
324 : request =
325 21 : CPLURLAddKVP(request.c_str(), pair[0].c_str(), pair[1].c_str());
326 : }
327 : }
328 : std::vector<std::string> pairs =
329 42 : Split(CPLGetXMLValue(psService, "GetCoverageExtra", ""), "&");
330 42 : for (unsigned int i = 0; i < pairs.size(); ++i)
331 : {
332 42 : std::vector<std::string> pair = Split(pairs[i].c_str(), "=");
333 21 : if (pair.size() > 1)
334 : {
335 : request =
336 21 : CPLURLAddKVP(request.c_str(), pair[0].c_str(), pair[1].c_str());
337 : }
338 : }
339 :
340 21 : CPLDebug("WCS", "Requesting %s", request.c_str());
341 42 : return request;
342 : }
343 :
344 : /************************************************************************/
345 : /* DescribeCoverageRequest() */
346 : /* */
347 : /************************************************************************/
348 :
349 7 : std::string WCSDataset201::DescribeCoverageRequest()
350 : {
351 7 : std::string request = CPLGetXMLValue(psService, "ServiceURL", "");
352 7 : request = CPLURLAddKVP(request.c_str(), "SERVICE", "WCS");
353 7 : request = CPLURLAddKVP(request.c_str(), "REQUEST", "DescribeCoverage");
354 14 : request = CPLURLAddKVP(request.c_str(), "VERSION",
355 14 : CPLGetXMLValue(psService, "Version", "2.0.1"));
356 14 : request = CPLURLAddKVP(request.c_str(), "COVERAGEID",
357 14 : CPLGetXMLValue(psService, "CoverageName", ""));
358 14 : std::string extra = CPLGetXMLValue(psService, "Parameters", "");
359 7 : if (extra != "")
360 : {
361 14 : std::vector<std::string> pairs = Split(extra.c_str(), "&");
362 14 : for (unsigned int i = 0; i < pairs.size(); ++i)
363 : {
364 7 : std::vector<std::string> pair = Split(pairs[i].c_str(), "=");
365 : request =
366 7 : CPLURLAddKVP(request.c_str(), pair[0].c_str(), pair[1].c_str());
367 : }
368 : }
369 7 : extra = CPLGetXMLValue(psService, "DescribeCoverageExtra", "");
370 7 : if (extra != "")
371 : {
372 0 : std::vector<std::string> pairs = Split(extra.c_str(), "&");
373 0 : for (unsigned int i = 0; i < pairs.size(); ++i)
374 : {
375 0 : std::vector<std::string> pair = Split(pairs[i].c_str(), "=");
376 : request =
377 0 : CPLURLAddKVP(request.c_str(), pair[0].c_str(), pair[1].c_str());
378 : }
379 : }
380 7 : CPLDebug("WCS", "Requesting %s", request.c_str());
381 14 : return request;
382 : }
383 :
384 : /************************************************************************/
385 : /* GridOffsets() */
386 : /* */
387 : /************************************************************************/
388 :
389 21 : bool WCSDataset201::GridOffsets(CPLXMLNode *grid, const std::string &subtype,
390 : bool swap_grid_axis,
391 : std::vector<double> &origin,
392 : std::vector<std::vector<double>> &offset,
393 : std::vector<std::string> axes, char ***metadata)
394 : {
395 : // todo: use domain_index
396 :
397 : // origin position, center of cell
398 21 : CPLXMLNode *point = CPLGetXMLNode(grid, "origin.Point.pos");
399 42 : origin = Flist(
400 63 : Split(CPLGetXMLValue(point, nullptr, ""), " ", axis_order_swap), 0, 2);
401 :
402 : // offsets = coefficients of affine transformation from cell coords to
403 : // CRS coords, (1,2) and (4,5)
404 :
405 21 : if (subtype == "RectifiedGrid")
406 : {
407 :
408 : // for rectified grid the geo transform is from origin and offsetVectors
409 18 : int i = 0;
410 126 : for (CPLXMLNode *node = grid->psChild; node != nullptr;
411 108 : node = node->psNext)
412 : {
413 126 : if (node->eType != CXT_Element ||
414 90 : !EQUAL(node->pszValue, "offsetVector"))
415 : {
416 90 : continue;
417 : }
418 36 : offset.push_back(Flist(
419 72 : Split(CPLGetXMLValue(node, nullptr, ""), " ", axis_order_swap),
420 : 0, 2));
421 36 : if (i == 1)
422 : {
423 18 : break;
424 : }
425 18 : i++;
426 : }
427 18 : if (offset.size() < 2)
428 : {
429 : // error or not?
430 0 : offset.push_back(std::vector<double>{1, 0}); // x
431 0 : offset.push_back(std::vector<double>{0, 1}); // y
432 : }
433 : // if axis_order_swap
434 : // the offset order should be swapped
435 : // Rasdaman does it
436 : // MapServer and GeoServer not
437 18 : if (swap_grid_axis)
438 : {
439 3 : std::swap(offset[0], offset[1]);
440 : }
441 : }
442 : else
443 : { // if (coverage_type == "ReferenceableGrid"(ByVector)) {
444 :
445 : // for vector referenceable grid the geo transform is from
446 : // offsetVector, coefficients, gridAxesSpanned, sequenceRule
447 : // in generalGridAxis.GeneralGridAxis
448 33 : for (CPLXMLNode *node = grid->psChild; node != nullptr;
449 30 : node = node->psNext)
450 : {
451 30 : CPLXMLNode *axis = CPLGetXMLNode(node, "GeneralGridAxis");
452 30 : if (!axis)
453 : {
454 21 : continue;
455 : }
456 9 : std::string spanned = CPLGetXMLValue(axis, "gridAxesSpanned", "");
457 9 : int index = IndexOf(spanned, axes);
458 9 : if (index == -1)
459 : {
460 0 : CPLError(CE_Failure, CPLE_AppDefined,
461 : "This is not a rectilinear grid(?).");
462 0 : return false;
463 : }
464 9 : std::string coeffs = CPLGetXMLValue(axis, "coefficients", "");
465 9 : if (coeffs != "")
466 : {
467 6 : *metadata = CSLSetNameValue(
468 : *metadata,
469 6 : CPLString().Printf("DIMENSION_%i_COEFFS", index).c_str(),
470 : coeffs.c_str());
471 : }
472 : std::string order =
473 9 : CPLGetXMLValue(axis, "sequenceRule.axisOrder", "");
474 9 : std::string rule = CPLGetXMLValue(axis, "sequenceRule", "");
475 9 : if (!(order == "+1" && rule == "Linear"))
476 : {
477 0 : CPLError(CE_Failure, CPLE_AppDefined,
478 : "Grids with sequence rule '%s' and axis order '%s' "
479 : "are not supported.",
480 : rule.c_str(), order.c_str());
481 0 : return false;
482 : }
483 9 : CPLXMLNode *offset_node = CPLGetXMLNode(axis, "offsetVector");
484 9 : if (offset_node)
485 : {
486 9 : offset.push_back(
487 18 : Flist(Split(CPLGetXMLValue(offset_node, nullptr, ""), " ",
488 9 : axis_order_swap),
489 : 0, 2));
490 : }
491 : else
492 : {
493 0 : CPLError(CE_Failure, CPLE_AppDefined,
494 : "Missing offset vector in grid axis.");
495 0 : return false;
496 : }
497 : }
498 : // todo: make sure offset order is the same as the axes order but see
499 : // above
500 : }
501 21 : if (origin.size() < 2 || offset.size() < 2)
502 : {
503 0 : CPLError(CE_Failure, CPLE_AppDefined,
504 : "Could not parse origin or offset vectors from grid.");
505 0 : return false;
506 : }
507 21 : return true;
508 : }
509 :
510 : /************************************************************************/
511 : /* GetSubdataset() */
512 : /* */
513 : /************************************************************************/
514 :
515 0 : std::string WCSDataset201::GetSubdataset(const std::string &coverage)
516 : {
517 0 : char **metadata = GDALPamDataset::GetMetadata("SUBDATASETS");
518 0 : std::string subdataset;
519 0 : if (metadata != nullptr)
520 : {
521 0 : for (int i = 0; metadata[i] != nullptr; ++i)
522 : {
523 : char *key;
524 0 : std::string url = CPLParseNameValue(metadata[i], &key);
525 0 : if (key != nullptr && strstr(key, "SUBDATASET_") &&
526 0 : strstr(key, "_NAME"))
527 : {
528 0 : if (coverage == CPLURLGetValue(url.c_str(), "coverageId"))
529 : {
530 0 : subdataset = key;
531 0 : subdataset.erase(subdataset.find("_NAME"), 5);
532 0 : CPLFree(key);
533 0 : break;
534 : }
535 : }
536 0 : CPLFree(key);
537 : }
538 : }
539 0 : return subdataset;
540 : }
541 :
542 : /************************************************************************/
543 : /* SetFormat() */
544 : /* */
545 : /************************************************************************/
546 :
547 21 : bool WCSDataset201::SetFormat(CPLXMLNode *coverage)
548 : {
549 : // set the Format value in service,
550 : // unless it is set by the user
551 42 : std::string format = CPLGetXMLValue(psService, "Format", "");
552 :
553 : // todo: check the value against list of supported formats?
554 21 : if (format != "")
555 : {
556 14 : return true;
557 : }
558 :
559 : /* We will prefer anything that sounds like TIFF, otherwise */
560 : /* falling back to the first supported format. Should we */
561 : /* consider preferring the nativeFormat if available? */
562 :
563 7 : char **metadata = GDALPamDataset::GetMetadata(nullptr);
564 : const char *value =
565 7 : CSLFetchNameValue(metadata, "WCS_GLOBAL#formatSupported");
566 7 : if (value == nullptr)
567 : {
568 0 : format = CPLGetXMLValue(coverage, "ServiceParameters.nativeFormat", "");
569 : }
570 : else
571 : {
572 14 : std::vector<std::string> format_list = Split(value, ",");
573 32 : for (unsigned j = 0; j < format_list.size(); ++j)
574 : {
575 32 : if (CPLString(format_list[j]).ifind("tiff") != std::string::npos)
576 : {
577 7 : format = format_list[j];
578 7 : break;
579 : }
580 : }
581 7 : if (format == "" && format_list.size() > 0)
582 : {
583 0 : format = format_list[0];
584 : }
585 : }
586 7 : if (format != "")
587 : {
588 7 : CPLSetXMLValue(psService, "Format", format.c_str());
589 7 : bServiceDirty = true;
590 7 : return true;
591 : }
592 : else
593 : {
594 0 : return false;
595 : }
596 : }
597 :
598 : /************************************************************************/
599 : /* ParseGridFunction() */
600 : /* */
601 : /************************************************************************/
602 :
603 21 : bool WCSDataset201::ParseGridFunction(CPLXMLNode *coverage,
604 : std::vector<int> &axisOrder)
605 : {
606 : CPLXMLNode *function =
607 21 : CPLGetXMLNode(coverage, "coverageFunction.GridFunction");
608 21 : if (function)
609 : {
610 15 : std::string path = "sequenceRule";
611 15 : std::string sequenceRule = CPLGetXMLValue(function, path.c_str(), "");
612 15 : path += ".axisOrder";
613 : axisOrder =
614 15 : Ilist(Split(CPLGetXMLValue(function, path.c_str(), ""), " "));
615 : // for now require simple
616 15 : if (sequenceRule != "Linear")
617 : {
618 0 : CPLError(CE_Failure, CPLE_AppDefined,
619 : "Can't handle '%s' coverages.", sequenceRule.c_str());
620 0 : return false;
621 : }
622 : }
623 21 : return true;
624 : }
625 :
626 : /************************************************************************/
627 : /* ParseRange() */
628 : /* */
629 : /************************************************************************/
630 :
631 21 : int WCSDataset201::ParseRange(CPLXMLNode *coverage,
632 : const std::string &range_subset, char ***metadata)
633 : {
634 21 : int fields = 0;
635 : // Default is to include all (types permitting?)
636 : // Can also be controlled with Range parameter
637 :
638 : // The contents of a rangeType is a swe:DataRecord
639 21 : const char *path = "rangeType.DataRecord";
640 21 : CPLXMLNode *record = CPLGetXMLNode(coverage, path);
641 21 : if (!record)
642 : {
643 0 : CPLError(CE_Failure, CPLE_AppDefined,
644 : "Attributes are not defined in a DataRecord, giving up.");
645 0 : return 0;
646 : }
647 :
648 : // mapserver does not like field names, it wants indexes
649 : // so we should be able to give those
650 :
651 : // if Range is set remove those not in it
652 42 : std::vector<std::string> range = Split(range_subset.c_str(), ",");
653 : // todo: add check for range subsetting profile existence in server metadata
654 : // here
655 21 : unsigned int range_index = 0; // index for reading from range
656 21 : bool in_band_range = false;
657 :
658 21 : unsigned int field_index = 1;
659 21 : std::vector<std::string> nodata_array;
660 :
661 63 : for (CPLXMLNode *field = record->psChild; field != nullptr;
662 42 : field = field->psNext)
663 : {
664 42 : if (field->eType != CXT_Element || !EQUAL(field->pszValue, "field"))
665 : {
666 0 : continue;
667 : }
668 42 : std::string fname = CPLGetXMLValue(field, "name", "");
669 42 : bool include = true;
670 :
671 42 : if (range.size() > 0)
672 : {
673 0 : include = false;
674 0 : if (range_index < range.size())
675 : {
676 0 : std::string current_range = range[range_index];
677 0 : std::string fname_test;
678 :
679 0 : if (atoi(current_range.c_str()) != 0)
680 : {
681 0 : fname_test = CPLString().Printf("%i", field_index);
682 : }
683 : else
684 : {
685 0 : fname_test = fname;
686 : }
687 :
688 0 : if (current_range == "*")
689 : {
690 0 : include = true;
691 : }
692 0 : else if (current_range == fname_test)
693 : {
694 0 : include = true;
695 0 : range_index += 1;
696 : }
697 0 : else if (current_range.find(fname_test + ":") !=
698 : std::string::npos)
699 : {
700 0 : include = true;
701 0 : in_band_range = true;
702 : }
703 0 : else if (current_range.find(":" + fname_test) !=
704 : std::string::npos)
705 : {
706 0 : include = true;
707 0 : in_band_range = false;
708 0 : range_index += 1;
709 : }
710 0 : else if (in_band_range)
711 : {
712 0 : include = true;
713 : }
714 : }
715 : }
716 :
717 42 : if (include)
718 : {
719 : const std::string key =
720 84 : CPLString().Printf("FIELD_%i_", field_index);
721 42 : *metadata = CSLSetNameValue(*metadata, (key + "NAME").c_str(),
722 : fname.c_str());
723 :
724 : std::string nodata =
725 84 : CPLGetXMLValue(field, "Quantity.nilValues.NilValue", "");
726 42 : if (nodata != "")
727 : {
728 0 : *metadata = CSLSetNameValue(*metadata, (key + "NODATA").c_str(),
729 : nodata.c_str());
730 : }
731 :
732 : std::string descr =
733 84 : CPLGetXMLValue(field, "Quantity.description", "");
734 42 : if (descr != "")
735 : {
736 18 : *metadata = CSLSetNameValue(*metadata, (key + "DESCR").c_str(),
737 : descr.c_str());
738 : }
739 :
740 42 : path = "Quantity.constraint.AllowedValues.interval";
741 42 : std::string interval = CPLGetXMLValue(field, path, "");
742 42 : if (interval != "")
743 : {
744 54 : *metadata = CSLSetNameValue(
745 54 : *metadata, (key + "INTERVAL").c_str(), interval.c_str());
746 : }
747 :
748 42 : nodata_array.push_back(std::move(nodata));
749 42 : fields += 1;
750 : }
751 :
752 42 : field_index += 1;
753 : }
754 :
755 21 : if (fields == 0)
756 : {
757 0 : CPLError(CE_Failure, CPLE_AppDefined,
758 : "No data fields found (bad Range?).");
759 : }
760 : else
761 : {
762 : // todo: default to the first one?
763 21 : bServiceDirty = CPLUpdateXML(psService, "NoDataValue",
764 59 : Join(nodata_array, ",").c_str()) ||
765 17 : bServiceDirty;
766 : }
767 :
768 21 : return fields;
769 : }
770 :
771 : /************************************************************************/
772 : /* ExtractGridInfo() */
773 : /* */
774 : /* Collect info about grid from describe coverage for WCS 2.0. */
775 : /* */
776 : /************************************************************************/
777 :
778 21 : bool WCSDataset201::ExtractGridInfo()
779 : {
780 : // this is for checking what's in service and for filling in empty slots in
781 : // it if the service file can be considered ready for use, this could be
782 : // skipped
783 :
784 21 : CPLXMLNode *coverage = CPLGetXMLNode(psService, "CoverageDescription");
785 :
786 21 : if (coverage == nullptr)
787 : {
788 0 : CPLError(CE_Failure, CPLE_AppDefined,
789 : "CoverageDescription missing from service.");
790 0 : return false;
791 : }
792 :
793 42 : std::string subtype = CoverageSubtype(coverage);
794 :
795 : // get CRS from boundedBy.Envelope and set the native flag to true
796 : // below we may set the CRS again but that won't be native (however, non
797 : // native CRS is not yet supported) also axis order swap is set
798 42 : std::string path = "boundedBy.Envelope";
799 21 : CPLXMLNode *envelope = CPLGetXMLNode(coverage, path.c_str());
800 21 : if (envelope == nullptr)
801 : {
802 0 : path = "boundedBy.EnvelopeWithTimePeriod";
803 0 : envelope = CPLGetXMLNode(coverage, path.c_str());
804 0 : if (envelope == nullptr)
805 : {
806 0 : CPLError(CE_Failure, CPLE_AppDefined, "Missing boundedBy.Envelope");
807 0 : return false;
808 : }
809 : }
810 42 : std::vector<std::string> bbox = ParseBoundingBox(envelope);
811 21 : if (!SetCRS(ParseCRS(envelope), true) || bbox.size() < 2)
812 : {
813 0 : return false;
814 : }
815 :
816 : // has the user set the domain?
817 : std::vector<std::string> domain =
818 42 : Split(CPLGetXMLValue(psService, "Domain", ""), ",");
819 :
820 : // names of axes
821 : std::vector<std::string> axes =
822 0 : Split(CPLGetXMLValue(coverage, (path + ".axisLabels").c_str(), ""), " ",
823 42 : axis_order_swap);
824 : std::vector<std::string> uoms =
825 0 : Split(CPLGetXMLValue(coverage, (path + ".uomLabels").c_str(), ""), " ",
826 42 : axis_order_swap);
827 :
828 21 : if (axes.size() < 2)
829 : {
830 0 : CPLError(CE_Failure, CPLE_AppDefined,
831 : "The coverage has less than 2 dimensions or no axisLabels.");
832 0 : return false;
833 : }
834 :
835 42 : std::vector<int> domain_indexes = IndexOf(domain, axes);
836 21 : if (Contains(domain_indexes, -1))
837 : {
838 0 : CPLError(CE_Failure, CPLE_AppDefined,
839 : "Axis in given domain does not exist in coverage.");
840 0 : return false;
841 : }
842 21 : if (domain_indexes.size() == 0)
843 : { // default is the first two
844 7 : domain_indexes.push_back(0);
845 7 : domain_indexes.push_back(1);
846 : }
847 21 : if (domain.size() == 0)
848 : {
849 7 : domain.push_back(axes[0]);
850 7 : domain.push_back(axes[1]);
851 7 : CPLSetXMLValue(psService, "Domain", Join(domain, ",").c_str());
852 7 : bServiceDirty = true;
853 : }
854 :
855 : // GridFunction (is optional)
856 : // We support only linear grid functions.
857 : // axisOrder determines how data is arranged in the grid <order><axis
858 : // number> specifically: +2 +1 => swap grid envelope and the order of the
859 : // offsets
860 42 : std::vector<int> axisOrder;
861 21 : if (!ParseGridFunction(coverage, axisOrder))
862 : {
863 0 : return false;
864 : }
865 :
866 21 : const char *md_domain = "";
867 63 : char **metadata = CSLDuplicate(
868 21 : GetMetadata(md_domain)); // coverage metadata to be added/updated
869 :
870 21 : metadata = CSLSetNameValue(metadata, "DOMAIN", Join(domain, ",").c_str());
871 :
872 : // add coverage metadata: GeoServer TimeDomain
873 :
874 : CPLXMLNode *timedomain =
875 21 : CPLGetXMLNode(coverage, "metadata.Extension.TimeDomain");
876 21 : if (timedomain)
877 : {
878 0 : std::vector<std::string> timePositions;
879 : // "//timePosition"
880 0 : for (CPLXMLNode *node = timedomain->psChild; node != nullptr;
881 0 : node = node->psNext)
882 : {
883 0 : if (node->eType != CXT_Element ||
884 0 : strcmp(node->pszValue, "TimeInstant") != 0)
885 : {
886 0 : continue;
887 : }
888 0 : for (CPLXMLNode *node2 = node->psChild; node2 != nullptr;
889 0 : node2 = node2->psNext)
890 : {
891 0 : if (node2->eType != CXT_Element ||
892 0 : strcmp(node2->pszValue, "timePosition") != 0)
893 : {
894 0 : continue;
895 : }
896 0 : timePositions.push_back(CPLGetXMLValue(node2, "", ""));
897 : }
898 : }
899 0 : metadata = CSLSetNameValue(metadata, "TimeDomain",
900 0 : Join(timePositions, ",").c_str());
901 : }
902 :
903 : // dimension metadata
904 :
905 : std::vector<std::string> slow =
906 42 : Split(bbox[0].c_str(), " ", axis_order_swap);
907 : std::vector<std::string> shigh =
908 42 : Split(bbox[1].c_str(), " ", axis_order_swap);
909 35 : bServiceDirty = CPLUpdateXML(psService, "Low", Join(slow, ",").c_str()) ||
910 14 : bServiceDirty;
911 35 : bServiceDirty = CPLUpdateXML(psService, "High", Join(shigh, ",").c_str()) ||
912 14 : bServiceDirty;
913 21 : if (slow.size() < 2 || shigh.size() < 2)
914 : {
915 0 : CPLError(CE_Failure, CPLE_AppDefined,
916 : "The coverage has less than 2 dimensions.");
917 0 : CSLDestroy(metadata);
918 0 : return false;
919 : }
920 : // todo: if our x,y domain is not the first two? use domain_indexes?
921 42 : std::vector<double> low = Flist(slow, 0, 2);
922 42 : std::vector<double> high = Flist(shigh, 0, 2);
923 42 : std::vector<double> env;
924 21 : env.insert(env.end(), low.begin(), low.begin() + 2);
925 21 : env.insert(env.end(), high.begin(), high.begin() + 2);
926 :
927 66 : for (unsigned int i = 0; i < axes.size(); ++i)
928 : {
929 90 : const std::string key = CPLString().Printf("DIMENSION_%i_", i);
930 45 : metadata =
931 45 : CSLSetNameValue(metadata, (key + "AXIS").c_str(), axes[i].c_str());
932 45 : if (i < uoms.size())
933 : {
934 39 : metadata = CSLSetNameValue(metadata, (key + "UOM").c_str(),
935 39 : uoms[i].c_str());
936 : }
937 45 : if (i < 2)
938 : {
939 42 : metadata = CSLSetNameValue(
940 84 : metadata, (key + "INTERVAL").c_str(),
941 84 : CPLString().Printf("%.15g,%.15g", low[i], high[i]));
942 : }
943 3 : else if (i < slow.size() && i < shigh.size())
944 : {
945 3 : metadata = CSLSetNameValue(
946 6 : metadata, (key + "INTERVAL").c_str(),
947 6 : CPLString().Printf("%s,%s", slow[i].c_str(), shigh[i].c_str()));
948 : }
949 0 : else if (i < bbox.size())
950 : {
951 0 : metadata = CSLSetNameValue(metadata, (key + "INTERVAL").c_str(),
952 0 : bbox[i].c_str());
953 : }
954 : }
955 :
956 : // domainSet
957 : // requirement 23: the srsName here _shall_ be the same as in boundedBy
958 : // => we ignore it
959 : // the CRS of this dataset is from boundedBy (unless it is overridden)
960 : // this is the size of this dataset
961 : // this gives the geotransform of this dataset (unless there is CRS
962 : // override)
963 :
964 21 : CPLXMLNode *grid = GetGridNode(coverage, subtype);
965 21 : if (!grid)
966 : {
967 0 : CSLDestroy(metadata);
968 0 : return false;
969 : }
970 :
971 : //
972 21 : bool swap_grid_axis = false;
973 30 : if (axisOrder.size() >= 2 && axisOrder[domain_indexes[0]] == 2 &&
974 9 : axisOrder[domain_indexes[1]] == 1)
975 : {
976 9 : swap_grid_axis = !CPLGetXMLBoolean(psService, "NoGridAxisSwap");
977 : }
978 21 : path = "limits.GridEnvelope";
979 : std::vector<std::vector<int>> size =
980 42 : ParseGridEnvelope(CPLGetXMLNode(grid, path.c_str()), swap_grid_axis);
981 42 : std::vector<int> grid_size;
982 21 : if (size.size() < 2)
983 : {
984 0 : CPLError(CE_Failure, CPLE_AppDefined, "Can't parse the grid envelope.");
985 0 : CSLDestroy(metadata);
986 0 : return false;
987 : }
988 :
989 21 : grid_size.push_back(size[1][domain_indexes[0]] -
990 21 : size[0][domain_indexes[0]] + 1);
991 21 : grid_size.push_back(size[1][domain_indexes[1]] -
992 21 : size[0][domain_indexes[1]] + 1);
993 :
994 21 : path = "axisLabels";
995 : bool swap_grid_axis_labels =
996 21 : swap_grid_axis || CPLGetXMLBoolean(psService, "GridAxisLabelSwap");
997 : std::vector<std::string> grid_axes = Split(
998 42 : CPLGetXMLValue(grid, path.c_str(), ""), " ", swap_grid_axis_labels);
999 : // autocorrect MapServer thing
1000 24 : if (grid_axes.size() >= 2 && grid_axes[0] == "lat" &&
1001 3 : grid_axes[1] == "long")
1002 : {
1003 3 : grid_axes[0] = "long";
1004 3 : grid_axes[1] = "lat";
1005 : }
1006 21 : bServiceDirty =
1007 35 : CPLUpdateXML(psService, "GridAxes", Join(grid_axes, ",").c_str()) ||
1008 14 : bServiceDirty;
1009 :
1010 42 : std::vector<double> origin;
1011 42 : std::vector<std::vector<double>> offsets;
1012 21 : if (!GridOffsets(grid, subtype, swap_grid_axis, origin, offsets, axes,
1013 : &metadata))
1014 : {
1015 0 : CSLDestroy(metadata);
1016 0 : return false;
1017 : }
1018 :
1019 21 : SetGeometry(grid_size, origin, offsets);
1020 :
1021 : // subsetting and dimension to bands
1022 42 : std::vector<std::string> dimensions;
1023 42 : std::string range;
1024 42 : std::vector<std::vector<std::string>> others;
1025 21 : ParseParameters(psService, dimensions, range, others);
1026 :
1027 : // it is ok to have trimming or even slicing for x/y, it just affects our
1028 : // bounding box but that is a todo item todo: BoundGeometry(domain_trim) if
1029 : // domain_trim.size() > 0
1030 42 : std::vector<std::vector<double>> domain_trim;
1031 :
1032 : // are all dimensions that are not x/y domain sliced?
1033 : // if not, bands can't be defined, see below
1034 21 : bool dimensions_are_ok = true;
1035 66 : for (unsigned int i = 0; i < axes.size(); ++i)
1036 : {
1037 45 : std::vector<std::string> params;
1038 51 : for (unsigned int j = 0; j < dimensions.size(); ++j)
1039 : {
1040 9 : if (dimensions[j].find(axes[i] + "(") != std::string::npos)
1041 : {
1042 3 : params = Split(FromParenthesis(dimensions[j]).c_str(), ",");
1043 3 : break;
1044 : }
1045 : }
1046 45 : int domain_index = IndexOf(axes[i], domain);
1047 45 : if (domain_index != -1)
1048 : {
1049 42 : domain_trim.push_back(Flist(params, 0, 2));
1050 42 : continue;
1051 : }
1052 : // size == 1 => sliced
1053 3 : if (params.size() != 1)
1054 : {
1055 0 : dimensions_are_ok = false;
1056 : }
1057 : }
1058 : // todo: add metadata: note: no bands, you need to subset to get data
1059 :
1060 : // check for CRS override
1061 42 : std::string crs = CPLGetXMLValue(psService, "SRS", "");
1062 21 : if (crs != "" && crs != osCRS)
1063 : {
1064 0 : if (!SetCRS(crs, false))
1065 : {
1066 0 : CSLDestroy(metadata);
1067 0 : return false;
1068 : }
1069 : // todo: support CRS override, it requires warping the grid to the new
1070 : // CRS SetGeometry(grid_size, origin, offsets);
1071 0 : CPLError(CE_Failure, CPLE_AppDefined,
1072 : "CRS override not yet supported.");
1073 0 : CSLDestroy(metadata);
1074 0 : return false;
1075 : }
1076 :
1077 : // todo: ElevationDomain, DimensionDomain
1078 :
1079 : // rangeType
1080 :
1081 : // get the field metadata
1082 : // get the count of fields
1083 : // if Range is set in service that may limit the fields
1084 21 : int fields = ParseRange(coverage, range, &metadata);
1085 : // if fields is 0 an error message has been emitted
1086 : // but we let this go on since the user may be experimenting
1087 : // and she wants to see the resulting metadata and not just an error message
1088 : // situation is ~the same when bands == 0 when we exit here
1089 :
1090 : // todo: do this only if metadata is dirty
1091 21 : this->SetMetadata(metadata, md_domain);
1092 21 : CSLDestroy(metadata);
1093 21 : TrySaveXML();
1094 :
1095 : // determine the band count
1096 21 : int bands = 0;
1097 21 : if (dimensions_are_ok)
1098 : {
1099 21 : bands = fields;
1100 : }
1101 21 : bServiceDirty =
1102 35 : CPLUpdateXML(psService, "BandCount", CPLString().Printf("%d", bands)) ||
1103 14 : bServiceDirty;
1104 :
1105 : // set the Format value in service, unless it is set
1106 : // by the user, either through direct edit or options
1107 21 : if (!SetFormat(coverage))
1108 : {
1109 : // all attempts to find a format have failed...
1110 0 : CPLError(CE_Failure, CPLE_AppDefined,
1111 : "All attempts to find a format have failed, giving up.");
1112 0 : return false;
1113 : }
1114 :
1115 21 : return true;
1116 : }
|
