/[escript]/trunk/escript/src/DataFactory.cpp
ViewVC logotype

Contents of /trunk/escript/src/DataFactory.cpp

Parent Directory Parent Directory | Revision Log Revision Log


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

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26