39 #include <gdal_priv.h> 40 #include <gdal_utils.h> 41 #include <ogr_spatialref.h> 42 #include <gdalwarper.h> 54 #include <BESDapError.h> 56 #include "ScaleGrid.h" 58 #define DEBUG_KEY "geo" 61 using namespace libdap;
69 inline static int is_valid_lat(
const double lat)
71 return (lat >= -90 && lat <= 90);
74 inline static int is_valid_lon(
const double lon)
76 return (lon >= -180 && lon <= 180);
86 SizeBox get_size_box(Array *x, Array *y)
88 int src_x_size = x->dimension_size(x->dim_begin(),
true);
89 int src_y_size = y->dimension_size(y->dim_begin(),
true);
91 return SizeBox(src_x_size, src_y_size);
101 static bool same_as(
const double a,
const double b)
104 return fabs(a - b) <= numeric_limits<float>::epsilon();
113 bool monotonic_and_uniform(
const vector<double> &values,
double res)
115 vector<double>::size_type end_index = values.size() - 1;
116 for (vector<double>::size_type i = 0; i < end_index; ++i) {
119 if (!same_as((values[i+1] - values[i]), res))
137 vector<double> get_geotransform_data(Array *x, Array *y,
bool test_maps )
143 SizeBox size = get_size_box(x, y);
146 vector<double> y_values(size.y_size);
147 extract_double_array(y, y_values);
149 double res_y = (y_values[y_values.size()-1] - y_values[0]) / (y_values.size() -1);
151 if (test_maps && !monotonic_and_uniform(y_values, res_y)){
152 string msg =
"The grids maps/dimensions must be monotonic and uniform (" + y->name() +
").";
153 BESDEBUG(DEBUG_KEY,
"ERROR get_geotransform_data(): " << msg << endl);
154 throw BESError(msg,BES_SYNTAX_USER_ERROR,__FILE__,__LINE__);
157 vector<double> x_values(size.x_size);
158 extract_double_array(x, x_values);
160 double res_x = (x_values[x_values.size()-1] - x_values[0]) / (x_values.size() -1);
162 if (test_maps && !monotonic_and_uniform(x_values, res_x)){
163 string msg =
"The grids maps/dimensions must be monotonic and uniform (" + x->name() +
").";
164 BESDEBUG(DEBUG_KEY,
"ERROR get_geotransform_data(): " << msg << endl);
165 throw BESError(msg,BES_SYNTAX_USER_ERROR,__FILE__,__LINE__);
181 vector<double> geo_transform(6);
182 geo_transform[0] = x_values[0];
183 geo_transform[1] = res_x;
184 geo_transform[2] = 0;
185 geo_transform[3] = y_values[0];
186 geo_transform[4] = 0;
187 geo_transform[5] = res_y;
189 return geo_transform;
203 vector<GDAL_GCP> get_gcp_data(Array *x, Array *y,
int sample_x,
int sample_y)
205 SizeBox size = get_size_box(x, y);
208 vector<double> y_values(size.y_size);
209 extract_double_array(y, y_values);
212 vector<double> x_values(size.x_size);
213 extract_double_array(x, x_values);
222 unsigned long n_gcps = (size.x_size/sample_x) * (size.y_size/sample_y);
224 vector<GDAL_GCP> gcp_list(n_gcps);
225 GDALInitGCPs(n_gcps, &gcp_list[0]);
227 unsigned long count = 0;
228 for (
int i = 0; i < size.x_size; i += sample_x) {
231 for (
int j = 0; count < n_gcps && j < size.y_size; j += sample_y) {
232 gcp_list[count].dfGCPLine = j;
233 gcp_list[count].dfGCPPixel = i;
234 gcp_list[count].dfGCPX = x_values[i];
235 gcp_list[count].dfGCPY = y_values[j];
244 GDALDataType get_array_type(
const Array *a)
246 switch (const_cast<Array*>(a)->var()->type()) {
276 string msg =
"Cannot perform geo-spatial operations on " 277 + const_cast<Array*>(a)->var()->type_name() +
" data.";
278 BESDEBUG(DEBUG_KEY,
"ERROR get_array_type(): " << msg << endl);
279 throw BESError(msg,BES_SYNTAX_USER_ERROR,__FILE__,__LINE__);
293 template <
typename T>
294 static Array *transfer_values_helper(GDALRasterBand *band,
const unsigned long x,
const unsigned long y, Array *a)
297 vector<T> buf(x * y);
298 CPLErr error = band->RasterIO(GF_Read, 0, 0, x, y, &buf[0], x, y, get_array_type(a), 0, 0);
300 if (error != CPLE_None){
301 string msg = string(
"Could not extract data for array.") + CPLGetLastErrorMsg();
302 BESDEBUG(DEBUG_KEY,
"ERROR transfer_values_helper(): " << msg << endl);
303 throw BESError(msg,BES_SYNTAX_USER_ERROR,__FILE__,__LINE__);
305 a->set_value(buf, buf.size());
317 Array *build_array_from_gdal_dataset(GDALDataset *source,
const Array *dest)
320 GDALRasterBand *band = source->GetRasterBand(1);
321 unsigned long x = band->GetXSize();
322 unsigned long y = band->GetYSize();
325 Array *result =
new Array(
"result", const_cast<Array*>(dest)->var()->ptr_duplicate());
326 result->append_dim(y);
327 result->append_dim(x);
330 switch (result->var()->type()) {
332 return transfer_values_helper<dods_byte>(source->GetRasterBand(1), x, y, result);
335 return transfer_values_helper<dods_uint16>(source->GetRasterBand(1), x, y, result);
338 return transfer_values_helper<dods_int16>(source->GetRasterBand(1), x, y, result);
341 return transfer_values_helper<dods_uint32>(source->GetRasterBand(1), x, y, result);
344 return transfer_values_helper<dods_int32>(source->GetRasterBand(1), x, y, result);
347 return transfer_values_helper<dods_float32>(source->GetRasterBand(1), x, y, result);
350 return transfer_values_helper<dods_float64>(source->GetRasterBand(1), x, y, result);
353 return transfer_values_helper<dods_byte>(source->GetRasterBand(1), x, y, result);
356 return transfer_values_helper<dods_int8>(source->GetRasterBand(1), x, y, result);
361 string msg =
"The source array to a geo-function contained an unsupported numeric type.";
362 BESDEBUG(DEBUG_KEY,
"ERROR build_array_from_gdal_dataset(): " << msg << endl);
363 throw BESError(msg,BES_SYNTAX_USER_ERROR,__FILE__,__LINE__);
389 void build_maps_from_gdal_dataset(GDALDataset *dst, Array *x_map, Array *y_map,
bool name_maps )
392 vector<double> gt(6);
393 dst->GetGeoTransform(>[0]);
396 GDALRasterBand *band = dst->GetRasterBand(1);
399 unsigned long x = band->GetXSize();
402 x_map->append_dim(x,
"Longitude");
405 x_map->append_dim(x);
409 vector<dods_float32> x_map_vals(x);
410 dods_float32 *cur_x = &x_map_vals[0];
411 dods_float32 *prev_x = cur_x;
414 for (
unsigned long i = 1; i < x; ++i) {
417 *cur_x++ = *prev_x++ + gt[1];
420 x_map->set_value(&x_map_vals[0], x);
423 unsigned long y = band->GetYSize();
426 y_map->append_dim(y,
"Latitude");
429 y_map->append_dim(y);
433 vector<dods_float32> y_map_vals(y);
434 dods_float32 *cur_y = &y_map_vals[0];
435 dods_float32 *prev_y = cur_y;
438 for (
unsigned long i = 1; i < y; ++i) {
441 *cur_y++ = *prev_y++ + gt[5];
444 y_map->set_value(&y_map_vals[0], y);
447 void build_maps_from_gdal_dataset_3D(GDALDataset *dst, Array *t_map, Array *x_map, Array *y_map,
bool name_maps )
450 vector<double> gt(6);
451 dst->GetGeoTransform(>[0]);
452 BESDEBUG(DEBUG_KEY,
"build_maps_from_gdal_dataset_3D: Band number: " << dst->GetRasterCount() << endl);
454 for(
int j=1; j<=dst->GetRasterCount(); j++ ){
456 GDALRasterBand *band = dst->GetRasterBand(j);
459 unsigned long x = band->GetXSize();
462 x_map->append_dim(x,
"Longitude");
465 x_map->append_dim(x);
469 vector<dods_float32> x_map_vals(x);
470 dods_float32 *cur_x = &x_map_vals[0];
471 dods_float32 *prev_x = cur_x;
474 for (
unsigned long i = 1; i < x; ++i) {
477 *cur_x++ = *prev_x++ + gt[1];
480 x_map->set_value(&x_map_vals[0], x);
483 unsigned long y = band->GetYSize();
486 y_map->append_dim(y,
"Latitude");
489 y_map->append_dim(y);
493 vector<dods_float32> y_map_vals(y);
494 dods_float32 *cur_y = &y_map_vals[0];
495 dods_float32 *prev_y = cur_y;
498 for (
unsigned long i = 1; i < y; ++i) {
501 *cur_y++ = *prev_y++ + gt[5];
504 y_map->set_value(&y_map_vals[0], y);
516 double get_missing_data_value(
const Array *src)
518 Array *a = const_cast<Array*>(src);
521 string mv_attr = a->get_attr_table().get_attr(
"missing_value");
522 if (mv_attr.empty()) mv_attr = a->get_attr_table().get_attr(
"_FillValue");
524 double missing_data = numeric_limits<double>::quiet_NaN();
525 if (!mv_attr.empty()) {
527 missing_data = strtod(mv_attr.c_str(), &endptr);
528 if (missing_data == 0.0 && endptr == mv_attr.c_str())
529 missing_data = numeric_limits<double>::quiet_NaN();
541 Array::Dim_iter get_x_dim(
const libdap::Array *src)
543 Array *a = const_cast<Array*>(src);
544 int numDims = a->dimensions();
547 ss <<
"Ouch! Retrieving the 'x' dimension for the array ";
548 a->print_decl(ss,
"",
false,
true,
true);
549 ss <<
" FAILED Because it has less than 2 dimensions" << endl;
550 BESDEBUG(DEBUG_KEY, ss.str());
551 throw BESError(ss.str(),BES_SYNTAX_USER_ERROR,__FILE__,__LINE__);
553 Array::Dim_iter start = a->dim_begin();
554 Array::Dim_iter xDim = start + numDims - 1;
565 Array::Dim_iter get_y_dim(
const libdap::Array *src)
567 Array *a = const_cast<Array*>(src);
568 int numDims = a->dimensions();
571 ss <<
"Ouch! Retrieving the 'y' dimension for the array ";
572 a->print_decl(ss,
"",
false,
true,
true);
573 ss <<
" FAILED Because it has less than 2 dimensions" << endl;
574 BESDEBUG(DEBUG_KEY, ss.str());
575 throw BESError(ss.str(),BES_SYNTAX_USER_ERROR,__FILE__,__LINE__);
577 Array::Dim_iter start = a->dim_begin();
578 Array::Dim_iter yDim = start + numDims - 2;
593 bool array_is_effectively_2D(
const libdap::Array *src)
595 Array *a = const_cast<Array*>(src);
596 int numDims = a->dimensions();
597 if (numDims == 2)
return true;
598 if (numDims < 2)
return false;
600 Array::Dim_iter xDim = get_x_dim(a);
601 for (Array::Dim_iter thisDim = a->dim_begin(); thisDim < xDim; thisDim++) {
602 unsigned long size = a->dimension_size(thisDim,
true);
621 void read_band_data(
const Array *src, GDALRasterBand* band)
623 Array *a = const_cast<Array*>(src);
625 if (!array_is_effectively_2D(src)) {
627 ss <<
"Cannot perform geo-spatial operations on an Array (";
628 ss << a->name() <<
") with " << a->dimensions() <<
" dimensions.";
629 ss <<
"Because the constrained shape of the array: ";
630 a->print_decl(ss,
"",
false,
true,
true);
631 ss <<
" is not a two-dimensional array." << endl;
632 BESDEBUG(DEBUG_KEY, ss.str());
633 throw BESError(ss.str(), BES_SYNTAX_USER_ERROR, __FILE__, __LINE__);
639 unsigned long x = a->dimension_size(get_x_dim(a),
true);
640 unsigned long y = a->dimension_size(get_y_dim(a),
true);
647 CPLErr error = band->RasterIO(GF_Write, 0, 0, x, y, a->get_buf(), x, y, get_array_type(a), 0, 0);
649 if (error != CPLE_None){
650 string msg =
"Could not load data for grid '" + a->name() +
"' msg: '" + CPLGetLastErrorMsg() +
"'";
651 BESDEBUG(DEBUG_KEY,
"ERROR read_band_data(): " << msg << endl);
652 throw BESError(msg,BES_SYNTAX_USER_ERROR,__FILE__,__LINE__);
666 void add_band_data(
const Array *src, GDALDataset* ds)
668 Array *a = const_cast<Array*>(src);
675 char **options = NULL;
677 oss << reinterpret_cast<unsigned long>(a->get_buf());
678 options = CSLSetNameValue(options,
"DATAPOINTER", oss.str().c_str());
680 CPLErr error = ds->AddBand(get_array_type(a), options);
684 if (error != CPLE_None){
685 string msg =
"Could not add data for grid '" + a->name() +
"': " + CPLGetLastErrorMsg();
686 BESDEBUG(DEBUG_KEY,
"ERROR add_band_data(): " << msg << endl);
687 throw BESError(msg,BES_INTERNAL_ERROR,__FILE__,__LINE__);
710 auto_ptr<GDALDataset> build_src_dataset(Array *data, Array *x, Array *y,
const string &srs)
712 GDALDriver *driver = GetGDALDriverManager()->GetDriverByName(
"MEM");
714 string msg = string(
"Could not get the Memory driver for GDAL: ") + CPLGetLastErrorMsg();
715 BESDEBUG(DEBUG_KEY,
"ERROR build_src_dataset(): " << msg << endl);
716 throw BESError(msg,BES_INTERNAL_ERROR,__FILE__,__LINE__);
719 SizeBox array_size = get_size_box(x, y);
722 auto_ptr<GDALDataset> ds(driver->Create(
"result", array_size.x_size, array_size.y_size,
723 1 , get_array_type(data), NULL ));
726 add_band_data(data, ds.get());
730 GDALRasterBand *band = ds->GetRasterBand(1);
732 string msg =
"Could not get the GDAL RasterBand for Array '" + data->name() +
"': " + CPLGetLastErrorMsg();
733 BESDEBUG(DEBUG_KEY,
"ERROR build_src_dataset(): " << msg << endl);
734 throw BESError(msg,BES_INTERNAL_ERROR,__FILE__,__LINE__);
738 double no_data = get_missing_data_value(data);
739 band->SetNoDataValue(no_data);
742 read_band_data(data, band);
745 vector<double> geo_transform = get_geotransform_data(x, y);
746 ds->SetGeoTransform(&geo_transform[0]);
748 OGRSpatialReference native_srs;
749 if (CE_None != native_srs.SetWellKnownGeogCS(srs.c_str())){
750 string msg =
"Could not set '" + srs +
"' as the dataset native CRS.";
751 BESDEBUG(DEBUG_KEY,
"ERROR build_src_dataset(): " << msg << endl);
752 throw BESError(msg,BES_SYNTAX_USER_ERROR,__FILE__,__LINE__);
758 char *pszSRS_WKT = NULL;
759 native_srs.exportToWkt( &pszSRS_WKT );
760 ds->SetProjection( pszSRS_WKT );
761 CPLFree( pszSRS_WKT );
776 auto_ptr<GDALDataset> scale_dataset(auto_ptr<GDALDataset> src,
const SizeBox &size,
const string &crs ,
777 const string &interp )
780 argv = CSLAddString(argv,
"-of");
781 argv = CSLAddString(argv,
"MEM");
783 argv = CSLAddString(argv,
"-outsize");
786 argv = CSLAddString(argv, oss.str().c_str());
789 argv = CSLAddString(argv, oss.str().c_str());
791 argv = CSLAddString(argv,
"-b");
792 argv = CSLAddString(argv,
"1");
794 argv = CSLAddString(argv,
"-r");
795 argv = CSLAddString(argv, interp.c_str());
798 argv = CSLAddString(argv,
"-a_srs");
799 argv = CSLAddString(argv, crs.c_str());
802 if (BESISDEBUG(DEBUG_KEY)) {
805 BESDEBUG(DEBUG_KEY,
"argv: " << *local++ << endl);
812 cerr <<
"argv: " << *local++ << endl;
816 GDALTranslateOptions *options = GDALTranslateOptionsNew(argv, NULL );
818 int usage_error = CE_None;
819 GDALDatasetH dst_handle = GDALTranslate(
"warped_dst", src.get(), options, &usage_error);
820 if (!dst_handle || usage_error != CE_None) {
821 GDALClose(dst_handle);
822 GDALTranslateOptionsFree(options);
823 string msg = string(
"Error calling GDAL translate: ") + CPLGetLastErrorMsg();
824 BESDEBUG(DEBUG_KEY,
"ERROR scale_dataset(): " << msg << endl);
825 throw BESError(msg, BES_INTERNAL_ERROR, __FILE__, __LINE__);
828 auto_ptr<GDALDataset> dst(static_cast<GDALDataset*>(dst_handle));
830 GDALTranslateOptionsFree(options);
836 auto_ptr<GDALDataset> scale_dataset_3D(auto_ptr<GDALDataset> src,
const SizeBox &size,
const string &crs ,
837 const string &interp )
840 argv = CSLAddString(argv,
"-of");
841 argv = CSLAddString(argv,
"MEM");
843 argv = CSLAddString(argv,
"-outsize");
846 argv = CSLAddString(argv, oss.str().c_str());
849 argv = CSLAddString(argv, oss.str().c_str());
852 int n_bands = src.get()->GetRasterCount();
853 for(
int i=0; i < n_bands; i++){
856 argv = CSLAddString(argv,
"-b");
857 argv = CSLAddString(argv, oss.str().c_str());
860 argv = CSLAddString(argv,
"-r");
861 argv = CSLAddString(argv, interp.c_str());
864 argv = CSLAddString(argv,
"-a_srs");
865 argv = CSLAddString(argv, crs.c_str());
868 if (BESISDEBUG(DEBUG_KEY)) {
871 BESDEBUG(DEBUG_KEY,
"argv: " << *local++ << endl);
878 cerr <<
"argv: " << *local++ << endl;
882 GDALTranslateOptions *options = GDALTranslateOptionsNew(argv, NULL );
883 int usage_error = CE_None;
884 GDALDatasetH dst_handle = GDALTranslate(
"warped_dst", src.get(), options, &usage_error);
885 if (!dst_handle || usage_error != CE_None) {
886 GDALClose(dst_handle);
887 GDALTranslateOptionsFree(options);
888 string msg = string(
"Error calling GDAL translate: ") + CPLGetLastErrorMsg();
889 BESDEBUG(DEBUG_KEY,
"ERROR scale_dataset(): " << msg << endl);
890 throw BESError(msg, BES_INTERNAL_ERROR, __FILE__, __LINE__);
893 auto_ptr<GDALDataset> dst(static_cast<GDALDataset*>(dst_handle));
895 GDALTranslateOptionsFree(options);
913 Grid *scale_dap_array(
const Array *data,
const Array *x,
const Array *y,
const SizeBox &size,
914 const string &crs,
const string &interp)
917 Array *d = const_cast<Array*>(data);
919 auto_ptr<GDALDataset> src = build_src_dataset(d, const_cast<Array*>(x), const_cast<Array*>(y));
922 auto_ptr<GDALDataset> dst = scale_dataset(src, size, crs, interp);
925 auto_ptr<Array> built_data(build_array_from_gdal_dataset(dst.get(), d));
927 auto_ptr<Array> built_lat(
new Array(y->name(),
new Float32(y->name())));
928 auto_ptr<Array> built_lon(
new Array(x->name(),
new Float32(x->name())));
930 build_maps_from_gdal_dataset(dst.get(), built_lon.get(), built_lat.get());
932 auto_ptr<Grid> result(
new Grid(d->name()));
933 result->set_array(built_data.release());
934 result->add_map(built_lat.release(),
false);
935 result->add_map(built_lon.release(),
false);
937 return result.release();
950 Grid *scale_dap_grid(
const Grid *g,
const SizeBox &size,
const string &crs,
const string &interp)
952 string func(__func__);
956 throw BESError(func+
"The Grid object is null.",BES_INTERNAL_ERROR,__FILE__,__LINE__);
959 Array *data = dynamic_cast<Array*>(const_cast<Grid*>(g)->array_var());
961 throw BESError(func+
"Unable to obtain data array from Grid '"+g->name()+
"'",BES_INTERNAL_ERROR,__FILE__,__LINE__);
964 Grid::Map_riter ritr = const_cast<Grid*>(g)->map_rbegin();
965 Array *x = dynamic_cast<Array*>(*ritr);
966 Array *y = dynamic_cast<Array*>(*++ritr);
969 throw BESError(func+
"Unable to obtain 2 Map arrays from Grid '"+g->name()+
"'",BES_INTERNAL_ERROR,__FILE__,__LINE__);
972 return scale_dap_array(data, x, y, size, crs, interp);
994 auto_ptr<GDALDataset> build_src_dataset_3D(Array *data, Array *t, Array *x, Array *y,
const string &srs)
996 Array *d = dynamic_cast<Array*>(data);
998 GDALDriver *driver = GetGDALDriverManager()->GetDriverByName(
"MEM");
1000 string msg = string(
"Could not get the Memory driver for GDAL: ") + CPLGetLastErrorMsg();
1001 BESDEBUG(DEBUG_KEY,
"ERROR build_src_dataset(): " << msg << endl);
1002 throw BESError(msg,BES_INTERNAL_ERROR,__FILE__,__LINE__);
1005 SizeBox array_size = get_size_box(x, y);
1006 int nBands = t->length();
1007 BESDEBUG(DEBUG_KEY,
"nBands = " << nBands << endl);
1008 int nBytes = data->prototype()->width();
1009 const int data_size = x->length() * y->length();
1010 unsigned int dsize = data_size * nBytes;
1012 auto_ptr<GDALDataset> ds(driver->Create(
"result", array_size.x_size, array_size.y_size, nBands, get_array_type(d),
1016 for(
int i=1; i<=nBands; i++){
1018 GDALRasterBand *band = ds->GetRasterBand(i);
1020 string msg =
"Could not get the GDAL RasterBand for Array '" + data->name() +
"': " + CPLGetLastErrorMsg();
1021 BESDEBUG(DEBUG_KEY,
"ERROR build_src_dataset(): " << msg << endl);
1022 throw BESError(msg,BES_INTERNAL_ERROR,__FILE__,__LINE__);
1025 double no_data = get_missing_data_value(data);
1026 band->SetNoDataValue(no_data);
1029 CPLErr error = band->RasterIO(GF_Write, 0, 0, x->length(), y->length(), data->get_buf() + dsize*(i-1), x->length(), y->length(), get_array_type(data), 0, 0);
1030 if (error != CPLE_None)
1031 throw Error(
"Could not write data for band: " + long_to_string(i) +
": " +
string(CPLGetLastErrorMsg()));
1036 vector<double> geo_transform = get_geotransform_data(x, y);
1037 ds->SetGeoTransform(&geo_transform[0]);
1039 OGRSpatialReference native_srs;
1040 if (CE_None != native_srs.SetWellKnownGeogCS(srs.c_str())){
1041 string msg =
"Could not set '" + srs +
"' as the dataset native CRS.";
1042 BESDEBUG(DEBUG_KEY,
"ERROR build_src_dataset(): " << msg << endl);
1043 throw BESError(msg,BES_SYNTAX_USER_ERROR,__FILE__,__LINE__);
1047 char *pszSRS_WKT = NULL;
1048 native_srs.exportToWkt( &pszSRS_WKT );
1049 ds->SetProjection( pszSRS_WKT );
1050 CPLFree( pszSRS_WKT );
1068 Grid *scale_dap_array_3D(
const Array *data,
const Array *t,
const Array *x,
const Array *y,
const SizeBox &size,
1069 const string &crs,
const string &interp)
1072 Array *d = const_cast<Array*>(data);
1074 auto_ptr<GDALDataset> src = build_src_dataset_3D(d, const_cast<Array*>(t), const_cast<Array*>(x), const_cast<Array*>(y));
1078 auto_ptr<GDALDataset> dst = scale_dataset_3D(src, size, crs, interp);
1081 auto_ptr<Array> built_data(build_array_from_gdal_dataset(dst.get(), d));
1082 auto_ptr<Array> built_time(
new Array(t->name(),
new Float32(t->name())));
1083 auto_ptr<Array> built_lat(
new Array(y->name(),
new Float32(y->name())));
1084 auto_ptr<Array> built_lon(
new Array(x->name(),
new Float32(x->name())));
1086 build_maps_from_gdal_dataset_3D(dst.get(), built_time.get(), built_lon.get(), built_lat.get());
1088 auto_ptr<Grid> result(
new Grid(d->name()));
1089 result->set_array(built_data.release());
1090 result->add_map(built_time.release(),
false);
1091 result->add_map(built_lat.release(),
false);
1092 result->add_map(built_lon.release(),
false);
1093 BESDEBUG(DEBUG_KEY,
"result length: " << result.release()->get_array()->dimensions() << endl);
1094 return result.release();
Abstract exception class for the BES with basic string message.