/[escript]/branches/split/weipa/src/DataVar.cpp
ViewVC logotype

Annotation of /branches/split/weipa/src/DataVar.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 2548 - (hide annotations)
Mon Jul 20 06:20:06 2009 UTC (9 years, 8 months ago) by jfenwick
Original Path: trunk/tools/libescriptreader/src/escriptreader/DataVar.cpp
File size: 22803 byte(s)
Updating copyright notices
1 caltinay 2183
2     /*******************************************************
3     *
4 jfenwick 2548 * Copyright (c) 2003-2009 by University of Queensland
5 caltinay 2183 * Earth Systems Science Computational Center (ESSCC)
6     * http://www.uq.edu.au/esscc
7     *
8     * Primary Business: Queensland, Australia
9     * Licensed under the Open Software License version 3.0
10     * http://www.opensource.org/licenses/osl-3.0.php
11     *
12     *******************************************************/
13    
14     //
15     // DataVar.cpp
16     //
17     #include <escriptreader/DataVar.h>
18     #include <escriptreader/ElementData.h>
19     #include <escriptreader/MeshWithElements.h>
20     #include <netcdf.hh>
21     #if HAVE_SILO
22     #include <silo.h>
23     #endif
24    
25     using namespace std;
26    
27 caltinay 2187 namespace EscriptReader {
28    
29 caltinay 2183 enum {
30     NODE_CENTERED = 1,
31     ZONE_CENTERED = 2
32     };
33    
34     //
35     // Constructor
36     //
37     DataVar::DataVar(const string& name) :
38     varName(name), numSamples(0), rank(0), ptsPerSample(0), centering(0),
39     reorderedNumSamples(0), fullMesh(NULL)
40     {
41     }
42    
43     //
44     // Destructor
45     //
46     DataVar::~DataVar()
47     {
48     CoordArray::iterator it;
49     for (it = reorderedData.begin(); it != reorderedData.end(); it++)
50     delete[] *it;
51     for (it = rawData.begin(); it != rawData.end(); it++)
52     delete[] *it;
53     }
54    
55     //
56     // Copy constructor
57     //
58     DataVar::DataVar(const DataVar& d) :
59     varName(d.varName), numSamples(d.numSamples),
60     rank(d.rank), ptsPerSample(d.ptsPerSample), centering(d.centering),
61     funcSpace(d.funcSpace), shape(d.shape), sampleID(d.sampleID),
62     reorderedNumSamples(d.reorderedNumSamples), fullMesh(d.fullMesh)
63     {
64     CoordArray::const_iterator it;
65     for (it = d.rawData.begin(); it != d.rawData.end(); it++) {
66     float* c = new float[numSamples];
67     copy(*it, (*it)+numSamples, c);
68     rawData.push_back(c);
69     }
70     for (it = d.reorderedData.begin(); it != d.reorderedData.end(); it++) {
71     float* c = new float[reorderedNumSamples];
72     copy(*it, (*it)+reorderedNumSamples, c);
73     reorderedData.push_back(c);
74     }
75     }
76    
77     //
78     // Special constructor for mesh data
79     //
80     DataVar::DataVar(const string& name, const IntVec& data,
81     MeshWithElements* mesh) :
82     varName(name)
83     {
84     numSamples = data.size();
85    
86     float* c = new float[numSamples];
87     rawData.push_back(c);
88     IntVec::const_iterator it;
89     for (it=data.begin(); it != data.end(); it++)
90     *c++ = static_cast<float>(*it);
91    
92     rank = 0;
93     ptsPerSample = 1;
94     if (name.compare(0, 6, "Nodes_") == 0) {
95     funcSpace = FINLEY_NODES;
96     centering = NODE_CENTERED;
97     sampleID.insert(sampleID.end(), mesh->getNodeIDs().begin(),
98     mesh->getNodeIDs().end());
99     } else if (name.compare(0, 9, "Elements_") == 0) {
100     funcSpace = FINLEY_ELEMENTS;
101     centering = ZONE_CENTERED;
102     sampleID.insert(sampleID.end(), mesh->getElements()->getIDs().begin(),
103     mesh->getElements()->getIDs().end());
104     } else if (name.compare(0, 13, "FaceElements_") == 0) {
105     funcSpace = FINLEY_FACE_ELEMENTS;
106     centering = ZONE_CENTERED;
107     sampleID.insert(sampleID.end(),
108     mesh->getFaceElements()->getIDs().begin(),
109     mesh->getFaceElements()->getIDs().end());
110     } else if (name.compare(0, 16, "ContactElements_") == 0) {
111     funcSpace = FINLEY_CONTACT_ELEMENTS_1;
112     centering = ZONE_CENTERED;
113     sampleID.insert(sampleID.end(),
114     mesh->getContactElements()->getIDs().begin(),
115     mesh->getContactElements()->getIDs().end());
116     } else if (name.compare(0, 7, "Points_") == 0) {
117     funcSpace = FINLEY_POINTS;
118     centering = NODE_CENTERED;
119     sampleID.insert(sampleID.end(), mesh->getPoints()->getIDs().begin(),
120     mesh->getPoints()->getIDs().end());
121     }
122    
123     shape.clear();
124     reorderedNumSamples = 0;
125     }
126    
127     //
128 caltinay 2196 // Appends raw data including IDs from rhs.
129 caltinay 2183 //
130     bool DataVar::append(const DataVar& rhs)
131     {
132     // check if variables are compatible
133     if (varName != rhs.varName || ptsPerSample != rhs.ptsPerSample ||
134     rank != rhs.rank || shape.size() != rhs.shape.size() ||
135     centering != rhs.centering)
136     return false;
137    
138     for (size_t i=0; i<shape.size(); i++)
139     if (shape[i] != rhs.shape[i])
140     return false;
141    
142     sampleID.insert(sampleID.end(), rhs.sampleID.begin(), rhs.sampleID.end());
143     for (size_t i=0; i<rawData.size(); i++) {
144     float* c = new float[numSamples+rhs.numSamples];
145     copy(rawData[i], rawData[i]+numSamples, c);
146     copy(rhs.rawData[i], rhs.rawData[i]+rhs.numSamples, c+numSamples);
147     delete[] rawData[i];
148     rawData[i] = c;
149     }
150     numSamples += rhs.numSamples;
151    
152     // invalidate previously reordered data
153     CoordArray::iterator it;
154     for (it = reorderedData.begin(); it != reorderedData.end(); it++)
155     delete[] *it;
156     reorderedData.clear();
157     reorderedNumSamples = 0;
158    
159     return true;
160     }
161    
162     //
163     // Returns a subset of the src array according to stride parameter.
164     // If samples consist of multiple values they are averaged beforehand.
165     // Used to separate (x0,y0,z0,x1,y1,z1,...) into (x0,x1,...), (y0,y1,...) and
166     // (z0,z1,...)
167     //
168     float* DataVar::averageData(const float* src, size_t stride)
169     {
170     float* res = new float[numSamples];
171    
172     if (ptsPerSample == 1) {
173     float* dest = res;
174     for (int i=0; i<numSamples; i++, src+=stride)
175     *dest++ = *src;
176     } else {
177     float* dest = res;
178     for (int i=0; i<numSamples; i++) {
179     double tmpVal = 0.0;
180     for (int j=0; j<ptsPerSample; j++, src+=stride)
181     tmpVal += *src;
182     *dest++ = (float)(tmpVal / ptsPerSample);
183     }
184     }
185     return res;
186     }
187    
188     //
189     // Reads scalar data (rank=0) from NetCDF file and stores the values
190     // after averaging.
191     //
192     void DataVar::readRank0Data(NcFile* ncfile)
193     {
194     shape.clear();
195     float* tempData = new float[ptsPerSample*numSamples];
196     NcVar* var = ncfile->get_var("data");
197     var->get(tempData, ptsPerSample, numSamples);
198    
199     float* c = averageData(tempData, 1);
200     rawData.push_back(c);
201    
202     delete[] tempData;
203     }
204    
205     //
206     // Reads vector data (rank=1) from NetCDF file and stores the components
207     // separately after averaging.
208     //
209     void DataVar::readRank1Data(NcFile* ncfile)
210     {
211     shape.clear();
212     NcDim* dim = ncfile->get_dim("d0");
213     shape.push_back(dim->size());
214    
215     float* tempData = new float[shape[0]*ptsPerSample*numSamples];
216     NcVar* var = ncfile->get_var("data");
217     var->get(tempData, shape[0], ptsPerSample, numSamples);
218    
219     for (int i=0; i<shape[0]; i++) {
220     const float* src = tempData;
221     src+=i;
222     float* c = averageData(src, shape[0]);
223     rawData.push_back(c);
224     }
225     delete[] tempData;
226     }
227    
228     //
229     // Like readRank1Data() but for tensor data (rank=2).
230     //
231     void DataVar::readRank2Data(NcFile* ncfile)
232     {
233     shape.clear();
234     NcDim* dim = ncfile->get_dim("d0");
235     shape.push_back(dim->size());
236     dim = ncfile->get_dim("d1");
237     shape.push_back(dim->size());
238    
239     float* tempData = new float[shape[0]*shape[1]*ptsPerSample*numSamples];
240     NcVar* var = ncfile->get_var("data");
241     var->get(tempData, shape[0], shape[1], ptsPerSample, numSamples);
242    
243     for (int i=0; i<shape[1]; i++) {
244     for (int j=0; j<shape[0]; j++) {
245     const float* src = tempData;
246     src+=i*shape[0]+j;
247     float* c = averageData(src, shape[0]*shape[1]);
248     rawData.push_back(c);
249     }
250     }
251     delete[] tempData;
252     }
253    
254     //
255     // Reads a NetCDF file in escript/finley format.
256     //
257 caltinay 2196 bool DataVar::readFromNc(const string& filename)
258 caltinay 2183 {
259     NcError ncerr(NcError::silent_nonfatal);
260 caltinay 2196 NcFile* input = new NcFile(filename.c_str());
261 caltinay 2183 if (!input->is_valid()) {
262 caltinay 2196 cerr << "Could not open input file " << filename << "." << endl;
263 caltinay 2183 delete input;
264     return false;
265     }
266    
267     NcDim* dim;
268     NcAtt* att;
269    
270     dim = input->get_dim("num_samples");
271     numSamples = dim->size();
272    
273     // if there are no data samples bail out
274     if (numSamples == 0) {
275     delete input;
276     return false;
277     }
278    
279     att = input->get_att("type_id");
280     int typeID = att->as_int(0);
281     if (typeID != 2) {
282     cerr << "WARNING: Only expanded data supported at the moment!" << endl;
283     delete input;
284     return false;
285     }
286    
287     att = input->get_att("rank");
288     rank = att->as_int(0);
289    
290     dim = input->get_dim("num_data_points_per_sample");
291     ptsPerSample = dim->size();
292    
293     att = input->get_att("function_space_type");
294     funcSpace = att->as_int(0);
295    
296     #ifdef _DEBUG
297     cout << varName << ":\t" << numSamples << " samples, "
298     << ptsPerSample << " pts/s, rank: " << rank << endl;
299     #endif
300    
301     sampleID.clear();
302     sampleID.insert(sampleID.end(), numSamples, 0);
303     NcVar* var = input->get_var("id");
304     var->get(&sampleID[0], numSamples);
305    
306     switch (rank) {
307     case 0:
308     readRank0Data(input);
309     break;
310     case 1:
311     readRank1Data(input);
312     break;
313     case 2:
314     readRank2Data(input);
315     break;
316     default:
317     cerr << "WARNING: Rank " << rank << " data is not supported!\n";
318     delete input;
319     return false;
320     }
321    
322     delete input;
323     return true;
324     }
325    
326     //
327     // Returns one of the mesh names provided by mainMesh that matches the
328     // data variable's function space type and reduced/unreduced state.
329     //
330     string DataVar::getMeshName(MeshWithElements* mainMesh) const
331     {
332     string name;
333    
334     switch (funcSpace) {
335     case FINLEY_REDUCED_NODES:
336     case FINLEY_REDUCED_DEGREES_OF_FREEDOM:
337     case FINLEY_REDUCED_ELEMENTS:
338     case FINLEY_ELEMENTS:
339     if (mainMesh->getElements()->reducedCount > 0) {
340     name = mainMesh->getElements()->reducedMesh->getName();
341     } else {
342     name = mainMesh->getElements()->fullMesh->getName();
343     }
344     break;
345    
346     case FINLEY_REDUCED_FACE_ELEMENTS:
347     case FINLEY_FACE_ELEMENTS:
348     if (mainMesh->getFaceElements()->reducedCount > 0) {
349     name = mainMesh->getFaceElements()->reducedMesh->getName();
350     } else {
351     name = mainMesh->getFaceElements()->fullMesh->getName();
352     }
353     break;
354    
355     case FINLEY_REDUCED_CONTACT_ELEMENTS_1:
356     case FINLEY_REDUCED_CONTACT_ELEMENTS_2:
357     case FINLEY_CONTACT_ELEMENTS_1:
358     case FINLEY_CONTACT_ELEMENTS_2:
359     if (mainMesh->getContactElements()->reducedCount > 0) {
360     name = mainMesh->getContactElements()->reducedMesh->getName();
361     } else {
362     name = mainMesh->getContactElements()->fullMesh->getName();
363     }
364     break;
365    
366     case FINLEY_NODES:
367     case FINLEY_DEGREES_OF_FREEDOM:
368     name = mainMesh->getElements()->fullMesh->getName();
369     break;
370    
371     case FINLEY_POINTS:
372     name = mainMesh->getPoints()->fullMesh->getName();
373     break;
374     }
375     return name;
376     }
377    
378     //
379     // Returns true if the data values are nodal, false if they are zonal.
380     //
381     bool DataVar::isNodeCentered() const
382     {
383     return (funcSpace == FINLEY_REDUCED_NODES ||
384     funcSpace == FINLEY_REDUCED_DEGREES_OF_FREEDOM ||
385     funcSpace == FINLEY_NODES ||
386     funcSpace == FINLEY_DEGREES_OF_FREEDOM ||
387     funcSpace == FINLEY_POINTS);
388     }
389    
390     //
391     // Filters and reorders the raw sample values according to the IDs provided
392     // in 'requiredIDs'. This is used to have data arrays ordered according to
393     // the underlying mesh (i.e. DataID[i]==MeshNodeID[i])
394     //
395     void DataVar::reorderSamples(const IndexMap& id2idxMap,
396     const IntVec& requiredIDs)
397     {
398     CoordArray::iterator it;
399     for (it = reorderedData.begin(); it != reorderedData.end(); it++)
400     delete[] *it;
401     reorderedData.clear();
402    
403     buildIndexMap();
404     for (size_t i=0; i < rawData.size(); i++) {
405     float* c = new float[reorderedNumSamples];
406     reorderedData.push_back(c);
407     const float* src = rawData[i];
408     IntVec::const_iterator idIt = requiredIDs.begin();
409     for (; idIt != requiredIDs.end(); idIt++) {
410     size_t srcIdx = sampleID2idx.find(*idIt)->second;
411     size_t destIdx = id2idxMap.find(*idIt)->second;
412     c[destIdx] = src[srcIdx];
413     }
414     }
415     }
416    
417     //
418     // For zonal data this method reorders the values according to the indices
419     // given in reorderArray. This is used to move ghost zones to the end of
420     // the arrays which conforms to Silo's expected format.
421     // Nodal data is not changed by this method.
422     //
423     void DataVar::handleGhostZones(const IntVec& reorderArray)
424     {
425     if (centering == NODE_CENTERED)
426     return;
427    
428     vector<float*>::iterator it;
429     for (it = reorderedData.begin(); it!=reorderedData.end(); it++) {
430     float* original = *it;
431     float* reordered = new float[reorderedNumSamples];
432     float* arrIt = reordered;
433     IntVec::const_iterator idxIt;
434     for (idxIt=reorderArray.begin(); idxIt!=reorderArray.end(); idxIt++)
435     *arrIt++ = original[*idxIt];
436    
437     delete[] *it;
438     *it = reordered;
439     }
440     }
441    
442     //
443     // Makes the top-level mesh known to this data variable. The mesh is used
444     // to reorder and filter the samples and inquire the number of ghost zones.
445     //
446     bool DataVar::setMesh(MeshWithElements* mesh)
447     {
448     if (fullMesh == mesh)
449     return true;
450    
451     const IndexMap* id2idxMap;
452     const IntVec* reqIDs;
453     const IntVec* reorderArray = NULL;
454    
455     switch (funcSpace) {
456     case FINLEY_REDUCED_NODES:
457     case FINLEY_REDUCED_DEGREES_OF_FREEDOM:
458     {
459     centering = NODE_CENTERED;
460     ElementData* cells = mesh->getElements();
461     if (cells->reducedCount > 0) {
462     if (cells->getReducedGhostCount())
463     reorderArray = &cells->reducedIndexArray;
464     siloMeshName = cells->reducedMesh->getFullSiloName();
465     id2idxMap = &cells->reducedMesh->getIndexMap();
466     reqIDs = &cells->reducedMesh->getNodeIDs();
467     reorderedNumSamples = cells->reducedMesh->getNumNodes();
468     } else {
469     if (cells->getGhostCount())
470     reorderArray = &cells->indexArray;
471     siloMeshName = cells->fullMesh->getFullSiloName();
472     id2idxMap = &cells->fullMesh->getIndexMap();
473     reqIDs = &cells->fullMesh->getNodeIDs();
474     reorderedNumSamples = cells->fullMesh->getNumNodes();
475     }
476     }
477     break;
478    
479     case FINLEY_NODES:
480     case FINLEY_DEGREES_OF_FREEDOM:
481     {
482     centering = NODE_CENTERED;
483     ElementData* cells = mesh->getElements();
484     if (cells->getGhostCount())
485     reorderArray = &cells->indexArray;
486     siloMeshName = cells->fullMesh->getFullSiloName();
487     id2idxMap = &cells->fullMesh->getIndexMap();
488     reqIDs = &cells->fullMesh->getNodeIDs();
489     reorderedNumSamples = cells->fullMesh->getNumNodes();
490     }
491     break;
492    
493     case FINLEY_REDUCED_ELEMENTS:
494     case FINLEY_ELEMENTS:
495     {
496     centering = ZONE_CENTERED;
497     ElementData* cells = mesh->getElements();
498     id2idxMap = &cells->ID2idx;
499     reqIDs = &cells->getIDs();
500     if (cells->reducedCount > 0) {
501     if (cells->getReducedGhostCount())
502     reorderArray = &cells->reducedIndexArray;
503     reorderedNumSamples = cells->reducedCount;
504     siloMeshName = cells->reducedMesh->getFullSiloName();
505     } else {
506     if (cells->getGhostCount())
507     reorderArray = &cells->indexArray;
508     reorderedNumSamples = cells->count;
509     siloMeshName = cells->fullMesh->getFullSiloName();
510     }
511     }
512     break;
513    
514     case FINLEY_REDUCED_FACE_ELEMENTS:
515     case FINLEY_FACE_ELEMENTS:
516     {
517     centering = ZONE_CENTERED;
518     ElementData* cells = mesh->getFaceElements();
519     id2idxMap = &cells->ID2idx;
520     reqIDs = &cells->getIDs();
521     if (cells->reducedCount > 0) {
522     if (cells->getReducedGhostCount())
523     reorderArray = &cells->reducedIndexArray;
524     reorderedNumSamples = cells->reducedCount;
525     siloMeshName = cells->reducedMesh->getFullSiloName();
526     } else {
527     if (cells->getGhostCount())
528     reorderArray = &cells->indexArray;
529     reorderedNumSamples = cells->count;
530     siloMeshName = cells->fullMesh->getFullSiloName();
531     }
532     }
533     break;
534    
535     case FINLEY_REDUCED_CONTACT_ELEMENTS_1:
536     case FINLEY_REDUCED_CONTACT_ELEMENTS_2:
537     case FINLEY_CONTACT_ELEMENTS_1:
538     case FINLEY_CONTACT_ELEMENTS_2:
539     {
540     centering = ZONE_CENTERED;
541     ElementData* cells = mesh->getContactElements();
542     id2idxMap = &cells->ID2idx;
543     reqIDs = &cells->getIDs();
544     if (cells->reducedCount > 0) {
545     if (cells->getReducedGhostCount())
546     reorderArray = &cells->reducedIndexArray;
547     reorderedNumSamples = cells->reducedCount;
548     siloMeshName = cells->reducedMesh->getFullSiloName();
549     } else {
550     if (cells->getGhostCount())
551     reorderArray = &cells->indexArray;
552     reorderedNumSamples = cells->count;
553     siloMeshName = cells->fullMesh->getFullSiloName();
554     }
555     }
556     break;
557    
558     case FINLEY_POINTS:
559     {
560     centering = NODE_CENTERED;
561     ElementData* cells = mesh->getPoints();
562     if (cells->getGhostCount())
563     reorderArray = &cells->indexArray;
564     siloMeshName = cells->fullMesh->getFullSiloName();
565     id2idxMap = &cells->ID2idx;
566     reqIDs = &cells->getIDs();
567     reorderedNumSamples = cells->count;
568     }
569     break;
570    
571     default:
572     cerr << "Unknown function space type " << funcSpace << "!\n";
573     return false;
574     }
575    
576     if (reorderedNumSamples > numSamples) {
577     cerr << "WARNING: " << varName << " has " << numSamples
578     << " instead of " << reorderedNumSamples << " samples!" << endl;
579     return false;
580     }
581    
582     fullMesh = mesh;
583    
584     reorderSamples(*id2idxMap, *reqIDs);
585     if (reorderArray)
586     handleGhostZones(*reorderArray);
587     return true;
588     }
589    
590     ///////////////////////////////
591     // SILO related methods follow
592     ///////////////////////////////
593    
594     //
595     // If the data is tensor data then the components of the tensor are stored
596     // separately in the Silo file. This method then returns a string that
597     // contains the proper Silo expression to put the tensor together again.
598     // For non-tensor data this method returns an empty string.
599     //
600     string DataVar::getTensorDef() const
601     {
602     if (rank < 2)
603     return string();
604    
605     /// Format string for Silo 2x2 tensor
606     const string tensor2DefFmt =
607     "{{ <%sa_00>, <%sa_01> },"
608     " { <%sa_10>, <%sa_11> }}";
609    
610     /// Format string for Silo 3x3 tensor
611     const string tensor3DefFmt =
612     "{{ <%sa_00>, <%sa_01>, <%sa_02> },"
613     " { <%sa_10>, <%sa_11>, <%sa_12> },"
614     " { <%sa_20>, <%sa_21>, <%sa_22> }}";
615    
616     string tensorDef;
617     string tensorDir = varName+string("_comps/");
618     if (shape[1] == 3) {
619     char* tDef = new char[tensor3DefFmt.length()+9*tensorDir.length()];
620     sprintf(tDef, tensor3DefFmt.c_str(),
621     tensorDir.c_str(), tensorDir.c_str(), tensorDir.c_str(),
622     tensorDir.c_str(), tensorDir.c_str(), tensorDir.c_str(),
623     tensorDir.c_str(), tensorDir.c_str(), tensorDir.c_str());
624     tensorDef = tDef;
625     delete[] tDef;
626     } else {
627     char* tDef = new char[tensor2DefFmt.length()+4*tensorDir.length()];
628     sprintf(tDef, tensor2DefFmt.c_str(),
629     tensorDir.c_str(), tensorDir.c_str(),
630     tensorDir.c_str(), tensorDir.c_str(),
631     tensorDir.c_str(), tensorDir.c_str());
632     tensorDef = tDef;
633     delete[] tDef;
634     }
635     return tensorDef;
636     }
637    
638     //
639     // Writes the data to given Silo file under the virtual path provided.
640     // The corresponding mesh must have been written already and made known
641     // to this variable by a call to setMesh().
642     //
643     bool DataVar::writeToSilo(DBfile* dbfile, const string& siloPath)
644     {
645     #if HAVE_SILO
646     if (numSamples == 0)
647     return true;
648    
649     // have to set mesh first
650     if (!fullMesh)
651     return false;
652    
653     int ret;
654    
655     if (siloPath != "") {
656     ret = DBSetDir(dbfile, siloPath.c_str());
657     if (ret != 0)
658     return false;
659     }
660    
661     char* meshName = (char*)siloMeshName.c_str();
662     int dcenter = (centering == NODE_CENTERED ? DB_NODECENT : DB_ZONECENT);
663    
664     if (rank == 0) {
665     ret = DBPutUcdvar1(dbfile, varName.c_str(), meshName, reorderedData[0],
666     reorderedNumSamples, NULL, 0, DB_FLOAT, dcenter, NULL);
667     }
668     else if (rank == 1) {
669     const string comps[3] = {
670     varName+string("_x"), varName+string("_y"), varName+string("_z")
671     };
672     const char* varnames[3] = {
673     comps[0].c_str(), comps[1].c_str(), comps[2].c_str()
674     };
675    
676     ret = DBPutUcdvar(dbfile, varName.c_str(), meshName, shape[0],
677     (char**)varnames, &reorderedData[0], reorderedNumSamples, NULL,
678     0, DB_FLOAT, dcenter, NULL);
679     }
680     else {
681     string tensorDir = varName+string("_comps/");
682     ret = DBMkdir(dbfile, tensorDir.c_str());
683     if (ret == 0) {
684     int one = 1;
685     DBoptlist* optList = DBMakeOptlist(1);
686     DBAddOption(optList, DBOPT_HIDE_FROM_GUI, &one);
687    
688     for (int i=0; i<shape[1]; i++) {
689     for (int j=0; j<shape[0]; j++) {
690     ostringstream varname;
691     varname << tensorDir << "a_" << i << j;
692     ret = DBPutUcdvar1(dbfile, varname.str().c_str(), meshName,
693     reorderedData[i*shape[0]+j], reorderedNumSamples,
694     NULL, 0, DB_FLOAT, dcenter, optList);
695     if (ret != 0) break;
696     }
697     if (ret != 0) break;
698     }
699     DBFreeOptlist(optList);
700     } // ret==0
701     } // rank
702    
703     DBSetDir(dbfile, "/");
704     return (ret == 0);
705    
706     #else // !HAVE_SILO
707     return false;
708     #endif
709     }
710    
711 caltinay 2187 } // namespace EscriptReader
712    

  ViewVC Help
Powered by ViewVC 1.1.26