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

Annotation of /trunk/escript/src/DataTagged.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1796 - (hide annotations)
Wed Sep 17 01:45:46 2008 UTC (11 years, 6 months ago) by jfenwick
File size: 34058 byte(s)
Merged noarrayview branch onto trunk.


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