bes Updated for version 3.20.13
ugrid_restrict.cc
1// -*- mode: c++; c-basic-offset:4 -*-
2
3// This file is part of libdap, A C++ implementation of the OPeNDAP Data
4// Access Protocol.
5
6// Copyright (c) 2002,2003,2011,2012 OPeNDAP, Inc.
7// Authors: Nathan Potter <ndp@opendap.org>
8// James Gallagher <jgallagher@opendap.org>
9// Scott Moe <smeest1@gmail.com>
10// Bill Howe <billhowe@cs.washington.edu>
11//
12// This library is free software; you can redistribute it and/or
13// modify it under the terms of the GNU Lesser General Public
14// License as published by the Free Software Foundation; either
15// version 2.1 of the License, or (at your option) any later version.
16//
17// This library is distributed in the hope that it will be useful,
18// but WITHOUT ANY WARRANTY; without even the implied warranty of
19// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
20// Lesser General Public License for more details.
21//
22// You should have received a copy of the GNU Lesser General Public
23// License along with this library; if not, write to the Free Software
24// Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
25//
26// You can contact OPeNDAP, Inc. at PO Box 112, Saunderstown, RI. 02874-0112.
27
28// NOTE: This file is built only when the gridfields library is linked with
29// the netcdf_handler (i.e., the handler's build is configured using the
30// --with-gridfields=... option to the 'configure' script).
31
32#include "config.h"
33
34#include <limits.h>
35
36#include <cstdlib> // used by strtod()
37#include <cerrno>
38#include <cmath>
39#include <iostream>
40#include <sstream>
41//#include <cxxabi.h>
42
43#include <curl/curl.h>
44
45#include <libdap/BaseType.h>
46#include <libdap/Int32.h>
47#include <libdap/Str.h>
48#include <libdap/Array.h>
49#include <libdap/Structure.h>
50#include <libdap/Error.h>
51#include <libdap/util.h>
52#include <libdap/escaping.h>
53
54#include "BESDebug.h"
55#include "BESError.h"
56#include "BESStopWatch.h"
57#include <libdap/util.h>
58
59#include "ugrid_utils.h"
60#include "MeshDataVariable.h"
61#include "TwoDMeshTopology.h"
62#include "NDimensionalArray.h"
63#include <gridfields/GFError.h>
64
65#include "ugrid_restrict.h"
66
67#ifdef NDEBUG
68#undef BESDEBUG
69#define BESDEBUG( x, y )
70#endif
71
72using namespace std;
73using namespace libdap;
74
75namespace ugrid {
76
84struct UgridRestrictArgs {
94 locationType dimension;
95
102 vector<libdap::Array *> rangeVars;
103
107 string filterExpression;
108};
109
115static void addRangeVar(DDS *dds, libdap::Array *rangeVar, map<string, vector<MeshDataVariable *> *> *rangeVariables)
116{
117 MeshDataVariable *mdv = new MeshDataVariable();
118 mdv->init(rangeVar);
119 string meshVarName = mdv->getMeshName();
120
121 BaseType *meshVar = dds->var(meshVarName);
122
123 if (meshVar == 0) {
124 string msg = "The range variable '" + mdv->getName() + "' references the mesh variable '" + meshVarName
125 + "' which cannot be located in this dataset.";
126 throw Error(malformed_expr, msg);
127 }
128
129 // Get the rangeVariable vector for this mesh name from the map.
130 vector<MeshDataVariable *> *requestedRangeVarsForMesh;
131 map<string, vector<MeshDataVariable *> *>::iterator mit = rangeVariables->find(meshVarName);
132 if (mit == rangeVariables->end()) {
133 // Not there? Make a new one.
134 BESDEBUG("ugrid",
135 "addRangeVar() - MeshTopology object for '" << meshVarName <<"' does not exist. Getting a new one... " << endl);
136
137 requestedRangeVarsForMesh = new vector<MeshDataVariable *>();
138 (*rangeVariables)[meshVarName] = requestedRangeVarsForMesh;
139 }
140 else {
141 // Sweet! Found it....
142 BESDEBUG("ugrid",
143 "addRangeVar() - MeshTopology object for '" << meshVarName <<"' exists. Retrieving... " << endl);
144 requestedRangeVarsForMesh = mit->second;
145 }
146
147 requestedRangeVarsForMesh->push_back(mdv);
148}
149
150
151string usage(string fnc){
152
153 string usage = fnc+"(rangeVariable:string, [rangeVariable:string, ... ] condition:string)";
154
155 return usage;
156}
157
161static UgridRestrictArgs processUgrArgs(string func_name, locationType dimension, int argc, BaseType *argv[])
162{
163
164 BESDEBUG("ugrid", "processUgrArgs() - BEGIN" << endl);
165
166 UgridRestrictArgs args;
167 args.dimension = dimension;
168 BESDEBUG("ugrid", "args.dimension: " << libdap::long_to_string(args.dimension) << endl);
169
170 args.rangeVars = vector<libdap::Array *>();
171
172 // Check number of arguments;
173 if (argc < 2)
174 throw Error(malformed_expr,
175 "Wrong number of arguments to ugrid restrict function: " + usage(func_name) + " was passed " + long_to_string(argc)
176 + " argument(s)");
177
178 BaseType * bt;
179
180#if 0
181 // ---------------------------------------------
182 // Process the first arg, which is the rank of the Restriction Clause
183 bt = argv[0];
184 if (bt->type() != dods_int32_c)
185 throw Error(malformed_expr,
186 "Wrong type for first argument, expected DAP Int32. " + ugrSyntax + " was passed a/an " + bt->type_name());
187
188 args.dimension = (locationType) dynamic_cast<Int32&>(*argv[0]).value();
189 BESDEBUG("ugrid", "args.dimension: " << libdap::long_to_string(args.dimension) << endl);
190#endif
191
192 // ---------------------------------------------
193 // Process the last argument, the relational/filter expression used to restrict the ugrid content.
194 bt = argv[argc - 1];
195 if (bt->type() != dods_str_c)
196 throw Error(malformed_expr,
197 "Wrong type for third argument, expected DAP String. " + usage(func_name) + " was passed a/an "
198 + bt->type_name());
199
200 args.filterExpression = dynamic_cast<Str&>(*bt).value();
201 BESDEBUG("ugrid", "args.filterExpression: '" << args.filterExpression << "' (AS RECEIVED)" << endl);
202
203 args.filterExpression = www2id(args.filterExpression);
204
205 BESDEBUG("ugrid", "args.filterExpression: '" << args.filterExpression << "' (URL DECODED)" << endl);
206
207 // --------------------------------------------------
208 // Process the range variables selected by the user.
209 // We know that argc>=3, because we checked so the
210 // following loop will try to find at least one rangeVar,
211 // and it won't try to process the first or last members
212 // of argv.
213 for (int i = 0; i < (argc - 1); i++) {
214 bt = argv[i];
215 if (bt->type() != dods_array_c)
216 throw Error(malformed_expr,
217 "Wrong type for second argument, expected DAP Array. " + usage(func_name) + " was passed a/an "
218 + bt->type_name());
219
220 libdap::Array *newRangeVar = dynamic_cast<libdap::Array*>(bt);
221 if (newRangeVar == 0) {
222 throw Error(malformed_expr,
223 "Wrong type for second argument. " + usage(func_name) + " was passed a/an " + bt->type_name());
224 }
225 args.rangeVars.push_back(newRangeVar);
226 }
227
228 BESDEBUG("ugrid", "processUgrArgs() - END" << endl);
229
230 return args;
231}
232
233static string arrayState(libdap::Array *dapArray, string indent)
234{
235
236 libdap::Array::Dim_iter thisDim;
237 stringstream arrayState;
238 arrayState << endl;
239 arrayState << indent << "dapArray: '" << dapArray->name() << "' ";
240 arrayState << "type: '" << dapArray->type_name() << "' ";
241
242 arrayState << "read(): " << dapArray->read() << " ";
243 arrayState << "read_p(): " << dapArray->read_p() << " ";
244 arrayState << endl;
245
246 for (thisDim = dapArray->dim_begin(); thisDim != dapArray->dim_end(); thisDim++) {
247
248 arrayState << indent << " dim: '" << dapArray->dimension_name(thisDim) << "' ";
249 arrayState << indent << " start: " << dapArray->dimension_start(thisDim) << " ";
250 arrayState << indent << " stride: " << dapArray->dimension_stride(thisDim) << " ";
251 arrayState << indent << " stop: " << dapArray->dimension_stop(thisDim) << " ";
252 arrayState << endl;
253 }
254 return arrayState.str();
255}
256
257static void copyUsingSubsetIndex(libdap::Array *sourceArray, vector<unsigned int> *subsetIndex, void *result)
258{
259 BESDEBUG("ugrid", "ugrid::copyUsingIndex() - BEGIN" << endl);
260
261 switch (sourceArray->var()->type()) {
262 case dods_byte_c:
263 sourceArray->value(subsetIndex, (dods_byte *) result);
264 break;
265 case dods_uint16_c:
266 sourceArray->value(subsetIndex, (dods_uint16 *) result);
267 break;
268 case dods_int16_c:
269 sourceArray->value(subsetIndex, (dods_int16 *) result);
270 break;
271 case dods_uint32_c:
272 sourceArray->value(subsetIndex, (dods_uint32 *) result);
273 break;
274 case dods_int32_c:
275 sourceArray->value(subsetIndex, (dods_int32 *) result);
276 break;
277
278 case dods_float32_c:
279 sourceArray->value(subsetIndex, (dods_float32 *) result);
280 break;
281 case dods_float64_c:
282 sourceArray->value(subsetIndex, (dods_float64 *) result);
283 break;
284
285 default:
286 throw InternalErr(__FILE__, __LINE__, "ugrid::hgr5::copyUsingSubsetIndex() - Unknown DAP type encountered.");
287 }
288
289 BESDEBUG("ugrid", "ugrid::copyUsingIndex() - END" << endl);
290}
291
292// This is only used for the ugrid2 BESDEBUG lines.
293static string vectorToString(vector<unsigned int> *index)
294{
295 BESDEBUG("ugrid", "indexToString() - BEGIN"<< endl);
296 BESDEBUG("ugrid", "indexToString() - index.size(): " << libdap::long_to_string(index->size()) << endl);
297 stringstream s;
298 s << "[";
299
300 for (unsigned int i = 0; i < index->size(); ++i) {
301 s << ((i > 0) ? ", " : " ");
302 s << (*index)[i];
303 }
304 s << "]";
305 BESDEBUG("ugrid", "indexToString() - END"<< endl);
306
307 return s.str();
308}
309
313static void rDAWorker(MeshDataVariable *mdv, libdap::Array::Dim_iter thisDim, vector<unsigned int> *slab_subset_index,
314 NDimensionalArray *results)
315{
316 libdap::Array *dapArray = mdv->getDapArray();
317
318 // For real data, this output is huge. jhrg 4/15/15
319 BESDEBUG("ugrid2",
320 "rDAWorker() - slab_subset_index" << vectorToString(slab_subset_index) << " size: " << slab_subset_index->size() << endl);
321
322 // The locationCoordinateDimension is the dimension of the array that is associated with the ugrid "rank" - e.g. it is the
323 // dimension that ties the variable to the 'nodes' (rank 0) or 'edges' (rank 1) or 'faces' (rank 2) of the ugrid.
324 libdap::Array::Dim_iter locationCoordinateDim = mdv->getLocationCoordinateDimension();
325
326 BESDEBUG("ugrid", "rdaWorker() - thisDim: '" << dapArray->dimension_name(thisDim) << "'" << endl);
327 BESDEBUG("ugrid",
328 "rdaWorker() - locationCoordinateDim: '" << dapArray->dimension_name(locationCoordinateDim) << "'" << endl);
329 // locationCoordinateDim is, e.g., 'nodes'. jhrg 10/25/13
330 if (thisDim != locationCoordinateDim) {
331 unsigned int start = dapArray->dimension_start(thisDim, true);
332 unsigned int stride = dapArray->dimension_stride(thisDim, true);
333 unsigned int stop = dapArray->dimension_stop(thisDim, true);
334
335 BESDEBUG("ugrid", "rdaWorker() - array state: " << arrayState(dapArray, " "));
336
337 for (unsigned int dimIndex = start; dimIndex <= stop; dimIndex += stride) {
338 dapArray->add_constraint(thisDim, dimIndex, 1, dimIndex);
339 rDAWorker(mdv, thisDim + 1, slab_subset_index, results);
340 }
341
342 // Reset the constraint for this dimension.
343 dapArray->add_constraint(thisDim, start, stride, stop);
344 }
345 else {
346 BESDEBUG("ugrid", "rdaWorker() - Found locationCoordinateDim" << endl);
347
348 if ((thisDim + 1) != dapArray->dim_end()) {
349 string msg =
350 "rDAWorker() - The location coordinate dimension is not the last dimension in the array. Hyperslab subsetting of this dimension is not supported.";
351 BESDEBUG("ugrid", msg << endl);
352 throw Error(malformed_expr, msg);
353 }
354
355 BESDEBUG("ugrid", "rdaWorker() - Arrived At Last Dimension" << endl);
356
357 dapArray->set_read_p(false);
358
359 BESDEBUG("ugrid", "rdaWorker() - dap array: " << arrayState(dapArray, " "));
360
361 vector<unsigned int> lastDimHyperSlabLocation;
362 NDimensionalArray::retrieveLastDimHyperSlabLocationFromConstrainedArrray(dapArray, &lastDimHyperSlabLocation);
363
364 BESDEBUG("ugrid",
365 "rdaWorker() - lastDimHyperSlabLocation: " << NDimensionalArray::vectorToIndices(&lastDimHyperSlabLocation) << endl);
366
367 // unused. 4/7/14 jhrg. unsigned int elementCount;
368
369 void *slab;
370 //results->getLastDimensionHyperSlab(&lastDimHyperSlabLocation, &slab, &elementCount);
371 results->getNextLastDimensionHyperSlab(&slab);
372
373 dapArray->read();
374
375 copyUsingSubsetIndex(dapArray, slab_subset_index, slab);
376 }
377}
378
386static libdap::Array *restrictRangeVariableByOneDHyperSlab(MeshDataVariable *mdv,
387 vector<unsigned int> *slab_subset_index)
388{
389
390 long restrictedSlabSize = slab_subset_index->size();
391
392 BESDEBUG("ugrid2",
393 "restrictRangeVariableByOneDHyperSlab() - slab_subset_index"<< vectorToString(slab_subset_index) << " size: " << libdap::long_to_string(restrictedSlabSize) << endl);
394
395 libdap::Array *sourceDapArray = mdv->getDapArray();
396
397 BESDEBUG("ugrid",
398 "restrictDapArrayByOneDHyperSlab() - locationCoordinateDim: '" << sourceDapArray->dimension_name(mdv->getLocationCoordinateDimension()) << "'" << endl);
399
400 // We want the manipulate the Array's Dimensions so that only a single dimensioned slab of the location coordinate dimension
401 // is read at a time. We need to cache the original constrained dimensions so that we can build the correct collection of
402 // location coordinate dimension slabs.
403
404 // Now we need to compute the shape of the final subset result array from the source range variable array and the slab subset.
405 // What's the shape of the source array with any constraint applied?
406 vector<unsigned int> resultArrayShape(sourceDapArray->dimensions(true));
407 NDimensionalArray::computeConstrainedShape(sourceDapArray, &resultArrayShape);
408
409 stringstream msg;
410 for (unsigned int i = 0; i < resultArrayShape.size(); i++) {
411 msg << "[" << resultArrayShape[i] << "]";
412 }
413 BESDEBUG("ugrid", "restrictDapArrayByOneDHyperSlab() - Constrained source array shape" << msg.str() << endl);
414 msg.str("");
415
416 // Now, we know that the result array has a last dimension slab size determined by the slab_subset_index (whichwas made by the ugrid sub-setting),
417 // so we make the result array shape reflect that
418
419 resultArrayShape[sourceDapArray->dimensions(true) - 1] = restrictedSlabSize; // jhrg 10/25/13 resultArrayShape.last() = ...
420 libdap::Type dapType = sourceDapArray->var()->type();
421
422 BESDEBUG("ugrid",
423 "restrictDapArrayByOneDHyperSlab() - UGrid restricted HyperSlab has "<< restrictedSlabSize << " elements." << endl);
424 BESDEBUG("ugrid",
425 "restrictDapArrayByOneDHyperSlab() - Array is of type '"<< libdap::type_name(dapType) << "'" << endl);
426
427 for (unsigned int i = 0; i < resultArrayShape.size(); i++) {
428 msg << "[" << resultArrayShape[i] << "]";
429 }
430 BESDEBUG("ugrid", "restrictDapArrayByOneDHyperSlab() - resultArrayShape" << msg.str() << endl);
431 msg.str(std::string());
432
433 // Now we make a new NDimensionalArray instance that we will use to hold the results.
434 NDimensionalArray *result = new NDimensionalArray(&resultArrayShape, dapType);
435
436 // And we pass that along with other stuff into the recursive rDAWorker that's going to go get all the stuff
437 rDAWorker(mdv, sourceDapArray->dim_begin(), slab_subset_index, result);
438
439 // And now that the recursion we grab have the NDimensionalArray cough up the rteuslt as a libdap::Array
440 libdap::Array *resultDapArray = result->getArray(sourceDapArray);
441
442 // Delete the NDimensionalArray
443 delete result;
444
445 // Return the result as libdap:Array
446 return resultDapArray;
447}
448
468void ugrid_restrict(string func_name, locationType location, int argc, BaseType *argv[], DDS &dds, BaseType **btpp)
469{
470 try { // This top level try block is used to catch gridfields library errors.
471
472 BESStopWatch sw;
473 if (BESDebug::IsSet(TIMING_LOG_KEY)) sw.start("ugrid::ugrid_restrict()", "[function_invocation]");
474
475 BESDEBUG("ugrid", "ugrid_restrict() - BEGIN" << endl);
476
477 string info = string("<?xml version=\"1.0\" encoding=\"UTF-8\"?>\n")
478 + "<function name=\"ugrid_restrict\" version=\"0.1\">\n" + "Server function for Unstructured grid operations.\n"
479 + "usage: " + usage(func_name) + "\n" + "</function>";
480
481 if (argc == 0) {
482 Str *response = new Str("info");
483 response->set_value(info);
484 *btpp = response;
485 return;
486 }
487
488 // Process and QC the arguments
489 UgridRestrictArgs args = processUgrArgs(func_name, location, argc, argv);
490
491 // Each range variable is associated with a "mesh" i.e. a mesh topology variable. Since there may be more than one mesh in a
492 // dataset, and the user may request more than one range variable for each mesh we need to sift through the list of requested
493 // range variables and organize them by mesh topology variable name.
494 map<string, vector<MeshDataVariable *> *> *meshToRangeVarsMap = new map<string, vector<MeshDataVariable *> *>();
495
496 // For every Range variable in the arguments list, locate it and ingest it.
497 vector<libdap::Array *>::iterator it;
498 for (it = args.rangeVars.begin(); it != args.rangeVars.end(); ++it) {
499 addRangeVar(&dds, *it, meshToRangeVarsMap);
500 }
501 BESDEBUG("ugrid", "ugrid_restrict() - The user requested "<< args.rangeVars.size() << " range data variables." << endl);
502 BESDEBUG("ugrid",
503 "ugrid_restrict() - The user's request referenced "<< meshToRangeVarsMap->size() << " mesh topology variables." << endl);
504
505 // ----------------------------------
506 // OK, so up to this point we have not read any data from the data set, but we have QC'd the inputs and verified that
507 // it looks like the request is consistent with the semantics of the dataset.
508 // Now it's time to read some data and pack it into the GridFields library...
509
510 // TODO This returns a single structure but it would make better sense to the
511 // world if it could return a vector of objects and have them appear at the
512 // top level of the DDS.
513 // FIXME fix the names of the variables in the mesh_topology attributes
514 // If the server side function can be made to return a DDS or a collection of BaseType's then the
515 // names won't change and the original mesh_topology variable and it's metadata will be valid
516 Structure *dapResult = new Structure("ugr_result_unwrap");
517
518 // Now we need to grab an top level metadata (attriubutes) and copy them into the dapResult Structure
519 // Add any global attributes to the netcdf file
520#if 1
521 AttrTable &globals = dds.get_attr_table();
522 BESDEBUG("ugrid", "ugrid_restrict() - Copying Global Attributes" << endl << globals << endl);
523
524 AttrTable *newDatasetAttrTable = new AttrTable(globals);
525 dapResult->set_attr_table(*newDatasetAttrTable);
526 BESDEBUG("ugrid", "ugrid_restrict() - Result Structure attrs: " << endl << dapResult->get_attr_table() << endl);
527
528
529#endif
530
531 // Since we only want each ugrid structure to appear in the results one time (cause otherwise we might be trying to add
532 // the same variables with the same names to the result multiple times.) we grind on this by iterating over the
533 // names of the mesh topology names.
534 map<string, vector<MeshDataVariable *> *>::iterator mit;
535 for (mit = meshToRangeVarsMap->begin(); mit != meshToRangeVarsMap->end(); ++mit) {
536
537 string meshVariableName = mit->first;
538 vector<MeshDataVariable *> *requestedRangeVarsForMesh = mit->second;
539
540 // Building the restricted TwoDMeshTopology without adding any range variables and then converting the result
541 // Grid field to Dap Objects should return all of the Ugrid structural stuff - mesh variable, node coordinate variables,
542 // face and edge coordinate variables if present.
543 BESDEBUG("ugrid",
544 "ugrid_restrict() - Adding restricted mesh_topology structure for mesh '" << meshVariableName << "' to DAP response." << endl);
545
546 TwoDMeshTopology *tdmt = new TwoDMeshTopology();
547 tdmt->init(meshVariableName, &dds);
548
549 tdmt->buildBasicGfTopology();
550 tdmt->addIndexVariable(node);
551 tdmt->addIndexVariable(face);
552 tdmt->applyRestrictOperator(args.dimension, args.filterExpression);
553
554 // 3: because there are nodes (rank = 0), edges (rank = 1), and faces (rank = 2). jhrg 10/25/13
555 vector<vector<unsigned int> *> location_subset_indices(3);
556
557 long nodeResultSize = tdmt->getResultGridSize(node);
558 BESDEBUG("ugrid", "ugrid_restrict() - there are "<< nodeResultSize << " nodes in the subset." << endl);
559 vector<unsigned int> node_subset_index(nodeResultSize);
560 if (nodeResultSize > 0) {
561 tdmt->getResultIndex(node, node_subset_index.data());
562 }
563 location_subset_indices[node] = &node_subset_index;
564
565 BESDEBUG("ugrid2", "ugrid_restrict() - node_subset_index"<< vectorToString(&node_subset_index) << endl);
566
567 long faceResultSize = tdmt->getResultGridSize(face);
568 BESDEBUG("ugrid", "ugrid_restrict() - there are "<< faceResultSize << " faces in the subset." << endl);
569 vector<unsigned int> face_subset_index(faceResultSize);
570 if (faceResultSize > 0) {
571 tdmt->getResultIndex(face, face_subset_index.data());
572 }
573 location_subset_indices[face] = &face_subset_index;
574 BESDEBUG("ugrid2", "ugrid_restrict() - face_subset_index: "<< vectorToString(&face_subset_index) << endl);
575
576 // This gets all the stuff that's attached to the grid - which at this point does not include the range variables but does include the
577 // index variable. good enough for now but need to drop the index....
578 vector<BaseType *> dapResults;
579 tdmt->convertResultGridFieldStructureToDapObjects(&dapResults);
580
581 BESDEBUG("ugrid",
582 "ugrid_restrict() - Restriction of mesh_topology '"<< tdmt->getMeshVariable()->name() << "' structure completed." << endl);
583
584 // now that we have the mesh topology variable we are going to look at each of the requested
585 // range variables (aka MeshDataVariable instances) and we're going to subset that using the
586 // gridfields library and add its subset version to the results.
587 vector<MeshDataVariable *>::iterator rvit;
588 for (rvit = requestedRangeVarsForMesh->begin(); rvit != requestedRangeVarsForMesh->end(); rvit++) {
589 MeshDataVariable *mdv = *rvit;
590
591 BESDEBUG("ugrid",
592 "ugrid_restrict() - Processing MeshDataVariable '"<< mdv->getName() << "' associated with rank/location: "<< mdv->getGridLocation() << endl);
593
594 tdmt->setLocationCoordinateDimension(mdv);
595
600 libdap::Array *restrictedRangeVarArray = restrictRangeVariableByOneDHyperSlab(mdv,
601 location_subset_indices[mdv->getGridLocation()]);
602
603 BESDEBUG("ugrid",
604 "ugrid_restrict() - Adding resulting dapArray '"<< restrictedRangeVarArray->name() << "' to dapResults." << endl);
605
606 dapResults.push_back(restrictedRangeVarArray);
607 }
608
609 delete tdmt;
610
611 BESDEBUG("ugrid", "ugrid_restrict() - Adding GF::GridField results to DAP structure " << dapResult->name() << endl);
612
613 for (vector<BaseType *>::iterator i = dapResults.begin(); i != dapResults.end(); ++i) {
614 BESDEBUG("ugrid",
615 "ugrid_restrict() - Adding variable "<< (*i)->name() << " to DAP structure " << dapResult->name() << endl);
616 dapResult->add_var_nocopy(*i);
617 }
618 }
619
620 *btpp = dapResult;
621
622 // TODO replace with auto_ptr. jhrg 4/17/15
623
624 BESDEBUG("ugrid", "ugrid_restrict() - Releasing maps and vectors..." << endl);
625 for (mit = meshToRangeVarsMap->begin(); mit != meshToRangeVarsMap->end(); ++mit) {
626 vector<MeshDataVariable *> *requestedRangeVarsForMesh = mit->second;
627 vector<MeshDataVariable *>::iterator rvit;
628 for (rvit = requestedRangeVarsForMesh->begin(); rvit != requestedRangeVarsForMesh->end(); rvit++) {
629 MeshDataVariable *mdv = *rvit;
630 delete mdv;
631 }
632 delete requestedRangeVarsForMesh;
633 }
634 delete meshToRangeVarsMap;
635
636 BESDEBUG("ugrid", "ugrid_restrict() - END" << endl);
637 }
638 catch (GFError &gfe) {
639 throw BESError(gfe.get_message(), gfe.get_error_type(), gfe.get_file(), gfe.get_line());
640 }
641
642 return;
643}
644
645
646
647
648
667void ugnr(int argc, BaseType *argv[], DDS &dds, BaseType **btpp) {
668 ugrid_restrict("ugnr",node,argc,argv,dds,btpp);
669}
670
671
672
691void uger(int argc, BaseType *argv[], DDS &dds, BaseType **btpp) {
692 ugrid_restrict("uger",edge,argc,argv,dds,btpp);
693}
694
695
714void ugfr(int argc, BaseType *argv[], DDS &dds, BaseType **btpp) {
715 ugrid_restrict("ugfr",face,argc,argv,dds,btpp);
716}
717
718
719
720
721} // namespace ugrid_restrict
static bool IsSet(const std::string &flagName)
see if the debug context flagName is set to true
Definition: BESDebug.h:168
Base exception class for the BES with basic string message.
Definition: BESError.h:59
virtual bool start(std::string name)
Definition: BESStopWatch.cc:67
static long computeConstrainedShape(libdap::Array *a, vector< unsigned int > *shape)