00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036 #include "config.h"
00037
00038 static char rcsid[]not_used =
00039 { "$Id: ce_functions.cc 17299 2007-10-24 19:39:50Z jimg $"
00040 };
00041
00042 #include <errno.h>
00043 #include <limits.h>
00044 #include <math.h>
00045
00046 #include <iostream>
00047 #include <vector>
00048 #include <algorithm>
00049
00050
00051
00052 #include "BaseType.h"
00053 #include "Array.h"
00054 #include "Sequence.h"
00055 #include "Grid.h"
00056 #include "Error.h"
00057 #include "RValue.h"
00058
00059 #include "GSEClause.h"
00060 #include "GridGeoConstraint.h"
00061 #include "ArrayGeoConstraint.h"
00062
00063 #include "ce_functions.h"
00064 #include "gse_parser.h"
00065 #include "gse.tab.h"
00066 #include "debug.h"
00067 #include "util.h"
00068
00069
00070
00071
00072
00073 #ifdef WIN32
00074 #include <limits>
00075 double w32strtod(const char *, char **);
00076 #endif
00077
00078 using namespace std;
00079
00080 int gse_parse(void *arg);
00081 void gse_restart(FILE * in);
00082
00083
00084 void gse_switch_to_buffer(void *new_buffer);
00085 void gse_delete_buffer(void *buffer);
00086 void *gse_string(const char *yy_str);
00087
00088 namespace libdap {
00089
00091 inline bool double_eq(double lhs, double rhs, double epsilon = 1.0e-5)
00092 {
00093 if (lhs > rhs)
00094 return (lhs - rhs) < ((lhs + rhs) / epsilon);
00095 else
00096 return (rhs - lhs) < ((lhs + rhs) / epsilon);
00097 }
00098
00106 string extract_string_argument(BaseType * arg)
00107 {
00108 if (arg->type() != dods_str_c)
00109 throw Error(malformed_expr,
00110 "The function requires a DAP string argument.");
00111
00112 if (!arg->read_p())
00113 throw InternalErr(__FILE__, __LINE__,
00114 "The CE Evaluator built an argument list where some constants held no values.");
00115
00116 string s = dynamic_cast<Str&>(*arg).value();
00117
00118 DBG(cerr << "s: " << s << endl);
00119
00120 return s;
00121 }
00122 template<class T> static void set_array_using_double_helper(Array * a,
00123 double *src, int src_len)
00124 {
00125 T *values = new T[src_len];
00126 for (int i = 0; i < src_len; ++i)
00127 values[i] = (T) src[i];
00128 a->val2buf(values, true);
00129 delete[]values;
00130 }
00131
00149 void set_array_using_double(Array * dest, double *src, int src_len)
00150 {
00151
00152 if (dest->type() == dods_array_c && !dest->var()->is_simple_type() || dest->var()->type() == dods_str_c || dest->var()->type() == dods_url_c)
00153 throw InternalErr(__FILE__, __LINE__,
00154 "The function requires a DAP numeric-type array argument.");
00155
00156
00157
00158
00159
00160 if (dest->length() != src_len)
00161 throw InternalErr(__FILE__, __LINE__,
00162 "The source and destination array sizes don't match ("
00163 + long_to_string(src_len) + " versus "
00164 + long_to_string(dest->length()) + ").");
00165
00166
00167
00168
00169
00170 switch (dest->var()->type()) {
00171 case dods_byte_c:
00172 set_array_using_double_helper<dods_byte>(dest, src, src_len);
00173 break;
00174 case dods_uint16_c:
00175 set_array_using_double_helper<dods_uint16>(dest, src, src_len);
00176 break;
00177 case dods_int16_c:
00178 set_array_using_double_helper<dods_int16>(dest, src, src_len);
00179 break;
00180 case dods_uint32_c:
00181 set_array_using_double_helper<dods_uint32>(dest, src, src_len);
00182 break;
00183 case dods_int32_c:
00184 set_array_using_double_helper<dods_int32>(dest, src, src_len);
00185 break;
00186 case dods_float32_c:
00187 set_array_using_double_helper<dods_float32>(dest, src, src_len);
00188 break;
00189 case dods_float64_c:
00190 set_array_using_double_helper<dods_float64>(dest, src, src_len);
00191 break;
00192 default:
00193 throw InternalErr(__FILE__, __LINE__,
00194 "The argument list built by the CE parser contained an unsupported numeric type.");
00195 }
00196
00197
00198 dest->set_read_p(true);
00199 }
00200
00201 template<class T> static double *extract_double_array_helper(Array * a)
00202 {
00203 int length = a->length();
00204
00205 T *b = new T[length];
00206 a->value(b);
00207
00208 double *dest = new double[length];
00209 for (int i = 0; i < length; ++i)
00210 dest[i] = (double) b[i];
00211 delete[]b;
00212
00213 return dest;
00214 }
00215
00220 double *extract_double_array(Array * a)
00221 {
00222
00223 if (a->type() == dods_array_c && !a->var()->is_simple_type() || a->var()->type() == dods_str_c || a->var()->type() == dods_url_c)
00224 throw Error(malformed_expr,
00225 "The function requires a DAP numeric-type array argument.");
00226
00227 if (!a->read_p())
00228 throw InternalErr(__FILE__, __LINE__,
00229 string("The Array '") + a->name() +
00230 "'does not contain values.");
00231
00232
00233
00234
00235
00236 switch (a->var()->type()) {
00237 case dods_byte_c:
00238 return extract_double_array_helper<dods_byte>(a);
00239 case dods_uint16_c:
00240 return extract_double_array_helper<dods_uint16>(a);
00241 case dods_int16_c:
00242 return extract_double_array_helper<dods_int16>(a);
00243 case dods_uint32_c:
00244 return extract_double_array_helper<dods_uint32>(a);
00245 case dods_int32_c:
00246 return extract_double_array_helper<dods_int32>(a);
00247 case dods_float32_c:
00248 return extract_double_array_helper<dods_float32>(a);
00249 case dods_float64_c:
00250 return extract_double_array_helper<dods_float64>(a);
00251 default:
00252 throw InternalErr(__FILE__, __LINE__,
00253 "The argument list built by the CE parser contained an unsupported numeric type.");
00254 }
00255 }
00256
00264 double extract_double_value(BaseType * arg)
00265 {
00266
00267 if (!arg->is_simple_type() || arg->type() == dods_str_c || arg->type()
00268 == dods_url_c)
00269 throw Error(malformed_expr,
00270 "The function requires a DAP numeric-type argument.");
00271
00272 if (!arg->read_p())
00273 throw InternalErr(__FILE__, __LINE__,
00274 "The CE Evaluator built an argument list where some constants held no values.");
00275
00276
00277
00278
00279
00280 switch (arg->type()) {
00281 case dods_byte_c:
00282 return (double)(dynamic_cast<Byte&>(*arg).value());
00283 case dods_uint16_c:
00284 return (double)(dynamic_cast<UInt16&>(*arg).value());
00285 case dods_int16_c:
00286 return (double)(dynamic_cast<Int16&>(*arg).value());
00287 case dods_uint32_c:
00288 return (double)(dynamic_cast<UInt32&>(*arg).value());
00289 case dods_int32_c:
00290 return (double)(dynamic_cast<Int32&>(*arg).value());
00291 case dods_float32_c:
00292 return (double)(dynamic_cast<Float32&>(*arg).value());
00293 case dods_float64_c:
00294 return dynamic_cast<Float64&>(*arg).value();
00295 default:
00296 throw InternalErr(__FILE__, __LINE__,
00297 "The argument list built by the CE parser contained an unsupported numeric type.");
00298 }
00299 }
00300
00303 BaseType *function_version(int, BaseType *[], DDS &, const string &)
00304 {
00305 string
00306 xml_value =
00307 "<?xml version=\"1.0\" encoding=\"UTF-8\"?>\
00308 <functions>\
00309 <function name=\"version\" version=\"1.0\"/>\
00310 <function name=\"grid\" version=\"1.0\"/>\
00311 <function name=\"geogrid\" version=\"1.0b2\"/>\
00312 <function name=\"geoarray\" version=\"0.9b1\"/>\
00313 <function name=\"linear_scale\" version=\"1.0b1\"/>\
00314 </functions>";
00315
00316 Str *response = new Str("version");
00317
00318 response->set_value(xml_value);
00319
00320 return response;
00321 }
00322
00323 static void parse_gse_expression(gse_arg * arg, BaseType * expr)
00324 {
00325 gse_restart(0);
00326 void *cls = gse_string(extract_string_argument(expr).c_str());
00327
00328 bool status = gse_parse((void *) arg) == 0;
00329 gse_delete_buffer(cls);
00330 if (!status)
00331 throw Error(malformed_expr, "Error parsing grid selection.");
00332 }
00333
00334 static void apply_grid_selection_expr(Grid * grid, GSEClause * clause)
00335 {
00336
00337
00338 Grid::Map_iter map_i = grid->map_begin();
00339 while (map_i != grid->map_end() && (*map_i)->name() != clause->get_map_name())
00340 ++map_i;
00341
00342 if (map_i == grid->map_end())
00343 throw Error("The map vector '" + clause->get_map_name()
00344 + "' is not in the grid '" + grid->name() + "'.");
00345
00346
00347 Array::Dim_iter grid_dim = (grid->get_array()->dim_begin() + (map_i - grid->map_begin()));
00348
00349 Array *map = dynamic_cast < Array * >((*map_i));
00350 int start = max(map->dimension_start(map->dim_begin()), clause->get_start());
00351 int stop = min(map->dimension_stop(map->dim_begin()), clause->get_stop());
00352
00353 if (start > stop) {
00354 ostringstream msg;
00355 msg
00356 << "The expressions passed to grid() do not result in an inclusive \n"
00357 << "subset of '" << clause->get_map_name()
00358 << "'. The map's values range " << "from "
00359 << clause->get_map_min_value() << " to "
00360 << clause->get_map_max_value() << ".";
00361 throw Error(msg.str());
00362 }
00363
00364 DBG(cerr << "Setting constraint on " << map->name()
00365 << "[" << start << ":" << stop << "]" << endl);
00366
00367
00368 map->add_constraint(map->dim_begin(), start, 1, stop);
00369 grid->get_array()->add_constraint(grid_dim, start, 1, stop);
00370 }
00371
00372 static void apply_grid_selection_expressions(Grid * grid,
00373 vector < GSEClause * >clauses)
00374 {
00375 vector < GSEClause * >::iterator clause_i = clauses.begin();
00376 while (clause_i != clauses.end())
00377 apply_grid_selection_expr(grid, *clause_i++);
00378
00379 grid->set_read_p(false);
00380 }
00381
00418 BaseType *function_grid(int argc, BaseType * argv[], DDS &,
00419 const string & dataset)
00420 {
00421 DBG(cerr << "Entering function_grid..." << endl);
00422
00423 string
00424 info =
00425 "<?xml version=\"1.0\" encoding=\"UTF-8\"?>\
00426 <function name=\"grid\" version=\"1.0\"/>\
00427 The grid() function takes a grid variable and zero or more relational\
00428 expressions. Each relational expression is applied to the grid using the\
00429 server's constraint evaluator and the resulting grid is returned. The\
00430 expressions may use constants and the grid's map vectors but may not use\
00431 any other variables. Two forms of expression are provide: \"var relop const\"\
00432 and \"const relop var relop const\". For example: grid(sst, \"10<=TIME<20\")\
00433 and grid(sst, \"10<=TIME\", \"TIME<20\") are both legal and, in this case,\
00434 also equivalent.\
00435 </function>";
00436
00437 if (argc == 0) {
00438 Str *response = new Str("info");
00439 response->set_value(info);
00440 return response;
00441 }
00442
00443 Grid *original_grid = dynamic_cast < Grid * >(argv[0]);
00444 if (!original_grid)
00445 throw Error("The first argument to grid() must be a Grid variable!");
00446
00447
00448
00449 Grid *l_grid = dynamic_cast < Grid * >(original_grid->ptr_duplicate());
00450
00451 DBG(cerr << "grid: past initialization code" << endl);
00452
00453
00454
00455
00456
00457
00458
00459 Grid::Map_iter i = l_grid->map_begin();
00460 while (i != l_grid->map_end())
00461 (*i++)->set_send_p(true);
00462 l_grid->read(dataset);
00463
00464 DBG(cerr << "grid: past map read" << endl);
00465
00466
00467
00468
00469
00470 vector < GSEClause * > clauses;
00471 gse_arg *arg = new gse_arg(l_grid);
00472 for (int i = 1; i < argc; ++i) {
00473 parse_gse_expression(arg, argv[i]);
00474 clauses.push_back(arg->get_gsec());
00475 }
00476 delete arg;
00477 arg = 0;
00478
00479 apply_grid_selection_expressions(l_grid, clauses);
00480
00481 DBG(cerr << "grid: past gse application" << endl);
00482
00483 l_grid->get_array()->set_send_p(true);
00484
00485 l_grid->read(dataset);
00486
00487 return l_grid;
00488 }
00489
00525 BaseType *function_geogrid(int argc, BaseType * argv[], DDS &,
00526 const string & dataset)
00527 {
00528 string
00529 info =
00530 "<?xml version=\"1.0\" encoding=\"UTF-8\"?>\
00531 <function name=\"geogrid\" version=\"1.0b2\"/>\
00532 geogrid() applies a constraint given in latitude and longitude to a\
00533 DAP Grid variable. The arguments to the function are:\
00534 geogrid(<grid variable>, <upper latitude>, <left longitude>,\
00535 <lower latitude>, <right longitude> [selection expressions - see grid()])\
00536 geogrid(\"version\") returns the version of the function.\
00537 The function will always return a single Grid variable whose values\
00538 completely cover the given region, although there may be cases when\
00539 some additional data is also returned. If the longitude values 'wrap\
00540 around' the right edge of the data, then the function will make two\
00541 requests and return those joined together as a single Grid.\
00542 </function>";
00543
00544 if (argc == 0) {
00545 Str *response = new Str("version");
00546 response->set_value(info);
00547 return response;
00548 }
00549
00550 if (argc < 5)
00551 throw Error("Wrong number of arguments to geogrid(). See geogrid() for more information.");
00552
00553 Grid *l_grid = dynamic_cast < Grid * >(argv[0]->ptr_duplicate());
00554 if (!l_grid)
00555 throw Error("The first argument to geogrid() must be a Grid variable!");
00556
00557
00558
00559
00560
00561
00562
00563
00564
00565
00566
00567 Grid::Map_iter i = l_grid->map_begin();
00568 while (i != l_grid->map_end())
00569 (*i++)->set_send_p(true);
00570 l_grid->read(dataset);
00571
00572
00573
00574 l_grid->get_array()->set_read_p(false);DBG(cerr << "geogrid: past map read" << endl);
00575
00576
00577
00578 if (argc > 5) {
00579
00580
00581 vector < GSEClause * > clauses;
00582 gse_arg *arg = new gse_arg(l_grid);
00583 for (int i = 5; i < argc; ++i) {
00584 parse_gse_expression(arg, argv[i]);
00585 clauses.push_back(arg->get_gsec());
00586 }
00587 delete arg;
00588 arg = 0;
00589
00590 apply_grid_selection_expressions(l_grid, clauses);
00591 }
00592
00593 try {
00594
00595
00596 GridGeoConstraint gc(l_grid, dataset);
00597
00598
00599
00600 double top = extract_double_value(argv[1]);
00601 double left = extract_double_value(argv[2]);
00602 double bottom = extract_double_value(argv[3]);
00603 double right = extract_double_value(argv[4]);
00604 gc.set_bounding_box(left, top, right, bottom);
00605 DBG(cerr << "geogrid: past bounding box set" << endl);
00606
00607
00608 gc.apply_constraint_to_data();
00609 DBG(cerr << "geogrid: past apply constraint" << endl);
00610
00611
00612
00613 return gc.get_constrained_grid();
00614 }
00615 catch (Error &e) {
00616 throw e;
00617 }
00618 catch (exception & e) {
00619 throw
00620 InternalErr(string
00621 ("A C++ exception was thrown from inside geogrid(): ")
00622 + e.what());
00623 }
00624 }
00625
00626
00627
00628
00629
00630
00631 static double string_to_double(const char *val)
00632 {
00633 char *ptr;
00634 errno = 0;
00635
00636
00637 #ifdef WIN32
00638 double v = w32strtod(val, &ptr);
00639 #else
00640 double v = strtod(val, &ptr);
00641 #endif
00642
00643 if ((v == 0.0 && (val == ptr || errno == HUGE_VAL || errno == ERANGE))
00644 || *ptr != '\0') {
00645 throw Error(string("Could not convert the string '") + val + "' to a double.");
00646 }
00647
00648 double abs_val = fabs(v);
00649 if (abs_val > DODS_DBL_MAX || (abs_val != 0.0 && abs_val < DODS_DBL_MIN))
00650 throw Error(string("Could not convert the string '") + val + "' to a double.");
00651
00652 return v;
00653 }
00654
00655 static string remove_quotes(const string &s)
00656 {
00657 if (s[0] == '\"' && s[s.length()-1] == '\"')
00658 return s.substr(1, s.length() - 2);
00659 else
00660 return s;
00661 }
00662
00666 static double get_attribute_double_value(BaseType *var,
00667 vector<string> &attributes)
00668 {
00669 AttrTable &attr = var->get_attr_table();
00670 string attribute_value = "";
00671 string values = "";
00672 vector<string>::iterator i = attributes.begin();
00673 while (attribute_value == "" && i != attributes.end()) {
00674 values += *i;
00675 if (!values.empty())
00676 values += ", ";
00677 attribute_value = attr.get_attr(*i++);
00678 }
00679
00680
00681
00682 if (attribute_value.empty()) {
00683 if (var->type() == dods_grid_c)
00684 return get_attribute_double_value(dynamic_cast<Grid&>(*var).get_array(), attributes);
00685 else
00686 throw Error(string("No COARDS '") + values.substr(0, values.length() - 2)
00687 + "' attribute was found for the variable '"
00688 + var->name() + "'.");
00689 }
00690
00691 return string_to_double(remove_quotes(attribute_value).c_str());
00692 }
00693
00694 static double get_attribute_double_value(BaseType *var, const string &attribute)
00695 {
00696 AttrTable &attr = var->get_attr_table();
00697 string attribute_value = attr.get_attr(attribute);
00698
00699
00700
00701 if (attribute_value.empty()) {
00702 if (var->type() == dods_grid_c)
00703 return get_attribute_double_value(dynamic_cast<Grid&>(*var).get_array(), attribute);
00704 else
00705 throw Error(string("No COARDS '") + attribute
00706 + "' attribute was found for the variable '"
00707 + var->name() + "'.");
00708 }
00709
00710 return string_to_double(remove_quotes(attribute_value).c_str());
00711 }
00712
00713 static double get_y_intercept(BaseType *var)
00714 {
00715 vector<string> attributes;
00716 attributes.push_back("add_offset");
00717 attributes.push_back("add_off");
00718 return get_attribute_double_value(var, attributes);
00719 }
00720
00721 static double get_slope(BaseType *var)
00722 {
00723 return get_attribute_double_value(var, "scale_factor");
00724 }
00725
00726 static double get_missing_value(BaseType *var)
00727 {
00728 return get_attribute_double_value(var, "missing_value");
00729 }
00730
00743 BaseType * function_linear_scale(int argc, BaseType * argv[], DDS &,
00744 const string & dataset)
00745 {
00746 string
00747 info =
00748 "<?xml version=\"1.0\" encoding=\"UTF-8\"?>\
00749 <function name=\"linear_scale\" version=\"1.0b1\">\
00750 The linear_scale() function applies the familiar y=mx+b equation to data.\
00751 If only the name of a variable is given, the function looks for the COARDS\
00752 scale_factor, add_offset and missing_value attributes. In the equation,\
00753 'm' is scale_factor, 'b' is add_offset and data values which match\
00754 missing_value are not scaled. If add_offset cannot be found, it defaults to\
00755 zero; if missing_value cannot be found, the test for it is not performed.\
00756 The caller can also provide values for these and, if provided, those values\
00757 are used and the dataset's attributes are ignored. If only scale_factor is\
00758 provided, the defaults for 'b' and the missing value flag are used.\
00759 Similarly, if only 'm' and 'b' are provided the missing value test is not\
00760 performed. The values are always returned in a Float64 array.\
00761 </function>";
00762
00763 if (argc == 0) {
00764 Str *response = new Str("info");
00765 response->set_value(info);
00766 return response;
00767 }
00768
00769
00770 DBG(cerr << "argc = " << argc << endl);
00771 if (!(argc == 1 || argc == 3 || argc == 4))
00772 throw Error("Wrong number of arguments to linear_scale(). See linear_scale() for more information");
00773
00774
00775 bool use_missing = false;
00776 double m, b, missing = 0.0;
00777 if (argc == 4) {
00778 m = extract_double_value(argv[1]);
00779 b = extract_double_value(argv[2]);
00780 missing = extract_double_value(argv[3]);
00781 use_missing = true;
00782 } else if (argc == 3) {
00783 m = extract_double_value(argv[1]);
00784 b = extract_double_value(argv[2]);
00785 use_missing = false;
00786 } else {
00787 m = get_slope(argv[0]);
00788
00789
00790
00791
00792 try {
00793 b = get_y_intercept(argv[0]);
00794 }
00795 catch (Error &e) {
00796 b = 0.0;
00797 }
00798
00799
00800
00801
00802 try {
00803 missing = get_missing_value(argv[0]);
00804 use_missing = true;
00805 }
00806 catch (Error &e) {
00807 use_missing = false;
00808 }
00809 }
00810
00811 DBG(cerr << "m: " << m << ", b: " << b << endl);DBG(cerr << "use_missing: " << use_missing << ", missing: " << missing << endl);
00812
00813
00814
00815 BaseType *dest = 0;
00816 double *data;
00817 if (argv[0]->type() == dods_grid_c) {
00818 Array &source = *dynamic_cast<Grid&>(*argv[0]).get_array();
00819 source.set_send_p(true);
00820 source.read(dataset);
00821 data = extract_double_array(&source);
00822 int length = source.length();
00823 int i = 0;
00824 while (i < length) {
00825 DBG2(cerr << "data[" << i << "]: " << data[i] << endl);
00826 if (!use_missing || !double_eq(data[i], missing))
00827 data[i] = data[i] * m + b;
00828 DBG2(cerr << " >> data[" << i << "]: " << data[i] << endl);
00829 ++i;
00830 }
00831
00832
00833 Float64 *temp_f = new Float64(source.name());
00834 source.add_var(temp_f);
00835 source.val2buf(static_cast<void*>(data), false);
00836 delete temp_f;
00837 dest = argv[0];
00838 } else if (argv[0]->is_vector_type()) {
00839 Array &source = dynamic_cast<Array&>(*argv[0]);
00840 source.set_send_p(true);
00841
00842
00843 if (source.get_parent() && source.get_parent()->type() == dods_grid_c)
00844 source.get_parent()->read(dataset);
00845 else
00846 source.read(dataset);
00847
00848 data = extract_double_array(&source);
00849 int length = source.length();
00850 int i = 0;
00851 while (i < length) {
00852 if (!use_missing || !double_eq(data[i], missing))
00853 data[i] = data[i] * m + b;
00854 ++i;
00855 }
00856
00857 Float64 *temp_f = new Float64(source.name());
00858 source.add_var(temp_f);
00859 source.val2buf(static_cast<void*>(data), false);
00860 delete temp_f;
00861
00862 dest = argv[0];
00863 } else if (argv[0]->is_simple_type() && !(argv[0]->type() == dods_str_c
00864 || argv[0]->type() == dods_url_c)) {
00865 double data = extract_double_value(argv[0]);
00866 if (!use_missing || !double_eq(data, missing))
00867 data = data * m + b;
00868
00869 dest = new Float64(argv[0]->name());
00870 dest->val2buf(static_cast<void*>(&data));
00871 } else {
00872 throw Error("The linear_scale() function works only for numeric Grids, Arrays and scalars.");
00873 }
00874
00875 return dest;
00876 }
00877
00894 BaseType * function_geoarray(int argc, BaseType * argv[], DDS &,
00895 const string & dataset)
00896 {
00897 string
00898 info =
00899 "<?xml version=\"1.0\" encoding=\"UTF-8\"?>\
00900 <function name=\"geoarray\" version=\"0.9b1\"/>\
00901 The geoarray() function supports two different sets of arguments:\
00902 geoarray(var,left,top,right,bottom)\
00903 geoarray(var,left,top,right,bottom,var_left,var_top,var_right,var_bottom)\
00904 In the first version 'var' is the target of the selection and 'left', 'top',\
00905 'right' and 'bottom' are the corners of a longitude-latitude box that defines\
00906 the selection. In the second version the 'var_left', ..., parameters give the\
00907 longitude and latitude extent of the entire array. The projection and datum are\
00908 assumed to be Plat-Carre and WGS84.\
00909 </function>";
00910
00911 if (argc == 0) {
00912 Str *response = new Str("version");
00913 response->set_value(info);
00914 return response;
00915 }
00916
00917 DBG(cerr << "argc = " << argc << endl);
00918 if (!(argc == 5 || argc == 9 || argc == 11))
00919 throw Error("Wrong number of arguments to geoarray(). See geoarray() for more information.");
00920
00921
00922 Array *l_array = dynamic_cast < Array * >(argv[0]->ptr_duplicate());
00923 if (!l_array)
00924 throw Error("The first argument to geoarray() must be an Array variable!");
00925
00926 try {
00927
00928
00929 double bb_top = extract_double_value(argv[1]);
00930 double bb_left = extract_double_value(argv[2]);
00931 double bb_bottom = extract_double_value(argv[3]);
00932 double bb_right = extract_double_value(argv[4]);
00933
00934 ArrayGeoConstraint *agc = 0;
00935 switch (argc) {
00936 case 5:
00937 agc = new ArrayGeoConstraint(l_array, dataset);
00938 break;
00939 case 9: {
00940 double var_top = extract_double_value(argv[5]);
00941 double var_left = extract_double_value(argv[6]);
00942 double var_bottom = extract_double_value(argv[7]);
00943 double var_right = extract_double_value(argv[8]);
00944 agc = new ArrayGeoConstraint(l_array, dataset,
00945 var_left, var_top, var_right, var_bottom);
00946 break;
00947 }
00948 case 11: {
00949 double var_top = extract_double_value(argv[5]);
00950 double var_left = extract_double_value(argv[6]);
00951 double var_bottom = extract_double_value(argv[7]);
00952 double var_right = extract_double_value(argv[8]);
00953 string projection = extract_string_argument(argv[9]);
00954 string datum = extract_string_argument(argv[10]);
00955 agc = new ArrayGeoConstraint(l_array, dataset,
00956 var_left, var_top, var_right, var_bottom,
00957 projection, datum);
00958 break;
00959 }
00960 default:
00961 throw InternalErr(__FILE__, __LINE__, "Wrong number of args to geoarray.");
00962 }
00963
00964 agc->set_bounding_box(bb_left, bb_top, bb_right, bb_bottom);
00965
00966
00967 agc->apply_constraint_to_data();
00968
00969
00970
00971 return agc->get_constrained_array();
00972 }
00973 catch (Error & e) {
00974 throw e;
00975 }
00976 catch (exception & e) {
00977 throw
00978 InternalErr(string
00979 ("A C++ exception was thrown from inside geoarray(): ")
00980 + e.what());
00981
00982 }
00983 }
00984
00985 void register_functions(ConstraintEvaluator & ce)
00986 {
00987 ce.add_function("grid", function_grid);
00988 ce.add_function("geogrid", function_geogrid);
00989 ce.add_function("linear_scale", function_linear_scale);
00990 ce.add_function("geoarray", function_geoarray);
00991 ce.add_function("version", function_version);
00992 }
00993
00994 }