/[escript]/branches/arrexp_2137_win/escript/src/DataTagged.cpp
ViewVC logotype

Annotation of /branches/arrexp_2137_win/escript/src/DataTagged.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1802 - (hide annotations)
Tue Sep 23 01:03:29 2008 UTC (10 years, 8 months ago) by jfenwick
Original Path: trunk/escript/src/DataTagged.cpp
File size: 34685 byte(s)
Added canTag methods to FunctionSpace and AbstractDomain (and its 
offspring).
This checks to see if the domain supports tags for the given type of 
function space.

Constructors for DataTagged now throw exceptions if you attempt to make 
a DataTagged with a FunctionSpace which does not support tags.
To allow the default constructor to work, NullDomain has a single 
functioncode which "supports" tagging.

Fixed a bug in DataTagged::toString and DataTypes::pointToString.

Added FunctionSpace::getListOfTagsSTL.

algorithm(DataTagged, BinaryFunction) in DataAlgorithm now only 
processes tags known to be in use.
This fixes mantis issue #0000186.

Added comment to Data.h intro warning about holding references if the 
underlying DataAbstract changes.

_python_ unit tests have been updated to test TaggedData with invalid 
FunctionSpaces and to give the correct answers to Lsup etc.


1 jgs 121
2 ksteube 1312 /* $Id$ */
3 jgs 82
4 ksteube 1312 /*******************************************************
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 jgs 474 #include "DataTagged.h"
17 phornby 1628 #include "esysUtils/esys_malloc.h"
18 jgs 496
19 jgs 474 #include "DataConstant.h"
20     #include "DataException.h"
21 gross 1023 #ifdef USE_NETCDF
22 ksteube 1312 #include <netcdfcpp.h>
23 gross 1023 #endif
24 jfenwick 1796 #include "DataMaths.h"
25 jgs 82
26     using namespace std;
27    
28     namespace escript {
29    
30 jgs 121 DataTagged::DataTagged()
31 jfenwick 1796 : DataAbstract(FunctionSpace(),DataTypes::scalarShape)
32 jgs 82 {
33 jgs 509 // default constructor
34    
35 jgs 82 // create a scalar default value
36 jgs 151 m_data.resize(1,0.,1);
37 jfenwick 1796 /* DataArrayView temp(m_data,DataTypes::ShapeType());
38     setPointDataView(temp);*/
39 jgs 82 }
40    
41 jfenwick 1796 // DataTagged::DataTagged(const TagListType& tagKeys,
42     // const ValueListType& values,
43     // const DataArrayView& defaultValue,
44     // const FunctionSpace& what)
45     // : DataAbstract(what)
46     // {
47     // // constructor
48     //
49     // // initialise the array of data values
50     // // the default value is always the first item in the values list
51     // int len = defaultValue.noValues();
52     // m_data.resize(len,0.,len);
53     // for (int i=0; i<defaultValue.noValues(); i++) {
54     // m_data[i]=defaultValue.getData(i);
55     // }
56     //
57     // // create the data view
58     // DataArrayView temp(m_data,defaultValue.getShape());
59     // setPointDataView(temp);
60     //
61     // // add remaining tags and values
62     // addTaggedValues(tagKeys,values);
63     // }
64 jgs 509
65 jgs 119 DataTagged::DataTagged(const FunctionSpace& what,
66 jfenwick 1796 const DataTypes::ShapeType &shape,
67 jgs 119 const int tags[],
68 jgs 500 const ValueType& data)
69 jfenwick 1796 : DataAbstract(what,shape)
70 jgs 119 {
71 jgs 509 // alternative constructor
72     // not unit_tested tested yet
73 jfenwick 1796 // It is not explicitly unit tested yet, but it is called from DataFactory
74 jgs 509
75 jfenwick 1802 if (!what.canTag())
76     {
77     throw DataException("Programming error - DataTag created with a non-taggable FunctionSpace.");
78     }
79 jgs 121 // copy the data
80 jgs 119 m_data=data;
81 jgs 121
82 jgs 119 // create the view of the data
83 jfenwick 1796 // DataArrayView tempView(m_data,shape);
84     // setPointDataView(tempView);
85 jgs 121
86 jfenwick 1796 // we can't rely on the tag array to give us the number of tags so
87     // use the data we have been passed
88     int valsize=DataTypes::noValues(shape);
89     int ntags=data.size()/valsize;
90    
91 jgs 119 // create the tag lookup map
92 jfenwick 1796 // we assume that the first value and first tag are the default value so we skip
93     for (int i=1;i<ntags;++i)
94     {
95     m_offsetLookup.insert(DataMapType::value_type(tags[i],i*valsize));
96 jgs 119 }
97     }
98    
99 woo409 757 DataTagged::DataTagged(const FunctionSpace& what,
100 jfenwick 1796 const DataTypes::ShapeType &shape,
101 woo409 757 const TagListType& tags,
102     const ValueType& data)
103 jfenwick 1796 : DataAbstract(what,shape)
104 woo409 757 {
105     // alternative constructor
106    
107 jfenwick 1802 if (!what.canTag())
108     {
109     throw DataException("Programming error - DataTag created with a non-taggable FunctionSpace.");
110     }
111    
112 woo409 757 // copy the data
113     m_data=data;
114    
115     // create the view of the data
116 jfenwick 1796 // DataArrayView tempView(m_data,shape);
117     // setPointDataView(tempView);
118 woo409 757
119     // create the tag lookup map
120 jfenwick 1796
121     // for (int sampleNo=0; sampleNo<getNumSamples(); sampleNo++) {
122     // m_offsetLookup.insert(DataMapType::value_type(sampleNo,tags[sampleNo]));
123     // }
124    
125     // The above code looks like it will create a map the wrong way around
126    
127     int valsize=DataTypes::noValues(shape);
128     int npoints=(data.size()/valsize)-1;
129     int ntags=tags.size();
130     if (ntags>npoints)
131     { // This throw is not unit tested yet
132     throw DataException("Programming error - Too many tags for the supplied values.");
133 woo409 757 }
134 jfenwick 1796
135     // create the tag lookup map
136     // we assume that the first value is the default value so we skip it (hence the i+1 below)
137     for (int i=0;i<ntags;++i)
138     {
139     m_offsetLookup.insert(DataMapType::value_type(tags[i],(i+1)*valsize));
140     }
141 woo409 757 }
142    
143    
144 jgs 102 DataTagged::DataTagged(const DataTagged& other)
145 jfenwick 1796 : DataAbstract(other.getFunctionSpace(),other.getShape()),
146 jgs 82 m_data(other.m_data),
147     m_offsetLookup(other.m_offsetLookup)
148     {
149 jgs 509 // copy constructor
150    
151 jgs 119 // create the data view
152 jfenwick 1796 // DataArrayView temp(m_data,other.getPointDataView().getShape());
153     // setPointDataView(temp);
154 jgs 82 }
155    
156 jgs 102 DataTagged::DataTagged(const DataConstant& other)
157 jfenwick 1796 : DataAbstract(other.getFunctionSpace(),other.getShape())
158 jgs 82 {
159 jgs 509 // copy constructor
160    
161 jfenwick 1802 if (!other.getFunctionSpace().canTag())
162     {
163     throw DataException("Programming error - DataTag created with a non-taggable FunctionSpace.");
164     }
165    
166 jgs 121 // fill the default value with the constant value item from "other"
167 jfenwick 1796 // const DataArrayView& value=other.getPointDataView();
168     int len = other.getNoValues();
169 jgs 151 m_data.resize(len,0.,len);
170 jfenwick 1796 for (int i=0; i<len; i++) {
171     m_data[i]=other.getVector()[i];
172 jgs 121 }
173    
174 jfenwick 1796 // // create the data view
175     // DataArrayView temp(m_data,value.getShape());
176     // setPointDataView(temp);
177 jgs 82 }
178    
179 jfenwick 1796
180     // Create a new object by copying tags
181     DataTagged::DataTagged(const FunctionSpace& what,
182     const DataTypes::ShapeType& shape,
183     const DataTypes::ValueType& defaultvalue,
184     const DataTagged* tagsource)
185     : DataAbstract(what,shape)
186     {
187     // This constructor has not been unit tested yet
188    
189     if (defaultvalue.size()!=DataTypes::noValues(shape)) {
190     throw DataException("Programming error - defaultvalue does not match supplied shape.");
191     }
192    
193    
194 jfenwick 1802 if (!what.canTag())
195     {
196     throw DataException("Programming error - DataTag created with a non-taggable FunctionSpace.");
197     }
198    
199 jfenwick 1796 if (tagsource!=0)
200     {
201     int numtags=tagsource->getTagLookup().size();
202     // m_offsetLookup.reserve(tagsource.getTagLookup().size());
203     m_data.resize(defaultvalue.size(),0.); // since this is tagged data, we should have blocksize=1
204    
205     DataTagged::DataMapType::const_iterator i;
206     for (i=tagsource->getTagLookup().begin();i!=tagsource->getTagLookup().end();i++) {
207     addTag(i->first);
208     }
209     }
210     else
211     {
212     m_data.resize(defaultvalue.size());
213     }
214    
215    
216    
217     // need to set the default value ....
218     for (int i=0; i<defaultvalue.size(); i++) {
219     m_data[i]=defaultvalue[i];
220     }
221     }
222    
223 jfenwick 1799 DataAbstract*
224     DataTagged::deepCopy()
225     {
226     return new DataTagged(*this);
227     }
228 jfenwick 1796
229 jgs 82 DataAbstract*
230 jfenwick 1796 DataTagged::getSlice(const DataTypes::RegionType& region) const
231 jgs 82 {
232 jgs 513 return new DataTagged(*this, region);
233 jgs 82 }
234    
235 jgs 513 DataTagged::DataTagged(const DataTagged& other,
236 jfenwick 1796 const DataTypes::RegionType& region)
237     : DataAbstract(other.getFunctionSpace(),DataTypes::getResultSliceShape(region))
238 jgs 513 {
239     // slice constructor
240    
241     // get the shape of the slice to copy from other
242 jfenwick 1796 DataTypes::ShapeType regionShape(DataTypes::getResultSliceShape(region));
243     DataTypes::RegionLoopRangeType regionLoopRange=DataTypes::getSliceRegionLoopRange(region);
244 jgs 513
245     // allocate enough space in this for all values
246     // (need to add one to allow for the default value)
247 jfenwick 1796 int len = DataTypes::noValues(regionShape)*(other.m_offsetLookup.size()+1);
248 jgs 513 m_data.resize(len,0.0,len);
249    
250     // create the data view
251 jfenwick 1796 // DataArrayView temp(m_data,regionShape);
252     // setPointDataView(temp);
253 jgs 513
254     // copy the default value from other to this
255 jfenwick 1796 // getDefaultValue().copySlice(other.getDefaultValue(), regionLoopRange);
256     const DataTypes::ShapeType& otherShape=other.getShape();
257     const DataTypes::ValueType& otherData=other.getVector();
258     DataTypes::copySlice(getVector(),getShape(),getDefaultOffset(),otherData,otherShape,other.getDefaultOffset(), regionLoopRange);
259 jgs 513
260     // loop through the tag values copying these
261     DataMapType::const_iterator pos;
262 jfenwick 1796 DataTypes::ValueType::size_type tagOffset=getNoValues();
263 jgs 513 for (pos=other.m_offsetLookup.begin();pos!=other.m_offsetLookup.end();pos++){
264 jfenwick 1796 // getPointDataView().copySlice(tagOffset,other.getPointDataView(),pos->second,regionLoopRange);
265     DataTypes::copySlice(m_data,getShape(),tagOffset,otherData, otherShape, pos->second, regionLoopRange);
266 jgs 513 m_offsetLookup.insert(DataMapType::value_type(pos->first,tagOffset));
267 jfenwick 1796 tagOffset+=getNoValues();
268 jgs 513 }
269     }
270    
271 jgs 82 void
272 jgs 513 DataTagged::setSlice(const DataAbstract* other,
273 jfenwick 1796 const DataTypes::RegionType& region)
274 jgs 82 {
275 jgs 513
276     // other must be another DataTagged object
277     // Data:setSlice implementation should ensure this
278     const DataTagged* otherTemp=dynamic_cast<const DataTagged*>(other);
279     if (otherTemp==0) {
280 jgs 82 throw DataException("Programming error - casting to DataTagged.");
281     }
282 jgs 121
283 jgs 513 // determine shape of the specified region
284 jfenwick 1796 DataTypes::ShapeType regionShape(DataTypes::getResultSliceShape(region));
285 jgs 513
286     // modify region specification as needed to match rank of this object
287 jfenwick 1796 DataTypes::RegionLoopRangeType regionLoopRange=DataTypes::getSliceRegionLoopRange(region);
288 jgs 513
289     // ensure rank/shape of this object is compatible with specified region
290 jfenwick 1796 if (getRank()!=region.size()) {
291 jgs 82 throw DataException("Error - Invalid slice region.");
292     }
293 jfenwick 1796 if (otherTemp->getRank()>0 && !DataTypes::checkShape(other->getShape(),regionShape)) {
294     throw DataException (DataTypes::createShapeErrorMessage(
295     "Error - Couldn't copy slice due to shape mismatch.",regionShape,other->getShape()));
296 jgs 82 }
297 jgs 121
298 jfenwick 1796 const DataTypes::ValueType& otherData=otherTemp->getVector();
299     const DataTypes::ShapeType& otherShape=otherTemp->getShape();
300 jgs 513 // copy slice from other default value to this default value
301 jfenwick 1796 // getDefaultValue().copySliceFrom(otherTemp->getDefaultValue(), regionLoopRange);
302     DataTypes::copySliceFrom(m_data,getShape(),getDefaultOffset(),otherData,otherShape,otherTemp->getDefaultOffset(),regionLoopRange);
303 jgs 121
304 jgs 541 // loop through tag values in other, adding any which aren't in this, using default value
305     DataMapType::const_iterator pos;
306     for (pos=otherTemp->m_offsetLookup.begin();pos!=otherTemp->m_offsetLookup.end();pos++) {
307     if (!isCurrentTag(pos->first)) {
308 jfenwick 1796 addTag(pos->first);
309 jgs 541 }
310     }
311    
312 jgs 513 // loop through the tag values copying slices from other to this
313 jgs 121 for (pos=m_offsetLookup.begin();pos!=m_offsetLookup.end();pos++) {
314 jfenwick 1796 // getDataPointByTag(pos->first).copySliceFrom(otherTemp->getDataPointByTag(pos->first), regionLoopRange);
315     DataTypes::copySliceFrom(m_data,getShape(),getOffsetForTag(pos->first),otherData, otherShape, otherTemp->getOffsetForTag(pos->first), regionLoopRange);
316    
317 jgs 121 }
318 jgs 513
319 jgs 82 }
320    
321 jgs 149 int
322     DataTagged::getTagNumber(int dpno)
323     {
324     //
325     // Get the number of samples and data-points per sample
326     int numSamples = getNumSamples();
327     int numDataPointsPerSample = getNumDPPSample();
328     int numDataPoints = numSamples * numDataPointsPerSample;
329    
330     if (numDataPointsPerSample==0) {
331     throw DataException("DataTagged::getTagNumber error: no data-points associated with this object.");
332     }
333    
334 jgs 496 if (dpno<0 || dpno>numDataPoints-1) {
335 jgs 149 throw DataException("DataTagged::getTagNumber error: invalid data-point number supplied.");
336     }
337    
338     //
339     // Determine the sample number which corresponds to this data-point number
340     int sampleNo = dpno / numDataPointsPerSample;
341    
342     //
343     // Determine the tag number which corresponds to this sample number
344     int tagNo = getFunctionSpace().getTagFromSampleNo(sampleNo);
345    
346     //
347     // return the tag number
348     return(tagNo);
349     }
350    
351 jfenwick 1796 // void
352     // DataTagged::setTaggedValues(const TagListType& tagKeys,
353     // const ValueListType& values)
354     // {
355     // addTaggedValues(tagKeys,values);
356     // }
357    
358 jgs 82 void
359 jfenwick 1796 DataTagged::setTaggedValue(int tagKey,
360     const DataTypes::ShapeType& pointshape,
361     const ValueType& value,
362     int dataOffset)
363 jgs 509 {
364 jfenwick 1796 if (!DataTypes::checkShape(getShape(), pointshape)) {
365     throw DataException(DataTypes::createShapeErrorMessage(
366     "Error - Cannot setTaggedValue due to shape mismatch.", pointshape,getShape()));
367     }
368     DataMapType::iterator pos(m_offsetLookup.find(tagKey));
369     if (pos==m_offsetLookup.end()) {
370     // tag couldn't be found so use addTaggedValue
371     addTaggedValue(tagKey,pointshape, value, dataOffset);
372     } else {
373     // copy the values into the data array at the offset determined by m_offsetLookup
374     int offset=pos->second;
375     for (int i=0; i<getNoValues(); i++) {
376     m_data[offset+i]=value[i+dataOffset];
377     }
378     }
379 jgs 509 }
380    
381 jfenwick 1796
382     /*
383 jgs 509 void
384 jgs 82 DataTagged::setTaggedValue(int tagKey,
385     const DataArrayView& value)
386     {
387 jgs 121 if (!getPointDataView().checkShape(value.getShape())) {
388 jfenwick 1796 throw DataException(DataTypes::createShapeErrorMessage(
389     "Error - Cannot setTaggedValue due to shape mismatch.", value.getShape(),getShape()));
390 jgs 121 }
391 jgs 82 DataMapType::iterator pos(m_offsetLookup.find(tagKey));
392     if (pos==m_offsetLookup.end()) {
393 jgs 121 // tag couldn't be found so use addTaggedValue
394 jgs 82 addTaggedValue(tagKey,value);
395     } else {
396 jgs 121 // copy the values into the data array at the offset determined by m_offsetLookup
397     int offset=pos->second;
398     for (int i=0; i<getPointDataView().noValues(); i++) {
399     m_data[offset+i]=value.getData(i);
400 jgs 82 }
401     }
402 jfenwick 1796 }*/
403    
404     // void
405     // DataTagged::addTaggedValues(const TagListType& tagKeys,
406     // const ValueListType& values)
407     // {
408     // if (values.size()==0) {
409     // // copy the current default value for each of the tags
410     // TagListType::const_iterator iT;
411     // for (iT=tagKeys.begin();iT!=tagKeys.end();iT++) {
412     // // the point data view for DataTagged points at the default value
413     // addTaggedValue(*iT,getPointDataView());
414     // }
415     // } else if (values.size()==1 && tagKeys.size()>1) {
416     // // assume the one given value will be used for all tag values
417     // TagListType::const_iterator iT;
418     // for (iT=tagKeys.begin();iT!=tagKeys.end();iT++) {
419     // addTaggedValue(*iT,values[0]);
420     // }
421     // } else {
422     // if (tagKeys.size()!=values.size()) {
423     // stringstream temp;
424     // temp << "Error - (addTaggedValue) Number of tags: " << tagKeys.size()
425     // << " doesn't match number of values: " << values.size();
426     // throw DataException(temp.str());
427     // } else {
428     // unsigned int i;
429     // for (i=0;i<tagKeys.size();i++) {
430     // addTaggedValue(tagKeys[i],values[i]);
431     // }
432     // }
433     // }
434     // }
435    
436    
437     void
438     DataTagged::addTaggedValues(const TagListType& tagKeys,
439     const ValueBatchType& values,
440     const ShapeType& vShape)
441     {
442     DataTypes::ValueType t(values.size(),0);
443     for (int i=0;i<values.size();++i)
444     {
445     t[i]=values[i];
446     }
447     addTaggedValues(tagKeys,t,vShape);
448 jgs 82 }
449    
450 jfenwick 1796
451     // Note: The check to see if vShape==our shape is done in the addTaggedValue method
452 jgs 82 void
453 jgs 509 DataTagged::addTaggedValues(const TagListType& tagKeys,
454 jfenwick 1796 const ValueType& values,
455     const ShapeType& vShape)
456 jgs 509 {
457 jfenwick 1796 int n=getNoValues();
458     int numVals=values.size()/n;
459 jgs 509 if (values.size()==0) {
460     // copy the current default value for each of the tags
461     TagListType::const_iterator iT;
462     for (iT=tagKeys.begin();iT!=tagKeys.end();iT++) {
463     // the point data view for DataTagged points at the default value
464 jfenwick 1796 addTag(*iT);
465 jgs 509 }
466 jfenwick 1796 } else if (numVals==1 && tagKeys.size()>1) {
467 jgs 509 // assume the one given value will be used for all tag values
468     TagListType::const_iterator iT;
469     for (iT=tagKeys.begin();iT!=tagKeys.end();iT++) {
470 jfenwick 1796 addTaggedValue(*iT, vShape, values,0);
471 jgs 509 }
472     } else {
473 jfenwick 1796 if (tagKeys.size()!=numVals) {
474 jgs 509 stringstream temp;
475     temp << "Error - (addTaggedValue) Number of tags: " << tagKeys.size()
476     << " doesn't match number of values: " << values.size();
477     throw DataException(temp.str());
478     } else {
479 phornby 1455 unsigned int i;
480 jfenwick 1796 int offset=0;
481     for (i=0;i<tagKeys.size();i++ ,offset+=n) {
482     addTaggedValue(tagKeys[i],vShape,values,offset);
483 jgs 509 }
484     }
485     }
486     }
487    
488 jfenwick 1796
489    
490    
491 jgs 509 void
492 jgs 82 DataTagged::addTaggedValue(int tagKey,
493 jfenwick 1796 const DataTypes::ShapeType& pointshape,
494     const ValueType& value,
495     int dataOffset)
496 jgs 82 {
497 jfenwick 1796 if (!DataTypes::checkShape(getShape(), pointshape)) {
498     throw DataException(DataTypes::createShapeErrorMessage(
499     "Error - Cannot addTaggedValue due to shape mismatch.", pointshape,getShape()));
500 jgs 82 }
501 jgs 121 DataMapType::iterator pos(m_offsetLookup.find(tagKey));
502     if (pos!=m_offsetLookup.end()) {
503     // tag already exists so use setTaggedValue
504 jfenwick 1796 setTaggedValue(tagKey,pointshape, value, dataOffset);
505 jgs 121 } else {
506     // save the key and the location of its data in the lookup tab
507     m_offsetLookup.insert(DataMapType::value_type(tagKey,m_data.size()));
508     // add the data given in "value" at the end of m_data
509     // need to make a temp copy of m_data, resize m_data, then copy
510     // all the old values plus the value to be added back into m_data
511     ValueType m_data_temp(m_data);
512     int oldSize=m_data.size();
513 jfenwick 1796 int newSize=m_data.size()+getNoValues();
514 jgs 151 m_data.resize(newSize,0.,newSize);
515 jgs 121 for (int i=0;i<oldSize;i++) {
516     m_data[i]=m_data_temp[i];
517     }
518 jfenwick 1796 for (int i=0;i<getNoValues();i++) {
519     m_data[oldSize+i]=value[i+dataOffset];
520 jgs 121 }
521     }
522 jgs 82 }
523    
524 jfenwick 1796
525    
526    
527     // void
528     // DataTagged::addTaggedValue(int tagKey,
529     // const DataArrayView& value)
530     // {
531     // if (!getPointDataView().checkShape(value.getShape())) {
532     // throw DataException(DataTypes::createShapeErrorMessage(
533     // "Error - Cannot addTaggedValue due to shape mismatch.", value.getShape(),getShape()));
534     // }
535     // DataMapType::iterator pos(m_offsetLookup.find(tagKey));
536     // if (pos!=m_offsetLookup.end()) {
537     // // tag already exists so use setTaggedValue
538     // setTaggedValue(tagKey,value);
539     // } else {
540     // // save the key and the location of its data in the lookup tab
541     // m_offsetLookup.insert(DataMapType::value_type(tagKey,m_data.size()));
542     // // add the data given in "value" at the end of m_data
543     // // need to make a temp copy of m_data, resize m_data, then copy
544     // // all the old values plus the value to be added back into m_data
545     // ValueType m_data_temp(m_data);
546     // int oldSize=m_data.size();
547     // int newSize=m_data.size()+value.noValues();
548     // m_data.resize(newSize,0.,newSize);
549     // for (int i=0;i<oldSize;i++) {
550     // m_data[i]=m_data_temp[i];
551     // }
552     // for (int i=0;i<value.noValues();i++) {
553     // m_data[oldSize+i]=value.getData(i);
554     // }
555     // }
556     // }
557    
558    
559     void
560     DataTagged::addTag(int tagKey)
561     {
562     DataMapType::iterator pos(m_offsetLookup.find(tagKey));
563     if (pos!=m_offsetLookup.end()) {
564     // tag already exists so use setTaggedValue
565     // setTaggedValue(tagKey,value);
566     } else {
567     // save the key and the location of its data in the lookup tab
568     m_offsetLookup.insert(DataMapType::value_type(tagKey,m_data.size()));
569     // add the data given in "value" at the end of m_data
570     // need to make a temp copy of m_data, resize m_data, then copy
571     // all the old values plus the value to be added back into m_data
572     ValueType m_data_temp(m_data);
573     int oldSize=m_data.size();
574     int newSize=m_data.size()+getNoValues();
575     m_data.resize(newSize,0.,newSize);
576     for (int i=0;i<oldSize;i++) {
577     m_data[i]=m_data_temp[i];
578     }
579     for (int i=0;i<getNoValues();i++) {
580     m_data[oldSize+i]=m_data[m_defaultValueOffset+i];
581     }
582     }
583     }
584    
585    
586 jgs 82 double*
587     DataTagged::getSampleDataByTag(int tag)
588     {
589     DataMapType::iterator pos(m_offsetLookup.find(tag));
590     if (pos==m_offsetLookup.end()) {
591     // tag couldn't be found so return the default value
592     return &(m_data[0]);
593     } else {
594     // return the data-point corresponding to the given tag
595     return &(m_data[pos->second]);
596     }
597     }
598    
599     string
600     DataTagged::toString() const
601     {
602 jfenwick 1796 using namespace escript::DataTypes;
603     string empty="";
604 jgs 82 stringstream temp;
605     DataMapType::const_iterator i;
606     temp << "Tag(Default)" << endl;
607 jfenwick 1796 temp << pointToString(m_data,getShape(),getDefaultOffset(),empty) << endl;
608 jgs 82 // create a temporary view as the offset will be changed
609 jfenwick 1796 // DataArrayView tempView(getPointDataView().getData(), getPointDataView().getShape());
610 jgs 82 for (i=m_offsetLookup.begin();i!=m_offsetLookup.end();++i) {
611     temp << "Tag(" << i->first << ")" << endl;
612 jfenwick 1802 temp << pointToString(m_data,getShape(),i->second,empty) << endl;
613 jfenwick 1796 // tempView.setOffset(i->second);
614     // temp << tempView.toString() << endl;
615 jgs 82 }
616 ksteube 1312 return temp.str();
617 jgs 82 }
618    
619 jfenwick 1796 DataTypes::ValueType::size_type
620 jgs 496 DataTagged::getPointOffset(int sampleNo,
621     int dataPointNo) const
622 jgs 82 {
623 jgs 496 int tagKey=getFunctionSpace().getTagFromSampleNo(sampleNo);
624     DataMapType::const_iterator pos(m_offsetLookup.find(tagKey));
625 jfenwick 1796 DataTypes::ValueType::size_type offset=m_defaultValueOffset;
626 jgs 82 if (pos!=m_offsetLookup.end()) {
627     offset=pos->second;
628     }
629 jgs 496 return offset;
630 jgs 82 }
631    
632 jfenwick 1796 // DataArrayView
633     // DataTagged::getDataPointByTag(int tag) const
634     // {
635     // DataMapType::const_iterator pos(m_offsetLookup.find(tag));
636     // DataTypes::ValueType::size_type offset=m_defaultValueOffset;
637     // if (pos!=m_offsetLookup.end()) {
638     // offset=pos->second;
639     // }
640     // DataArrayView temp(getPointDataView());
641     // temp.setOffset(offset);
642     // return temp;
643     // }
644     //
645    
646    
647     DataTypes::ValueType::size_type
648     DataTagged::getOffsetForTag(int tag) const
649 jgs 82 {
650 jgs 496 DataMapType::const_iterator pos(m_offsetLookup.find(tag));
651 jfenwick 1796 DataTypes::ValueType::size_type offset=m_defaultValueOffset;
652 jgs 82 if (pos!=m_offsetLookup.end()) {
653     offset=pos->second;
654     }
655 jfenwick 1796 return offset;
656 jgs 82 }
657    
658 jfenwick 1796 DataTypes::ValueType::const_reference
659     DataTagged::getDataByTag(int tag, DataTypes::ValueType::size_type i) const
660 jgs 82 {
661 jfenwick 1796 DataMapType::const_iterator pos(m_offsetLookup.find(tag));
662     DataTypes::ValueType::size_type offset=m_defaultValueOffset;
663     if (pos!=m_offsetLookup.end()) {
664     offset=pos->second;
665     }
666     return m_data[offset+i];
667     /* DataArrayView temp(getPointDataView());
668     temp.setOffset(offset);
669     return temp.getData()[offset+i];*/
670 jgs 82 }
671    
672 jfenwick 1796
673     DataTypes::ValueType::reference
674     DataTagged::getDataByTag(int tag, DataTypes::ValueType::size_type i)
675 jgs 123 {
676 jfenwick 1796 DataMapType::const_iterator pos(m_offsetLookup.find(tag));
677     DataTypes::ValueType::size_type offset=m_defaultValueOffset;
678     if (pos!=m_offsetLookup.end()) {
679     offset=pos->second;
680     }
681     return m_data[offset+i];
682     /* DataArrayView temp(getPointDataView());
683     temp.setOffset(offset);
684     return temp.getData()[offset+i];*/
685 jgs 123 }
686    
687 jfenwick 1796
688    
689    
690    
691    
692     // DataArrayView
693     // DataTagged::getDataPoint(int sampleNo,
694     // int dataPointNo)
695     // {
696     // EsysAssert(validSampleNo(sampleNo),"(getDataPoint) Invalid sampleNo: " << sampleNo);
697     // int tagKey=getFunctionSpace().getTagFromSampleNo(sampleNo);
698     // return getDataPointByTag(tagKey);
699     // }
700    
701    
702 gross 580 void
703 ksteube 775 DataTagged::symmetric(DataAbstract* ev)
704     {
705     DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
706     if (temp_ev==0) {
707     throw DataException("Error - DataTagged::symmetric casting to DataTagged failed (probably a programming error).");
708     }
709     const DataTagged::DataMapType& thisLookup=getTagLookup();
710     DataTagged::DataMapType::const_iterator i;
711     DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
712 jfenwick 1796 ValueType& evVec=temp_ev->getVector();
713     const ShapeType& evShape=temp_ev->getShape();
714 ksteube 775 for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
715 jfenwick 1796 temp_ev->addTag(i->first);
716     DataTypes::ValueType::size_type offset=getOffsetForTag(i->first);
717     // DataArrayView thisView=getDataPointByTag(i->first);
718     // DataArrayView evView=temp_ev->getDataPointByTag(i->first);
719     DataTypes::ValueType::size_type evoffset=temp_ev->getOffsetForTag(i->first);
720    
721     // DataArrayView::symmetric(thisView,0,evView,0);
722     DataMaths::symmetric(m_data,getShape(),offset,evVec, evShape, evoffset);
723 ksteube 775 }
724 jfenwick 1796 // symmetric(m_data,getShape(),getDefaultOffset(),
725     DataMaths::symmetric(m_data,getShape(),getDefaultOffset(),evVec,evShape,temp_ev->getDefaultOffset());
726 ksteube 775 }
727 jfenwick 1796
728    
729 ksteube 775 void
730     DataTagged::nonsymmetric(DataAbstract* ev)
731     {
732     DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
733     if (temp_ev==0) {
734     throw DataException("Error - DataTagged::nonsymmetric casting to DataTagged failed (probably a programming error).");
735     }
736     const DataTagged::DataMapType& thisLookup=getTagLookup();
737     DataTagged::DataMapType::const_iterator i;
738     DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
739 jfenwick 1796 ValueType& evVec=temp_ev->getVector();
740     const ShapeType& evShape=temp_ev->getShape();
741 ksteube 775 for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
742 jfenwick 1796 temp_ev->addTag(i->first);
743     /* DataArrayView thisView=getDataPointByTag(i->first);
744     DataArrayView evView=temp_ev->getDataPointByTag(i->first);*/
745     DataTypes::ValueType::size_type offset=getOffsetForTag(i->first);
746     DataTypes::ValueType::size_type evoffset=temp_ev->getOffsetForTag(i->first);
747     DataMaths::nonsymmetric(m_data,getShape(),offset,evVec, evShape, evoffset);
748 ksteube 775 }
749 jfenwick 1796 DataMaths::nonsymmetric(m_data,getShape(),getDefaultOffset(),evVec,evShape,temp_ev->getDefaultOffset());
750 ksteube 775 }
751 jfenwick 1796
752    
753 ksteube 775 void
754 gross 800 DataTagged::trace(DataAbstract* ev, int axis_offset)
755 ksteube 775 {
756     DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
757     if (temp_ev==0) {
758 gross 800 throw DataException("Error - DataTagged::trace casting to DataTagged failed (probably a programming error).");
759 ksteube 775 }
760     const DataTagged::DataMapType& thisLookup=getTagLookup();
761     DataTagged::DataMapType::const_iterator i;
762     DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
763 jfenwick 1796 ValueType& evVec=temp_ev->getVector();
764     const ShapeType& evShape=temp_ev->getShape();
765 ksteube 775 for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
766 jfenwick 1796 temp_ev->addTag(i->first);
767     // DataArrayView thisView=getDataPointByTag(i->first);
768     // DataArrayView evView=temp_ev->getDataPointByTag(i->first);
769     DataTypes::ValueType::size_type offset=getOffsetForTag(i->first);
770     DataTypes::ValueType::size_type evoffset=temp_ev->getOffsetForTag(i->first);
771     DataMaths::trace(m_data,getShape(),offset,evVec, evShape, evoffset, axis_offset);
772 ksteube 775 }
773 jfenwick 1796 DataMaths::trace(m_data,getShape(),getDefaultOffset(),evVec,evShape,temp_ev->getDefaultOffset(),axis_offset);
774 ksteube 775 }
775 gross 800
776 ksteube 775 void
777     DataTagged::transpose(DataAbstract* ev, int axis_offset)
778     {
779     DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
780     if (temp_ev==0) {
781     throw DataException("Error - DataTagged::transpose casting to DataTagged failed (probably a programming error).");
782     }
783     const DataTagged::DataMapType& thisLookup=getTagLookup();
784     DataTagged::DataMapType::const_iterator i;
785     DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
786 jfenwick 1796 ValueType& evVec=temp_ev->getVector();
787     const ShapeType& evShape=temp_ev->getShape();
788 ksteube 775 for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
789 jfenwick 1796 temp_ev->addTag(i->first);
790     // DataArrayView thisView=getDataPointByTag(i->first);
791     // DataArrayView evView=temp_ev->getDataPointByTag(i->first);
792     DataTypes::ValueType::size_type offset=getOffsetForTag(i->first);
793     DataTypes::ValueType::size_type evoffset=temp_ev->getOffsetForTag(i->first);
794     DataMaths::transpose(m_data,getShape(),offset,evVec, evShape, evoffset, axis_offset);
795 ksteube 775 }
796 jfenwick 1796 DataMaths::transpose(m_data,getShape(),getDefaultOffset(),evVec,evShape,temp_ev->getDefaultOffset(),axis_offset);
797 ksteube 775 }
798 gross 800
799 ksteube 775 void
800 gross 804 DataTagged::swapaxes(DataAbstract* ev, int axis0, int axis1)
801 gross 800 {
802     DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
803     if (temp_ev==0) {
804 gross 804 throw DataException("Error - DataTagged::swapaxes casting to DataTagged failed (probably a programming error).");
805 gross 800 }
806     const DataTagged::DataMapType& thisLookup=getTagLookup();
807     DataTagged::DataMapType::const_iterator i;
808     DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
809 jfenwick 1796 ValueType& evVec=temp_ev->getVector();
810     const ShapeType& evShape=temp_ev->getShape();
811 gross 800 for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
812 jfenwick 1796 temp_ev->addTag(i->first);
813     /* DataArrayView thisView=getDataPointByTag(i->first);
814     DataArrayView evView=temp_ev->getDataPointByTag(i->first);*/
815     DataTypes::ValueType::size_type offset=getOffsetForTag(i->first);
816     DataTypes::ValueType::size_type evoffset=temp_ev->getOffsetForTag(i->first);
817     DataMaths::swapaxes(m_data,getShape(),offset,evVec, evShape, evoffset,axis0,axis1);
818 gross 800 }
819 jfenwick 1796 DataMaths::swapaxes(m_data,getShape(),getDefaultOffset(),evVec,evShape,temp_ev->getDefaultOffset(),axis0,axis1);
820 gross 800 }
821    
822     void
823 gross 580 DataTagged::eigenvalues(DataAbstract* ev)
824     {
825     DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
826     if (temp_ev==0) {
827     throw DataException("Error - DataTagged::eigenvalues casting to DataTagged failed (propably a programming error).");
828     }
829     const DataTagged::DataMapType& thisLookup=getTagLookup();
830     DataTagged::DataMapType::const_iterator i;
831     DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
832 jfenwick 1796 ValueType& evVec=temp_ev->getVector();
833     const ShapeType& evShape=temp_ev->getShape();
834 gross 580 for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
835 jfenwick 1796 temp_ev->addTag(i->first);
836     // DataArrayView thisView=getDataPointByTag(i->first);
837     // DataArrayView evView=temp_ev->getDataPointByTag(i->first);
838     DataTypes::ValueType::size_type offset=getOffsetForTag(i->first);
839     DataTypes::ValueType::size_type evoffset=temp_ev->getOffsetForTag(i->first);
840     DataMaths::eigenvalues(m_data,getShape(),offset,evVec, evShape, evoffset);
841 gross 580 }
842 jfenwick 1796 DataMaths::eigenvalues(m_data,getShape(),getDefaultOffset(),evVec, evShape, temp_ev->getDefaultOffset());
843 gross 580 }
844     void
845     DataTagged::eigenvalues_and_eigenvectors(DataAbstract* ev,DataAbstract* V,const double tol)
846     {
847     DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
848     if (temp_ev==0) {
849     throw DataException("Error - DataTagged::eigenvalues_and_eigenvectors casting to DataTagged failed (propably a programming error).");
850     }
851     DataTagged* temp_V=dynamic_cast<DataTagged*>(V);
852     if (temp_V==0) {
853     throw DataException("Error - DataTagged::eigenvalues_and_eigenvectors casting to DataTagged failed (propably a programming error).");
854     }
855     const DataTagged::DataMapType& thisLookup=getTagLookup();
856     DataTagged::DataMapType::const_iterator i;
857     DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
858 jfenwick 1796 ValueType& evVec=temp_ev->getVector();
859     const ShapeType& evShape=temp_ev->getShape();
860     ValueType& VVec=temp_V->getVector();
861     const ShapeType& VShape=temp_V->getShape();
862 gross 580 for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
863 jfenwick 1796 temp_ev->addTag(i->first);
864     temp_V->addTag(i->first);
865     /* DataArrayView thisView=getDataPointByTag(i->first);
866 gross 580 DataArrayView evView=temp_ev->getDataPointByTag(i->first);
867 jfenwick 1796 DataArrayView VView=temp_V->getDataPointByTag(i->first);*/
868     DataTypes::ValueType::size_type offset=getOffsetForTag(i->first);
869     DataTypes::ValueType::size_type evoffset=temp_ev->getOffsetForTag(i->first);
870     DataTypes::ValueType::size_type Voffset=temp_V->getOffsetForTag(i->first);
871     /* DataArrayView::eigenvalues_and_eigenvectors(thisView,0,evView,0,VView,0,tol);*/
872     DataMaths::eigenvalues_and_eigenvectors(m_data,getShape(),offset,evVec, evShape, evoffset,VVec,VShape,Voffset,tol);
873    
874 gross 580 }
875 jfenwick 1796 DataMaths::eigenvalues_and_eigenvectors(m_data,getShape(),getDefaultOffset(),evVec, evShape,
876     temp_ev->getDefaultOffset(),VVec,VShape,
877     temp_V->getDefaultOffset(), tol);
878 jgs 123
879 gross 580
880     }
881    
882 gross 950 void
883 gross 1118 DataTagged::setToZero(){
884 jfenwick 1796 DataTypes::ValueType::size_type n=m_data.size();
885 gross 1118 for (int i=0; i<n ;++i) m_data[i]=0.;
886     }
887    
888     void
889 gross 950 DataTagged::dump(const std::string fileName) const
890     {
891 gross 1023 #ifdef PASO_MPI
892 ksteube 1312 throw DataException("Error - DataTagged:: dump is not implemented for MPI yet.");
893 gross 1023 #endif
894     #ifdef USE_NETCDF
895 jfenwick 1796 const int ldims=DataTypes::maxRank+1;
896 gross 983 const NcDim* ncdims[ldims];
897     NcVar *var, *tags_var;
898 jfenwick 1796 int rank = getRank();
899 gross 983 int type= getFunctionSpace().getTypeCode();
900     int ndims =0;
901     long dims[ldims];
902 gross 1141 const double* d_ptr=&(m_data[0]);
903 jfenwick 1796 DataTypes::ShapeType shape = getShape();
904 gross 983
905     // netCDF error handler
906     NcError err(NcError::verbose_nonfatal);
907     // Create the file.
908     NcFile dataFile(fileName.c_str(), NcFile::Replace);
909     // check if writing was successful
910     if (!dataFile.is_valid())
911     throw DataException("Error - DataTagged:: opening of netCDF file for output failed.");
912 gross 1141 if (!dataFile.add_att("type_id",1) )
913 gross 983 throw DataException("Error - DataTagged:: appending data type to netCDF file failed.");
914     if (!dataFile.add_att("rank",rank) )
915     throw DataException("Error - DataTagged:: appending rank attribute to netCDF file failed.");
916     if (!dataFile.add_att("function_space_type",type))
917     throw DataException("Error - DataTagged:: appending function space attribute to netCDF file failed.");
918     ndims=rank+1;
919     if ( rank >0 ) {
920     dims[0]=shape[0];
921     if (! (ncdims[0] = dataFile.add_dim("d0",shape[0])) )
922 ksteube 1800 throw DataException("Error - DataTagged:: appending ncdimension 0 to netCDF file failed.");
923 gross 983 }
924     if ( rank >1 ) {
925     dims[1]=shape[1];
926     if (! (ncdims[1] = dataFile.add_dim("d1",shape[1])) )
927 ksteube 1800 throw DataException("Error - DataTagged:: appending ncdimension 1 to netCDF file failed.");
928 gross 983 }
929     if ( rank >2 ) {
930     dims[2]=shape[2];
931     if (! (ncdims[2] = dataFile.add_dim("d2", shape[2])) )
932 ksteube 1800 throw DataException("Error - DataTagged:: appending ncdimension 2 to netCDF file failed.");
933 gross 983 }
934     if ( rank >3 ) {
935     dims[3]=shape[3];
936     if (! (ncdims[3] = dataFile.add_dim("d3", shape[3])) )
937 ksteube 1800 throw DataException("Error - DataTagged:: appending ncdimension 3 to netCDF file failed.");
938 gross 983 }
939     const DataTagged::DataMapType& thisLookup=getTagLookup();
940     DataTagged::DataMapType::const_iterator i;
941     DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
942     int ntags=1;
943     for (i=thisLookup.begin();i!=thisLookupEnd;i++) ntags++;
944 phornby 1628 int* tags =(int*) esysUtils::malloc(ntags*sizeof(int));
945 phornby 1020 int c=1;
946 gross 983 tags[0]=-1;
947     for (i=thisLookup.begin();i!=thisLookupEnd;i++) tags[c++]=i->first;
948     dims[rank]=ntags;
949     if (! (ncdims[rank] = dataFile.add_dim("num_tags", dims[rank])) )
950 phornby 1020 {
951 phornby 1628 esysUtils::free(tags);
952 phornby 1020 throw DataException("Error - DataTagged:: appending num_tags to netCDF file failed.");
953     }
954 gross 983 if (! ( tags_var = dataFile.add_var("tags", ncInt, ncdims[rank])) )
955 phornby 1020 {
956 phornby 1628 esysUtils::free(tags);
957 gross 983 throw DataException("Error - DataTagged:: appending tags to netCDF file failed.");
958 phornby 1020 }
959 gross 983 if (! (tags_var->put(tags,dims[rank])) )
960 phornby 1020 {
961 phornby 1628 esysUtils::free(tags);
962 gross 983 throw DataException("Error - DataTagged:: copy tags to netCDF buffer failed.");
963 phornby 1020 }
964 gross 983 if (! ( var = dataFile.add_var("data", ncDouble, ndims, ncdims)) )
965 phornby 1020 {
966 phornby 1628 esysUtils::free(tags);
967 gross 983 throw DataException("Error - DataTagged:: appending variable to netCDF file failed.");
968 phornby 1020 }
969 gross 1141 if (! (var->put(d_ptr,dims)) )
970 phornby 1020 {
971 phornby 1628 esysUtils::free(tags);
972 gross 983 throw DataException("Error - DataTagged:: copy data to netCDF buffer failed.");
973 phornby 1020 }
974 gross 1023 #else
975     throw DataException("Error - DataTagged:: dump is not configured with netCDF. Please contact your installation manager.");
976     #endif
977 gross 950 }
978 jfenwick 1796
979     DataTypes::ValueType&
980     DataTagged::getVector()
981     {
982     return m_data;
983     }
984    
985     const DataTypes::ValueType&
986     DataTagged::getVector() const
987     {
988     return m_data;
989     }
990    
991 jgs 82 } // end of namespace

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26