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

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

Parent Directory Parent Directory | Revision Log Revision Log


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

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26