# Contents of /branches/arrayview_from_1695_trunk/escript/src/DataFactory.cpp

Revision 1721 - (show annotations)
Fri Aug 22 00:39:32 2008 UTC (14 years, 7 months ago) by jfenwick
File size: 11683 byte(s)
Branch commit.

Fixed problems with copyFromNumArray
Removed version of setTaggedValueFromCPP() which required DataArrayView.
Updated tests.

 1 2 /* \$Id\$ */ 3 4 /******************************************************* 5 * 6 * Copyright 2003-2007 by ACceSS MNRF 7 * Copyright 2007 by University of Queensland 8 * 9 10 * Primary Business: Queensland, Australia 11 * Licensed under the Open Software License version 3.0 12 13 * 14 *******************************************************/ 15 16 #include "DataFactory.h" 17 #include "esysUtils/esys_malloc.h" 18 19 #include 20 #include 21 #include 22 #ifdef USE_NETCDF 23 #include 24 #endif 25 26 using namespace boost::python; 27 28 namespace escript { 29 30 Data 31 Scalar(double value, 32 const FunctionSpace& what, 33 bool expanded) 34 { 35 // 36 // an empty shape is a scalar 37 DataTypes::ShapeType shape; 38 return Data(value,shape,what,expanded); 39 } 40 41 Data 42 Vector(double value, 43 const FunctionSpace& what, 44 bool expanded) 45 { 46 DataTypes::ShapeType shape(1,what.getDomain().getDim()); 47 return Data(value,shape,what,expanded); 48 } 49 50 Data 51 Tensor(double value, 52 const FunctionSpace& what, 53 bool expanded) 54 { 55 DataTypes::ShapeType shape(2,what.getDomain().getDim()); 56 return Data(value,shape,what,expanded); 57 } 58 59 Data 60 Tensor3(double value, 61 const FunctionSpace& what, 62 bool expanded) 63 { 64 DataTypes::ShapeType shape(3,what.getDomain().getDim()); 65 return Data(value,shape,what,expanded); 66 } 67 68 Data 69 Tensor4(double value, 70 const FunctionSpace& what, 71 bool expanded) 72 { 73 DataTypes::ShapeType shape(4,what.getDomain().getDim()); 74 return Data(value,shape,what,expanded); 75 } 76 77 Data 78 load(const std::string fileName, 79 const AbstractDomain& domain) 80 { 81 #ifdef PASO_MPI 82 throw DataException("Error - load :: not implemented for MPI yet."); 83 #endif 84 85 #ifdef USE_NETCDF 86 NcAtt *type_att, *rank_att, *function_space_type_att; 87 // netCDF error handler 88 NcError err(NcError::silent_nonfatal); 89 // Create the file. 90 NcFile dataFile(fileName.c_str(), NcFile::ReadOnly); 91 if (!dataFile.is_valid()) 92 throw DataException("Error - load:: opening of netCDF file for input failed."); 93 /* recover function space */ 94 if (! (function_space_type_att=dataFile.get_att("function_space_type")) ) 95 throw DataException("Error - load:: cannot recover function_space_type attribute from escript netCDF file."); 96 int function_space_type = function_space_type_att->as_int(0); 97 delete function_space_type_att; 98 /* test if function space id is valid and create function space instance */ 99 if (! domain.isValidFunctionSpaceType(function_space_type) ) 100 throw DataException("Error - load:: function space type code in netCDF file is invalid for given domain."); 101 FunctionSpace function_space=FunctionSpace(domain, function_space_type); 102 /* recover rank */ 103 if (! (rank_att=dataFile.get_att("rank")) ) 104 throw DataException("Error - load:: cannot recover rank attribute from escript netCDF file."); 105 int rank = rank_att->as_int(0); 106 delete rank_att; 107 if (rank<0 || rank>DataTypes::maxRank) 108 throw DataException("Error - load:: rank in escript netCDF file is greater than maximum rank."); 109 /* recover type attribute */ 110 int type=-1; 111 if ((type_att=dataFile.get_att("type")) ) { 112 char* type_str = type_att->as_string(0); 113 if (strncmp(type_str, "constant", strlen("constant")) == 0 ) { 114 type =0; 115 } else if (strncmp(type_str, "tagged", strlen("tagged")) == 0 ) { 116 type =1; 117 } else if (strncmp(type_str, "expanded", strlen("expanded")) == 0 ) { 118 type =2; 119 } 120 esysUtils::free(type_str); 121 } else { 122 if (! (type_att=dataFile.get_att("type_id")) ) 123 throw DataException("Error - load:: cannot recover type attribute from escript netCDF file."); 124 type=type_att->as_int(0); 125 } 126 delete type_att; 127 128 /* recover dimension */ 129 int ndims=dataFile.num_dims(); 130 int ntags =0 , num_samples =0 , num_data_points_per_sample =0, d=0, len_data_point=1; 131 NcDim *d_dim, *tags_dim, *num_samples_dim, *num_data_points_per_sample_dim; 132 /* recover shape */ 133 DataTypes::ShapeType shape; 134 long dims[DataTypes::maxRank+2]; 135 if (rank>0) { 136 if (! (d_dim=dataFile.get_dim("d0")) ) 137 throw DataException("Error - load:: unable to recover d0 from netCDF file."); 138 d=d_dim->size(); 139 shape.push_back(d); 140 dims[0]=d; 141 len_data_point*=d; 142 } 143 if (rank>1) { 144 if (! (d_dim=dataFile.get_dim("d1")) ) 145 throw DataException("Error - load:: unable to recover d1 from netCDF file."); 146 d=d_dim->size(); 147 shape.push_back(d); 148 dims[1]=d; 149 len_data_point*=d; 150 } 151 if (rank>2) { 152 if (! (d_dim=dataFile.get_dim("d2")) ) 153 throw DataException("Error - load:: unable to recover d2 from netCDF file."); 154 d=d_dim->size(); 155 shape.push_back(d); 156 dims[2]=d; 157 len_data_point*=d; 158 } 159 if (rank>3) { 160 if (! (d_dim=dataFile.get_dim("d3")) ) 161 throw DataException("Error - load:: unable to recover d3 from netCDF file."); 162 d=d_dim->size(); 163 shape.push_back(d); 164 dims[3]=d; 165 len_data_point*=d; 166 } 167 /* recover stuff */ 168 Data out; 169 NcVar *var, *ids_var, *tags_var; 170 if (type == 0) { 171 /* constant data */ 172 if ( ! ( (ndims == rank && rank >0) || ( ndims ==1 && rank == 0 ) ) ) 173 throw DataException("Error - load:: illegal number of dimensions for constant data in netCDF file."); 174 if (rank == 0) { 175 if (! (d_dim=dataFile.get_dim("l")) ) 176 throw DataException("Error - load:: unable to recover d0 for scalar constant data in netCDF file."); 177 int d0 = d_dim->size(); 178 if (! d0 == 1) 179 throw DataException("Error - load:: d0 is expected to be one for scalar constant data in netCDF file."); 180 dims[0]=1; 181 } 182 out=Data(0,shape,function_space); 183 if (!(var = dataFile.get_var("data"))) 184 throw DataException("Error - load:: unable to find data in netCDF file."); 185 if (! var->get(&(out.getDataPoint(0,0).getData()[0]), dims) ) 186 throw DataException("Error - load:: unable to recover data from netCDF file."); 187 } else if (type == 1) { 188 /* tagged data */ 189 if ( ! (ndims == rank + 1) ) 190 throw DataException("Error - load:: illegal number of dimensions for tagged data in netCDF file."); 191 if (! (tags_dim=dataFile.get_dim("num_tags")) ) 192 throw DataException("Error - load:: unable to recover number of tags from netCDF file."); 193 ntags=tags_dim->size(); 194 dims[rank]=ntags; 195 int *tags = (int *) esysUtils::malloc(ntags*sizeof(int)); 196 if (! ( tags_var = dataFile.get_var("tags")) ) 197 { 198 esysUtils::free(tags); 199 throw DataException("Error - load:: unable to find tags in netCDF file."); 200 } 201 if (! tags_var->get(tags, ntags) ) 202 { 203 esysUtils::free(tags); 204 throw DataException("Error - load:: unable to recover tags from netCDF file."); 205 } 206 207 DataVector data(len_data_point * ntags, 0., len_data_point * ntags); 208 if (!(var = dataFile.get_var("data"))) 209 { 210 esysUtils::free(tags); 211 throw DataException("Error - load:: unable to find data in netCDF file."); 212 } 213 if (! var->get(&(data[0]), dims) ) 214 { 215 esysUtils::free(tags); 216 throw DataException("Error - load:: unable to recover data from netCDF file."); 217 } 218 out=Data(DataArrayView(data,shape,0),function_space); 219 for (int t=1; tsize(); 231 if ( ! (num_data_points_per_sample_dim = dataFile.get_dim("num_data_points_per_sample") ) ) 232 throw DataException("Error - load:: unable to recover number of data points per sample from netCDF file."); 233 num_data_points_per_sample=num_data_points_per_sample_dim->size(); 234 // check shape: 235 if ( ! (num_samples == function_space.getNumSamples() && num_data_points_per_sample == function_space.getNumDataPointsPerSample()) ) 236 throw DataException("Error - load:: data sample layout of file does not match data layout of function space."); 237 // get ids 238 if (! ( ids_var = dataFile.get_var("id")) ) 239 throw DataException("Error - load:: unable to find reference ids in netCDF file."); 240 const int* ids_p=function_space.borrowSampleReferenceIDs(); 241 int *ids_of_nc = (int *)esysUtils::malloc(num_samples*sizeof(int)); 242 if (! ids_var->get(ids_of_nc, (long) num_samples) ) 243 { 244 esysUtils::free(ids_of_nc); 245 throw DataException("Error - load:: unable to recover ids from netCDF file."); 246 } 247 // check order: 248 int failed=-1, local_failed=-1, i; 249 #pragma omp parallel private(local_failed) 250 { 251 local_failed=-1; 252 #pragma omp for private(i) schedule(static) 253 for (i=0;i < num_samples; ++i) { 254 if (ids_of_nc[i]!=ids_p[i]) local_failed=i; 255 } 256 #pragma omp critical 257 if (local_failed>=0) failed = local_failed; 258 } 259 /* if (failed>=0) 260 { 261 esysUtils::free(ids_of_nc); 262 throw DataException("Error - load:: data ordering in netCDF file does not match ordering of FunctionSpace."); 263 } */ 264 // get the data: 265 dims[rank]=num_data_points_per_sample; 266 dims[rank+1]=num_samples; 267 out=Data(0,shape,function_space,true); 268 if (!(var = dataFile.get_var("data"))) 269 { 270 esysUtils::free(ids_of_nc); 271 throw DataException("Error - load:: unable to find data in netCDF file."); 272 } 273 if (! var->get(&(out.getDataPoint(0,0).getData()[0]), dims) ) 274 { 275 esysUtils::free(ids_of_nc); 276 throw DataException("Error - load:: unable to recover data from netCDF file."); 277 } 278 if (failed>=0) { 279 try { 280 std::cout << "Information - load: start reordering data from netCDF file " << fileName << std::endl; 281 out.borrowData()->reorderByReferenceIDs(ids_of_nc); 282 } 283 catch (std::exception&) { 284 esysUtils::free(ids_of_nc); 285 throw DataException("Error - load:: unable to reorder data in netCDF file."); 286 } 287 } 288 } else { 289 throw DataException("Error - load:: unknown escript data type in netCDF file."); 290 } 291 return out; 292 #else 293 throw DataException("Error - load:: is not compiled with netCFD. Please contact your insatllation manager."); 294 #endif 295 } 296 297 bool 298 loadConfigured() 299 { 300 #ifdef USE_NETCDF 301 return true; 302 #else 303 return false; 304 #endif 305 } 306 307 Data 308 convertToData(const boost::python::object& value, 309 const FunctionSpace& what) 310 { 311 // first we try to extract a Data object from value 312 extract value_data(value); 313 if (value_data.check()) { 314 Data extracted_data=value_data(); 315 if (extracted_data.isEmpty()) { 316 return extracted_data; 317 } else { 318 return Data(extracted_data,what); 319 } 320 } else { 321 return Data(value,what); 322 } 323 } 324 325 } // end of namespace

## Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision

 ViewVC Help Powered by ViewVC 1.1.26