libdap++ Updated for version 3.8.2
|
00001 // -*- mode: c++; c-basic-offset:4 -*- 00002 00003 // This file is part of libdap, A C++ implementation of the OPeNDAP Data 00004 // Access Protocol. 00005 00006 // Copyright (c) 2002,2003 OPeNDAP, Inc. 00007 // Author: James Gallagher <jgallagher@opendap.org> 00008 // 00009 // This library is free software; you can redistribute it and/or 00010 // modify it under the terms of the GNU Lesser General Public 00011 // License as published by the Free Software Foundation; either 00012 // version 2.1 of the License, or (at your option) any later version. 00013 // 00014 // This library is distributed in the hope that it will be useful, 00015 // but WITHOUT ANY WARRANTY; without even the implied warranty of 00016 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU 00017 // Lesser General Public License for more details. 00018 // 00019 // You should have received a copy of the GNU Lesser General Public 00020 // License along with this library; if not, write to the Free Software 00021 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA 00022 // 00023 // You can contact OPeNDAP, Inc. at PO Box 112, Saunderstown, RI. 02874-0112. 00024 00025 // (c) COPYRIGHT URI/MIT 1999 00026 // Please read the full copyright statement in the file COPYRIGHT_URI. 00027 // 00028 // Authors: 00029 // jhrg,jimg James Gallagher <jgallagher@gso.uri.edu> 00030 00031 00032 // These functions are used by the CE evaluator 00033 // 00034 // 1/15/99 jhrg 00035 00036 #include "config.h" 00037 00038 static char rcsid[]not_used = 00039 { "$Id: ce_functions.cc 23551 2010-09-13 21:00:46Z jimg $" 00040 }; 00041 00042 #include <limits.h> 00043 00044 #include <cstdlib> // used by strtod() 00045 #include <cerrno> 00046 #include <cmath> 00047 #include <iostream> 00048 #include <vector> 00049 #include <algorithm> 00050 00051 //#define DODS_DEBUG 00052 00053 #include "BaseType.h" 00054 #include "Byte.h" 00055 #include "Int16.h" 00056 #include "UInt16.h" 00057 #include "Int32.h" 00058 #include "UInt32.h" 00059 #include "Float32.h" 00060 #include "Float64.h" 00061 #include "Str.h" 00062 #include "Url.h" 00063 #include "Array.h" 00064 #include "Structure.h" 00065 #include "Sequence.h" 00066 #include "Grid.h" 00067 #include "Error.h" 00068 #include "RValue.h" 00069 00070 #include "GSEClause.h" 00071 #include "GridGeoConstraint.h" 00072 #include "ArrayGeoConstraint.h" 00073 00074 #include "ce_functions.h" 00075 #include "gse_parser.h" 00076 #include "gse.tab.hh" 00077 #include "debug.h" 00078 #include "util.h" 00079 00080 // We wrapped VC++ 6.x strtod() to account for a short coming 00081 // in that function in regards to "NaN". I don't know if this 00082 // still applies in more recent versions of that product. 00083 // ROM - 12/2007 00084 #ifdef WIN32 00085 #include <limits> 00086 double w32strtod(const char *, char **); 00087 #endif 00088 00089 using namespace std; 00090 00091 int gse_parse(void *arg); 00092 void gse_restart(FILE * in); 00093 00094 // Glue routines declared in gse.lex 00095 void gse_switch_to_buffer(void *new_buffer); 00096 void gse_delete_buffer(void *buffer); 00097 void *gse_string(const char *yy_str); 00098 00099 namespace libdap { 00100 00102 inline bool double_eq(double lhs, double rhs, double epsilon = 1.0e-5) 00103 { 00104 if (lhs > rhs) 00105 return (lhs - rhs) < ((lhs + rhs) / epsilon); 00106 else 00107 return (rhs - lhs) < ((lhs + rhs) / epsilon); 00108 } 00109 00117 string extract_string_argument(BaseType * arg) 00118 { 00119 if (arg->type() != dods_str_c) 00120 throw Error(malformed_expr, 00121 "The function requires a DAP string argument."); 00122 00123 if (!arg->read_p()) 00124 throw InternalErr(__FILE__, __LINE__, 00125 "The CE Evaluator built an argument list where some constants held no values."); 00126 00127 string s = dynamic_cast<Str&>(*arg).value(); 00128 00129 DBG(cerr << "s: " << s << endl); 00130 00131 return s; 00132 } 00133 template<class T> static void set_array_using_double_helper(Array * a, 00134 double *src, int src_len) 00135 { 00136 T *values = new T[src_len]; 00137 for (int i = 0; i < src_len; ++i) 00138 values[i] = (T) src[i]; 00139 00140 #ifdef VAL2BUF 00141 a->val2buf(values, true); 00142 #else 00143 a->set_value(values, src_len); 00144 #endif 00145 00146 delete[]values; 00147 } 00148 00166 void set_array_using_double(Array * dest, double *src, int src_len) 00167 { 00168 // Simple types are Byte, ..., Float64, String and Url. 00169 if ((dest->type() == dods_array_c && !dest->var()->is_simple_type()) 00170 || dest->var()->type() == dods_str_c 00171 || dest->var()->type() == dods_url_c) 00172 throw InternalErr(__FILE__, __LINE__, 00173 "The function requires a DAP numeric-type array argument."); 00174 00175 // Test sizes. Note that Array::length() takes any constraint into account 00176 // when it returns the length. Even if this was removed, the 'helper' 00177 // function this uses calls Vector::val2buf() which uses Vector::width() 00178 // which in turn uses length(). 00179 if (dest->length() != src_len) 00180 throw InternalErr(__FILE__, __LINE__, 00181 "The source and destination array sizes don't match (" 00182 + long_to_string(src_len) + " versus " 00183 + long_to_string(dest->length()) + ")."); 00184 00185 // The types of arguments that the CE Parser will build for numeric 00186 // constants are limited to Uint32, Int32 and Float64. See ce_expr.y. 00187 // Expanded to work for any numeric type so it can be used for more than 00188 // just arguments. 00189 switch (dest->var()->type()) { 00190 case dods_byte_c: 00191 set_array_using_double_helper<dods_byte>(dest, src, src_len); 00192 break; 00193 case dods_uint16_c: 00194 set_array_using_double_helper<dods_uint16>(dest, src, src_len); 00195 break; 00196 case dods_int16_c: 00197 set_array_using_double_helper<dods_int16>(dest, src, src_len); 00198 break; 00199 case dods_uint32_c: 00200 set_array_using_double_helper<dods_uint32>(dest, src, src_len); 00201 break; 00202 case dods_int32_c: 00203 set_array_using_double_helper<dods_int32>(dest, src, src_len); 00204 break; 00205 case dods_float32_c: 00206 set_array_using_double_helper<dods_float32>(dest, src, src_len); 00207 break; 00208 case dods_float64_c: 00209 set_array_using_double_helper<dods_float64>(dest, src, src_len); 00210 break; 00211 default: 00212 throw InternalErr(__FILE__, __LINE__, 00213 "The argument list built by the CE parser contained an unsupported numeric type."); 00214 } 00215 00216 // Set the read_p property. 00217 dest->set_read_p(true); 00218 } 00219 00220 template<class T> static double *extract_double_array_helper(Array * a) 00221 { 00222 int length = a->length(); 00223 00224 T *b = new T[length]; 00225 a->value(b); 00226 00227 double *dest = new double[length]; 00228 for (int i = 0; i < length; ++i) 00229 dest[i] = (double) b[i]; 00230 delete[]b; 00231 00232 return dest; 00233 } 00234 00239 double *extract_double_array(Array * a) 00240 { 00241 // Simple types are Byte, ..., Float64, String and Url. 00242 if ((a->type() == dods_array_c && !a->var()->is_simple_type()) 00243 || a->var()->type() == dods_str_c || a->var()->type() == dods_url_c) 00244 throw Error(malformed_expr, 00245 "The function requires a DAP numeric-type array argument."); 00246 00247 if (!a->read_p()) 00248 throw InternalErr(__FILE__, __LINE__, 00249 string("The Array '") + a->name() + 00250 "'does not contain values."); 00251 00252 // The types of arguments that the CE Parser will build for numeric 00253 // constants are limited to Uint32, Int32 and Float64. See ce_expr.y. 00254 // Expanded to work for any numeric type so it can be used for more than 00255 // just arguments. 00256 switch (a->var()->type()) { 00257 case dods_byte_c: 00258 return extract_double_array_helper<dods_byte>(a); 00259 case dods_uint16_c: 00260 return extract_double_array_helper<dods_uint16>(a); 00261 case dods_int16_c: 00262 return extract_double_array_helper<dods_int16>(a); 00263 case dods_uint32_c: 00264 return extract_double_array_helper<dods_uint32>(a); 00265 case dods_int32_c: 00266 return extract_double_array_helper<dods_int32>(a); 00267 case dods_float32_c: 00268 return extract_double_array_helper<dods_float32>(a); 00269 case dods_float64_c: 00270 return extract_double_array_helper<dods_float64>(a); 00271 default: 00272 throw InternalErr(__FILE__, __LINE__, 00273 "The argument list built by the CE parser contained an unsupported numeric type."); 00274 } 00275 } 00276 00284 double extract_double_value(BaseType * arg) 00285 { 00286 // Simple types are Byte, ..., Float64, String and Url. 00287 if (!arg->is_simple_type() || arg->type() == dods_str_c || arg->type() 00288 == dods_url_c) 00289 throw Error(malformed_expr, 00290 "The function requires a DAP numeric-type argument."); 00291 00292 if (!arg->read_p()) 00293 throw InternalErr(__FILE__, __LINE__, 00294 "The CE Evaluator built an argument list where some constants held no values."); 00295 00296 // The types of arguments that the CE Parser will build for numeric 00297 // constants are limited to Uint32, Int32 and Float64. See ce_expr.y. 00298 // Expanded to work for any numeric type so it can be used for more than 00299 // just arguments. 00300 switch (arg->type()) { 00301 case dods_byte_c: 00302 return (double)(dynamic_cast<Byte&>(*arg).value()); 00303 case dods_uint16_c: 00304 return (double)(dynamic_cast<UInt16&>(*arg).value()); 00305 case dods_int16_c: 00306 return (double)(dynamic_cast<Int16&>(*arg).value()); 00307 case dods_uint32_c: 00308 return (double)(dynamic_cast<UInt32&>(*arg).value()); 00309 case dods_int32_c: 00310 return (double)(dynamic_cast<Int32&>(*arg).value()); 00311 case dods_float32_c: 00312 return (double)(dynamic_cast<Float32&>(*arg).value()); 00313 case dods_float64_c: 00314 return dynamic_cast<Float64&>(*arg).value(); 00315 default: 00316 throw InternalErr(__FILE__, __LINE__, 00317 "The argument list built by the CE parser contained an unsupported numeric type."); 00318 } 00319 } 00320 00323 void 00324 function_version(int, BaseType *[], DDS &, BaseType **btpp) 00325 { 00326 string 00327 xml_value = 00328 "<?xml version=\"1.0\" encoding=\"UTF-8\"?>\ 00329 <functions>\ 00330 <function name=\"geogrid\" version=\"1.2\"/>\ 00331 <function name=\"grid\" version=\"1.0\"/>\ 00332 <function name=\"linear_scale\" version=\"1.0b1\"/>\ 00333 <function name=\"version\" version=\"1.0\"/>\ 00334 </functions>"; 00335 00336 // <function name=\"geoarray\" version=\"0.9b1\"/> 00337 00338 Str *response = new Str("version"); 00339 00340 response->set_value(xml_value); 00341 *btpp = response; 00342 return; 00343 } 00344 00345 static void parse_gse_expression(gse_arg * arg, BaseType * expr) 00346 { 00347 gse_restart(0); // Restart the scanner. 00348 void *cls = gse_string(extract_string_argument(expr).c_str()); 00349 // gse_switch_to_buffer(cls); // Get set to scan the string. 00350 bool status = gse_parse((void *) arg) == 0; 00351 gse_delete_buffer(cls); 00352 if (!status) 00353 throw Error(malformed_expr, "Error parsing grid selection."); 00354 } 00355 00356 static void apply_grid_selection_expr(Grid * grid, GSEClause * clause) 00357 { 00358 // Basic plan: For each map, look at each clause and set start and stop 00359 // to be the intersection of the ranges in those clauses. 00360 Grid::Map_iter map_i = grid->map_begin(); 00361 while (map_i != grid->map_end() && (*map_i)->name() != clause->get_map_name()) 00362 ++map_i; 00363 00364 if (map_i == grid->map_end()) 00365 throw Error(malformed_expr,"The map vector '" + clause->get_map_name() 00366 + "' is not in the grid '" + grid->name() + "'."); 00367 00368 // Use pointer arith & the rule that map order must match array dim order 00369 Array::Dim_iter grid_dim = (grid->get_array()->dim_begin() + (map_i - grid->map_begin())); 00370 00371 Array *map = dynamic_cast < Array * >((*map_i)); 00372 if (!map) 00373 throw InternalErr(__FILE__, __LINE__, "Expected an Array"); 00374 int start = max(map->dimension_start(map->dim_begin()), clause->get_start()); 00375 int stop = min(map->dimension_stop(map->dim_begin()), clause->get_stop()); 00376 00377 if (start > stop) { 00378 ostringstream msg; 00379 msg 00380 << "The expressions passed to grid() do not result in an inclusive \n" 00381 << "subset of '" << clause->get_map_name() 00382 << "'. The map's values range " << "from " 00383 << clause->get_map_min_value() << " to " 00384 << clause->get_map_max_value() << "."; 00385 throw Error(malformed_expr,msg.str()); 00386 } 00387 00388 DBG(cerr << "Setting constraint on " << map->name() 00389 << "[" << start << ":" << stop << "]" << endl); 00390 00391 // Stride is always one. 00392 map->add_constraint(map->dim_begin(), start, 1, stop); 00393 grid->get_array()->add_constraint(grid_dim, start, 1, stop); 00394 } 00395 00396 static void apply_grid_selection_expressions(Grid * grid, 00397 vector < GSEClause * >clauses) 00398 { 00399 vector < GSEClause * >::iterator clause_i = clauses.begin(); 00400 while (clause_i != clauses.end()) 00401 apply_grid_selection_expr(grid, *clause_i++); 00402 00403 grid->set_read_p(false); 00404 } 00405 00442 void 00443 function_grid(int argc, BaseType * argv[], DDS &, BaseType **btpp) 00444 { 00445 DBG(cerr << "Entering function_grid..." << endl); 00446 00447 string info = 00448 string("<?xml version=\"1.0\" encoding=\"UTF-8\"?>\n") + 00449 "<function name=\"grid\" version=\"1.0\" href=\"http://docs.opendap.org/index.php/Server_Side_Processing_Functions#grid\">\n" + 00450 "</function>\n"; 00451 00452 if (argc == 0) { 00453 Str *response = new Str("info"); 00454 response->set_value(info); 00455 *btpp = response; 00456 return; 00457 } 00458 00459 Grid *original_grid = dynamic_cast < Grid * >(argv[0]); 00460 if (!original_grid) 00461 throw Error(malformed_expr,"The first argument to grid() must be a Grid variable!"); 00462 00463 // Duplicate the grid; DODSFilter::send_data() will delete the variable 00464 // after serializing it. 00465 Grid *l_grid = dynamic_cast < Grid * >(original_grid->ptr_duplicate()); 00466 if (!l_grid) 00467 throw InternalErr(__FILE__, __LINE__, "Expected a Grid."); 00468 00469 DBG(cerr << "grid: past initialization code" << endl); 00470 00471 // Read the maps. Do this before calling parse_gse_expression(). Avoid 00472 // reading the array until the constraints have been applied because it 00473 // might be really large. 00474 00475 // This version makes sure to set the send_p flags which is needed for 00476 // the hdf4 handler (and is what should be done in general). 00477 Grid::Map_iter i = l_grid->map_begin(); 00478 while (i != l_grid->map_end()) 00479 (*i++)->set_send_p(true); 00480 l_grid->read(); 00481 00482 DBG(cerr << "grid: past map read" << endl); 00483 00484 // argv[1..n] holds strings; each are little expressions to be parsed. 00485 // When each expression is parsed, the parser makes a new instance of 00486 // GSEClause. GSEClause checks to make sure the named map really exists 00487 // in the Grid and that the range of values given makes sense. 00488 vector < GSEClause * > clauses; 00489 gse_arg *arg = new gse_arg(l_grid); 00490 for (int i = 1; i < argc; ++i) { 00491 parse_gse_expression(arg, argv[i]); 00492 clauses.push_back(arg->get_gsec()); 00493 } 00494 delete arg; 00495 arg = 0; 00496 00497 apply_grid_selection_expressions(l_grid, clauses); 00498 00499 DBG(cerr << "grid: past gse application" << endl); 00500 00501 l_grid->get_array()->set_send_p(true); 00502 00503 l_grid->read(); 00504 00505 *btpp = l_grid; 00506 return; 00507 } 00508 00544 void 00545 function_geogrid(int argc, BaseType * argv[], DDS &, BaseType **btpp) 00546 { 00547 string info = 00548 string("<?xml version=\"1.0\" encoding=\"UTF-8\"?>\n") + 00549 "<function name=\"geogrid\" version=\"1.2\" href=\"http://docs.opendap.org/index.php/Server_Side_Processing_Functions#geogrid\">\n"+ 00550 "</function>"; 00551 00552 if (argc == 0) { 00553 Str *response = new Str("version"); 00554 response->set_value(info); 00555 *btpp = response; 00556 return ; 00557 } 00558 00559 // There are two main forms of this function, one that takes a Grid and one 00560 // that takes a Grid and two Arrays. The latter provides a way to explicitly 00561 // tell the function which maps contain lat and lon data. The remaining 00562 // arguments are the same for both versions, although that includes a 00563 // varying argument list. 00564 00565 // Look at the types of the first three arguments to determine which of the 00566 // two forms were used to call this function. 00567 Grid *l_grid = 0; 00568 if (argc < 1 || !(l_grid = dynamic_cast < Grid * >(argv[0]->ptr_duplicate()))) 00569 throw Error(malformed_expr,"The first argument to geogrid() must be a Grid variable!"); 00570 00571 // Both forms require at least this many args 00572 if (argc < 5) 00573 throw Error(malformed_expr,"Wrong number of arguments to geogrid() (expected at least 5 args). See geogrid() for more information."); 00574 00575 bool grid_lat_lon_form; 00576 Array *l_lat = 0; 00577 Array *l_lon = 0; 00578 if (!(l_lat = dynamic_cast < Array * >(argv[1]))) //->ptr_duplicate()))) 00579 grid_lat_lon_form = false; 00580 else if (!(l_lon = dynamic_cast < Array * >(argv[2]))) //->ptr_duplicate()))) 00581 throw Error(malformed_expr,"When using the Grid, Lat, Lon form of geogrid() both the lat and lon maps must be given (lon map missing)!"); 00582 else 00583 grid_lat_lon_form = true; 00584 00585 if (grid_lat_lon_form && argc < 7) 00586 throw Error(malformed_expr,"Wrong number of arguments to geogrid() (expected at least 7 args). See geogrid() for more information."); 00587 00588 #if 0 00589 Grid *l_grid = dynamic_cast < Grid * >(argv[0]->ptr_duplicate()); 00590 if (!l_grid) 00591 throw Error(malformed_expr,"The first argument to geogrid() must be a Grid variable!"); 00592 #endif 00593 // Read the maps. Do this before calling parse_gse_expression(). Avoid 00594 // reading the array until the constraints have been applied because it 00595 // might be really large. 00596 // 00597 // Trick: Some handlers build Grids from a combination of Array 00598 // variables and attributes. Those handlers (e.g., hdf4) use the send_p 00599 // property to determine which parts of the Grid to read *but they can 00600 // only read the maps from within Grid::read(), not the map's read()*. 00601 // Since the Grid's array does not have send_p set, it will not be read 00602 // by the call below to Grid::read(). 00603 Grid::Map_iter i = l_grid->map_begin(); 00604 while (i != l_grid->map_end()) 00605 (*i++)->set_send_p(true); 00606 00607 l_grid->read(); 00608 // Calling read() above sets the read_p flag for the entire grid; clear it 00609 // for the grid's array so that later on the code will be sure to read it 00610 // under all circumstances. 00611 l_grid->get_array()->set_read_p(false); 00612 DBG(cerr << "geogrid: past map read" << endl); 00613 00614 // Look for Grid Selection Expressions tacked onto the end of the BB 00615 // specification. If there are any, evaluate them before evaluating the BB. 00616 int min_arg_count = (grid_lat_lon_form) ? 7 : 5; 00617 if (argc > min_arg_count) { 00618 // argv[5..n] holds strings; each are little Grid Selection Expressions 00619 // to be parsed and evaluated. 00620 vector < GSEClause * > clauses; 00621 gse_arg *arg = new gse_arg(l_grid); 00622 for (int i = min_arg_count; i < argc; ++i) { 00623 parse_gse_expression(arg, argv[i]); 00624 clauses.push_back(arg->get_gsec()); 00625 } 00626 delete arg; 00627 arg = 0; 00628 00629 apply_grid_selection_expressions(l_grid, clauses); 00630 } 00631 00632 try { 00633 // Build a GeoConstraint object. If there are no longitude/latitude 00634 // maps then this constructor throws Error. 00635 GridGeoConstraint gc(l_grid); 00636 00637 // This sets the bounding box and modifies the maps to match the 00638 // notation of the box (0/359 or -180/179) 00639 int box_index_offset = (grid_lat_lon_form) ? 3 : 1; 00640 double top = extract_double_value(argv[box_index_offset]); 00641 double left = extract_double_value(argv[box_index_offset + 1]); 00642 double bottom = extract_double_value(argv[box_index_offset + 2]); 00643 double right = extract_double_value(argv[box_index_offset + 3]); 00644 gc.set_bounding_box(top, left, bottom, right); 00645 DBG(cerr << "geogrid: past bounding box set" << endl); 00646 00647 // This also reads all of the data into the grid variable 00648 gc.apply_constraint_to_data(); 00649 DBG(cerr << "geogrid: past apply constraint" << endl); 00650 00651 // In this function the l_grid pointer is the same as the pointer returned 00652 // by this call. The caller of the function must free the pointer. 00653 *btpp = gc.get_constrained_grid(); 00654 return; 00655 } 00656 catch (Error &e) { 00657 throw e; 00658 } 00659 catch (exception & e) { 00660 throw 00661 InternalErr(string 00662 ("A C++ exception was thrown from inside geogrid(): ") 00663 + e.what()); 00664 } 00665 } 00666 00667 // These static functions could be moved to a class that provides a more 00668 // general interface for COARDS/CF someday. Assume each BaseType comes bundled 00669 // with an attribute table. 00670 00671 // This was ripped from parser-util.cc 00672 static double string_to_double(const char *val) 00673 { 00674 char *ptr; 00675 errno = 0; 00676 // Clear previous value. 5/21/2001 jhrg 00677 00678 #ifdef WIN32 00679 double v = w32strtod(val, &ptr); 00680 #else 00681 double v = strtod(val, &ptr); 00682 #endif 00683 00684 if ((v == 0.0 && (val == ptr || errno == HUGE_VAL || errno == ERANGE)) 00685 || *ptr != '\0') { 00686 throw Error(malformed_expr,string("Could not convert the string '") + val + "' to a double."); 00687 } 00688 00689 double abs_val = fabs(v); 00690 if (abs_val > DODS_DBL_MAX || (abs_val != 0.0 && abs_val < DODS_DBL_MIN)) 00691 throw Error(malformed_expr,string("Could not convert the string '") + val + "' to a double."); 00692 00693 return v; 00694 } 00695 00705 static double get_attribute_double_value(BaseType *var, 00706 vector<string> &attributes) 00707 { 00708 // This code also builds a list of the attribute values that have been 00709 // passed in but not found so that an informative message can be returned. 00710 AttrTable &attr = var->get_attr_table(); 00711 string attribute_value = ""; 00712 string values = ""; 00713 vector<string>::iterator i = attributes.begin(); 00714 while (attribute_value == "" && i != attributes.end()) { 00715 values += *i; 00716 if (!values.empty()) 00717 values += ", "; 00718 attribute_value = attr.get_attr(*i++); 00719 } 00720 00721 // If the value string is empty, then look at the grid's array (if it's a 00722 // grid) or throw an Error. 00723 if (attribute_value.empty()) { 00724 if (var->type() == dods_grid_c) 00725 return get_attribute_double_value(dynamic_cast<Grid&>(*var).get_array(), attributes); 00726 else 00727 throw Error(malformed_expr,string("No COARDS/CF '") + values.substr(0, values.length() - 2) 00728 + "' attribute was found for the variable '" 00729 + var->name() + "'."); 00730 } 00731 00732 return string_to_double(remove_quotes(attribute_value).c_str()); 00733 } 00734 00735 static double get_attribute_double_value(BaseType *var, const string &attribute) 00736 { 00737 AttrTable &attr = var->get_attr_table(); 00738 string attribute_value = attr.get_attr(attribute); 00739 00740 // If the value string is empty, then look at the grid's array (if it's a 00741 // grid or throw an Error. 00742 if (attribute_value.empty()) { 00743 if (var->type() == dods_grid_c) 00744 return get_attribute_double_value(dynamic_cast<Grid&>(*var).get_array(), attribute); 00745 else 00746 throw Error(malformed_expr,string("No COARDS '") + attribute 00747 + "' attribute was found for the variable '" 00748 + var->name() + "'."); 00749 } 00750 00751 return string_to_double(remove_quotes(attribute_value).c_str()); 00752 } 00753 00754 static double get_y_intercept(BaseType *var) 00755 { 00756 vector<string> attributes; 00757 attributes.push_back("add_offset"); 00758 attributes.push_back("add_off"); 00759 return get_attribute_double_value(var, attributes); 00760 } 00761 00762 static double get_slope(BaseType *var) 00763 { 00764 return get_attribute_double_value(var, "scale_factor"); 00765 } 00766 00767 static double get_missing_value(BaseType *var) 00768 { 00769 return get_attribute_double_value(var, "missing_value"); 00770 } 00771 00784 void 00785 function_linear_scale(int argc, BaseType * argv[], DDS &, BaseType **btpp) 00786 { 00787 string info = 00788 string("<?xml version=\"1.0\" encoding=\"UTF-8\"?>\n") + 00789 "<function name=\"linear_scale\" version=\"1.0b1\" href=\"http://docs.opendap.org/index.php/Server_Side_Processing_Functions#linear_scale\">\n" + 00790 "</function>"; 00791 00792 if (argc == 0) { 00793 Str *response = new Str("info"); 00794 response->set_value(info); 00795 *btpp = response; 00796 return; 00797 } 00798 00799 // Check for 1 or 3 arguments: 1 --> use attributes; 3 --> m & b supplied 00800 DBG(cerr << "argc = " << argc << endl); 00801 if (!(argc == 1 || argc == 3 || argc == 4)) 00802 throw Error(malformed_expr,"Wrong number of arguments to linear_scale(). See linear_scale() for more information"); 00803 00804 // Get m & b 00805 bool use_missing = false; 00806 double m, b, missing = 0.0; 00807 if (argc == 4) { 00808 m = extract_double_value(argv[1]); 00809 b = extract_double_value(argv[2]); 00810 missing = extract_double_value(argv[3]); 00811 use_missing = true; 00812 } else if (argc == 3) { 00813 m = extract_double_value(argv[1]); 00814 b = extract_double_value(argv[2]); 00815 use_missing = false; 00816 } else { 00817 m = get_slope(argv[0]); 00818 00819 // This is really a hack; on a fair number of datasets, the y intercept 00820 // is not given and is assumed to be 0. Here the function looks and 00821 // catches the error if a y intercept is not found. 00822 try { 00823 b = get_y_intercept(argv[0]); 00824 } 00825 catch (Error &e) { 00826 b = 0.0; 00827 } 00828 00829 // This is not the best plan; the get_missing_value() function should 00830 // do something other than throw, but to do that would require mayor 00831 // surgery on get_attribute_double_value(). 00832 try { 00833 missing = get_missing_value(argv[0]); 00834 use_missing = true; 00835 } 00836 catch (Error &e) { 00837 use_missing = false; 00838 } 00839 } 00840 00841 DBG(cerr << "m: " << m << ", b: " << b << endl);DBG(cerr << "use_missing: " << use_missing << ", missing: " << missing << endl); 00842 00843 // Read the data, scale and return the result. Must replace the new data 00844 // in a constructor (i.e., Array part of a Grid). 00845 BaseType *dest = 0; 00846 double *data; 00847 if (argv[0]->type() == dods_grid_c) { 00848 Array &source = *dynamic_cast<Grid&>(*argv[0]).get_array(); 00849 source.set_send_p(true); 00850 source.read(); 00851 data = extract_double_array(&source); 00852 int length = source.length(); 00853 int i = 0; 00854 while (i < length) { 00855 DBG2(cerr << "data[" << i << "]: " << data[i] << endl); 00856 if (!use_missing || !double_eq(data[i], missing)) 00857 data[i] = data[i] * m + b; 00858 DBG2(cerr << " >> data[" << i << "]: " << data[i] << endl); 00859 ++i; 00860 } 00861 00862 // Vector::add_var will delete the existing 'template' variable 00863 Float64 *temp_f = new Float64(source.name()); 00864 source.add_var(temp_f); 00865 #ifdef VAL2BUF 00866 source.val2buf(static_cast<void*>(data), false); 00867 #else 00868 source.set_value(data, i); 00869 #endif 00870 delete [] data; // val2buf copies. 00871 delete temp_f; // add_var copies and then adds. 00872 dest = argv[0]; 00873 } else if (argv[0]->is_vector_type()) { 00874 Array &source = dynamic_cast<Array&>(*argv[0]); 00875 source.set_send_p(true); 00876 // If the array is really a map, make sure to read using the Grid 00877 // because of the HDF4 handler's odd behavior WRT dimensions. 00878 if (source.get_parent() && source.get_parent()->type() == dods_grid_c) 00879 source.get_parent()->read(); 00880 else 00881 source.read(); 00882 00883 data = extract_double_array(&source); 00884 int length = source.length(); 00885 int i = 0; 00886 while (i < length) { 00887 if (!use_missing || !double_eq(data[i], missing)) 00888 data[i] = data[i] * m + b; 00889 ++i; 00890 } 00891 00892 Float64 *temp_f = new Float64(source.name()); 00893 source.add_var(temp_f); 00894 00895 source.val2buf(static_cast<void*>(data), false); 00896 00897 delete [] data; // val2buf copies. 00898 delete temp_f; // add_var copies and then adds. 00899 00900 dest = argv[0]; 00901 } else if (argv[0]->is_simple_type() && !(argv[0]->type() == dods_str_c 00902 || argv[0]->type() == dods_url_c)) { 00903 double data = extract_double_value(argv[0]); 00904 if (!use_missing || !double_eq(data, missing)) 00905 data = data * m + b; 00906 00907 dest = new Float64(argv[0]->name()); 00908 00909 dest->val2buf(static_cast<void*>(&data)); 00910 00911 } else { 00912 throw Error(malformed_expr,"The linear_scale() function works only for numeric Grids, Arrays and scalars."); 00913 } 00914 00915 *btpp = dest; 00916 return; 00917 } 00918 00919 #if 0 00920 00936 void 00937 function_geoarray(int argc, BaseType * argv[], DDS &, BaseType **btpp) 00938 { 00939 string info = 00940 string("<?xml version=\"1.0\" encoding=\"UTF-8\"?>\n") + 00941 "<function name=\"geoarray\" version=\"0.9b1\" href=\"http://docs.opendap.org/index.php/Server_Side_Processing_Functions#geoarray\">\n" + 00942 "</function>"; 00943 00944 if (argc == 0) { 00945 Str *response = new Str("version"); 00946 response->set_value(info); 00947 *btpp = response; 00948 return; 00949 } 00950 00951 DBG(cerr << "argc = " << argc << endl); 00952 if (!(argc == 5 || argc == 9 || argc == 11)) 00953 throw Error(malformed_expr,"Wrong number of arguments to geoarray(). See geoarray() for more information."); 00954 00955 // Check the Array (and dup because the caller will free the variable). 00956 Array *l_array = dynamic_cast < Array * >(argv[0]->ptr_duplicate()); 00957 if (!l_array) 00958 throw Error(malformed_expr,"The first argument to geoarray() must be an Array variable!"); 00959 00960 try { 00961 00962 // Read the bounding box and variable extents from the params 00963 double bb_top = extract_double_value(argv[1]); 00964 double bb_left = extract_double_value(argv[2]); 00965 double bb_bottom = extract_double_value(argv[3]); 00966 double bb_right = extract_double_value(argv[4]); 00967 00968 switch (argc) { 00969 case 5: { 00970 ArrayGeoConstraint agc(l_array); 00971 00972 agc.set_bounding_box(bb_left, bb_top, bb_right, bb_bottom); 00973 // This also reads all of the data into the grid variable 00974 agc.apply_constraint_to_data(); 00975 *btpp = agc.get_constrained_array(); 00976 return; 00977 break; 00978 } 00979 case 9: { 00980 double var_top = extract_double_value(argv[5]); 00981 double var_left = extract_double_value(argv[6]); 00982 double var_bottom = extract_double_value(argv[7]); 00983 double var_right = extract_double_value(argv[8]); 00984 ArrayGeoConstraint agc (l_array, var_left, var_top, var_right, var_bottom); 00985 00986 agc.set_bounding_box(bb_left, bb_top, bb_right, bb_bottom); 00987 // This also reads all of the data into the grid variable 00988 agc.apply_constraint_to_data(); 00989 *btpp = agc.get_constrained_array(); 00990 return; 00991 break; 00992 } 00993 case 11: { 00994 double var_top = extract_double_value(argv[5]); 00995 double var_left = extract_double_value(argv[6]); 00996 double var_bottom = extract_double_value(argv[7]); 00997 double var_right = extract_double_value(argv[8]); 00998 string projection = extract_string_argument(argv[9]); 00999 string datum = extract_string_argument(argv[10]); 01000 ArrayGeoConstraint agc(l_array, 01001 var_left, var_top, var_right, var_bottom, 01002 projection, datum); 01003 01004 agc.set_bounding_box(bb_left, bb_top, bb_right, bb_bottom); 01005 // This also reads all of the data into the grid variable 01006 agc.apply_constraint_to_data(); 01007 *btpp = agc.get_constrained_array(); 01008 return; 01009 break; 01010 } 01011 default: 01012 throw InternalErr(__FILE__, __LINE__, "Wrong number of args to geoarray."); 01013 } 01014 } 01015 catch (Error & e) { 01016 throw e; 01017 } 01018 catch (exception & e) { 01019 throw 01020 InternalErr(string 01021 ("A C++ exception was thrown from inside geoarray(): ") 01022 + e.what()); 01023 01024 } 01025 01026 throw InternalErr(__FILE__, __LINE__, "Impossible condition in geoarray."); 01027 } 01028 #endif 01029 01030 void register_functions(ConstraintEvaluator & ce) 01031 { 01032 ce.add_function("grid", function_grid); 01033 ce.add_function("geogrid", function_geogrid); 01034 ce.add_function("linear_scale", function_linear_scale); 01035 #if 0 01036 ce.add_function("geoarray", function_geoarray); 01037 #endif 01038 ce.add_function("version", function_version); 01039 } 01040 01041 } // namespace libdap