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