/[escript]/trunk/ripley/src/Rectangle.cpp
ViewVC logotype

Diff of /trunk/ripley/src/Rectangle.cpp

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 4615 by caltinay, Mon Jan 13 05:05:33 2014 UTC revision 4618 by caltinay, Wed Jan 15 04:35:19 2014 UTC
# Line 169  bool Rectangle::operator==(const Abstrac Line 169  bool Rectangle::operator==(const Abstrac
169  }  }
170    
171  void Rectangle::readNcGrid(escript::Data& out, string filename, string varname,  void Rectangle::readNcGrid(escript::Data& out, string filename, string varname,
172              const GridParameters& params) const              const ReaderParameters& params) const
173  {  {
174  #ifdef USE_NETCDF  #ifdef USE_NETCDF
175      // check destination function space      // check destination function space
# Line 195  void Rectangle::readNcGrid(escript::Data Line 195  void Rectangle::readNcGrid(escript::Data
195      for (size_t i=0; i<params.multiplier.size(); i++)      for (size_t i=0; i<params.multiplier.size(); i++)
196          if (params.multiplier[i]<1)          if (params.multiplier[i]<1)
197              throw RipleyException("readNcGrid(): all multipliers must be positive");              throw RipleyException("readNcGrid(): all multipliers must be positive");
198        if (params.reverse.size() != 2)
199            throw RipleyException("readNcGrid(): argument 'reverse' must have 2 entries");
200    
201      // check file existence and size      // check file existence and size
202      NcFile f(filename.c_str(), NcFile::ReadOnly);      NcFile f(filename.c_str(), NcFile::ReadOnly);
# Line 232  void Rectangle::readNcGrid(escript::Data Line 234  void Rectangle::readNcGrid(escript::Data
234      // first coordinates in data object to write to      // first coordinates in data object to write to
235      const int first0 = max(0, params.first[0]-m_offset[0]);      const int first0 = max(0, params.first[0]-m_offset[0]);
236      const int first1 = max(0, params.first[1]-m_offset[1]);      const int first1 = max(0, params.first[1]-m_offset[1]);
237      // indices to first value in file      // indices to first value in file (not accounting for reverse yet)
238      const int idx0 = max(0, m_offset[0]-params.first[0]);      int idx0 = max(0, m_offset[0]-params.first[0]);
239      const int idx1 = max(0, m_offset[1]-params.first[1]);      int idx1 = max(0, m_offset[1]-params.first[1]);
240      // number of values to read      // number of values to read
241      const int num0 = min(params.numValues[0]-idx0, myN0-first0);      const int num0 = min(params.numValues[0]-idx0, myN0-first0);
242      const int num1 = min(params.numValues[1]-idx1, myN1-first1);      const int num1 = min(params.numValues[1]-idx1, myN1-first1);
243    
244        // make sure we read the right block if going backwards through file
245        if (params.reverse[0])
246            idx0 = edges[dims-1]-num0-idx0;
247        if (dims>1 && params.reverse[1])
248            idx1 = edges[dims-2]-num1-idx1;
249    
250      vector<double> values(num0*num1);      vector<double> values(num0*num1);
251      if (dims==2) {      if (dims==2) {
252          var->set_cur(idx1, idx0);          var->set_cur(idx1, idx0);
# Line 251  void Rectangle::readNcGrid(escript::Data Line 259  void Rectangle::readNcGrid(escript::Data
259      const int dpp = out.getNumDataPointsPerSample();      const int dpp = out.getNumDataPointsPerSample();
260      out.requireWrite();      out.requireWrite();
261    
262        // helpers for reversing
263        const int x0 = (params.reverse[0] ? num0-1 : 0);
264        const int x_mult = (params.reverse[0] ? -1 : 1);
265        const int y0 = (params.reverse[1] ? num1-1 : 0);
266        const int y_mult = (params.reverse[1] ? -1 : 1);
267    
268      for (index_t y=0; y<num1; y++) {      for (index_t y=0; y<num1; y++) {
269  #pragma omp parallel for  #pragma omp parallel for
270          for (index_t x=0; x<num0; x++) {          for (index_t x=0; x<num0; x++) {
271              const int baseIndex = first0+x*params.multiplier[0]              const int baseIndex = first0+x*params.multiplier[0]
272                                    +(first1+y*params.multiplier[1])*myN0;                                    +(first1+y*params.multiplier[1])*myN0;
273              const int srcIndex = y*num0+x;              const int srcIndex = (y0+y_mult*y)*num0+(x0+x_mult*x);
274              if (!isnan(values[srcIndex])) {              if (!isnan(values[srcIndex])) {
275                  for (index_t m1=0; m1<params.multiplier[1]; m1++) {                  for (index_t m1=0; m1<params.multiplier[1]; m1++) {
276                      for (index_t m0=0; m0<params.multiplier[0]; m0++) {                      for (index_t m0=0; m0<params.multiplier[0]; m0++) {
# Line 276  void Rectangle::readNcGrid(escript::Data Line 290  void Rectangle::readNcGrid(escript::Data
290  }  }
291    
292  void Rectangle::readBinaryGrid(escript::Data& out, string filename,  void Rectangle::readBinaryGrid(escript::Data& out, string filename,
293                                 const GridParameters& params) const                                 const ReaderParameters& params) const
294  {  {
295      // the mapping is not universally correct but should work on our      // the mapping is not universally correct but should work on our
296      // supported platforms      // supported platforms
# Line 297  void Rectangle::readBinaryGrid(escript:: Line 311  void Rectangle::readBinaryGrid(escript::
311    
312  template<typename ValueType>  template<typename ValueType>
313  void Rectangle::readBinaryGridImpl(escript::Data& out, const string& filename,  void Rectangle::readBinaryGridImpl(escript::Data& out, const string& filename,
314                                     const GridParameters& params) const                                     const ReaderParameters& params) const
315  {  {
316      // check destination function space      // check destination function space
317      int myN0, myN1;      int myN0, myN1;

Legend:
Removed from v.4615  
changed lines
  Added in v.4618

  ViewVC Help
Powered by ViewVC 1.1.26