/[escript]/trunk/tools/libescriptreader/src/escript2silo.cpp
ViewVC logotype

Contents of /trunk/tools/libescriptreader/src/escript2silo.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 2194 - (show annotations)
Tue Jan 6 01:23:37 2009 UTC (10 years, 7 months ago) by caltinay
File size: 5658 byte(s)
escript2silo: Allow different meshes for different timesteps.
escriptreader: Fixed segfaults for some error conditions.


1
2 /*******************************************************
3 *
4 * Copyright (c) 2003-2008 by University of Queensland
5 * 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 // escript2silo.cpp
16 //
17 #include <escriptreader/MPDataSet.h>
18 #include <escriptreader/MeshWithElements.h>
19 #include <escriptreader/DataVar.h>
20 #include <silo.h>
21 #include <netcdf.hh>
22 #include <cstring>
23 #include <fstream>
24
25 using namespace std;
26 using namespace EscriptReader;
27
28 string insertTimestep(const string& fString, int timeStep)
29 {
30 string s(fString);
31 size_t pos;
32 if ((pos = s.find("%")) != s.npos) {
33 size_t endpos = pos+1;
34 while (endpos<s.length() && s[endpos] != 'd')
35 endpos++;
36 string fmtStr = s.substr(pos, endpos-pos+1);
37 char ts[255];
38 snprintf(ts, 255, fmtStr.c_str(), timeStep);
39 s.replace(pos, endpos-pos+1, ts);
40 }
41 return s;
42 }
43
44 int main(int argc, char** argv)
45 {
46 // turn off for debugging purposes
47 bool writeMultiMesh = true;
48
49 // whether time-varying datasets should use the same mesh (from T=0)
50 // TODO: Add a command line option for this
51 bool writeMeshOnce = true;
52
53 if (argc < 2) {
54 cerr << "Usage: " << argv[0] << " <file.esd>" << endl;
55 return -1;
56 }
57
58 string esdFile(argv[1]);
59
60 ifstream in(esdFile.c_str());
61 if (!in.is_open()) {
62 cerr << "Could not open " << esdFile << "." << endl;
63 return -1;
64 }
65
66 // first line (header) should be "#escript datafile Vx.y"
67 char line[256];
68 in.getline(line, 256);
69 int major, minor;
70 if (sscanf(line, "#escript datafile V%d.%d", &major, &minor) != 2) {
71 cerr << esdFile << " is not a valid escript datafile." << endl;
72 in.close();
73 return -1;
74 }
75
76 int nParts=0, nTimesteps=1;
77 string meshFile;
78 StringVec varFiles;
79 StringVec varNames;
80
81 while (in.good()) {
82 in.getline(line, 256);
83 if (line[0] == '#' || strlen(line) == 0)
84 continue;
85 int iVal;
86 char sVal[256];
87 if (sscanf(line, "N=%d", &iVal) == 1)
88 nParts = iVal;
89 else if (sscanf(line, "T=%d", &iVal) == 1)
90 nTimesteps = iVal;
91 else if (sscanf(line, "M=%s", sVal) == 1)
92 meshFile = sVal;
93 else if (sscanf(line, "V=%s", sVal) == 1 && strchr(sVal, ':')) {
94 // split into filename and variable name
95 char* colon = strchr(sVal, ':');
96 *colon = 0;
97 varFiles.push_back(sVal);
98 varNames.push_back(colon+1);
99 } else {
100 cerr << esdFile << " is not a valid escript datafile." << endl;
101 in.close();
102 return -1;
103 }
104 }
105
106 in.close();
107
108 if (nParts < 1 || meshFile == "" || nTimesteps < 1) {
109 cerr << esdFile << " is not a valid escript datafile." << endl;
110 return -1;
111 }
112
113 cout << "Converting " << esdFile << "..." << endl;
114
115 MeshBlocks meshFromTzero;
116
117 // load and save all timesteps
118 for (int timeStep = 0; timeStep < nTimesteps; timeStep++) {
119 StringVec varFilesTS;
120 StringVec::const_iterator it;
121
122 // look for "%d" in filename and replace by timestep if found
123 for (it=varFiles.begin(); it!=varFiles.end(); it++) {
124 string v = insertTimestep(*it, timeStep);
125 if (nParts > 1)
126 v.append(".nc.%04d");
127 else
128 v.append(".nc");
129 varFilesTS.push_back(v);
130 }
131
132 if (nTimesteps > 1)
133 cout << "T = " << timeStep << endl;
134
135 MPDataSet* ds = new MPDataSet();
136
137 if (writeMeshOnce && timeStep > 0) {
138 if (!ds->load(meshFromTzero, meshFile, varFilesTS, varNames)) {
139 delete ds;
140 break;
141 }
142 } else {
143 string meshTS = insertTimestep(meshFile, timeStep);
144 if (nParts > 1)
145 meshTS.append(".nc.%04d");
146 else
147 meshTS.append(".nc");
148
149 if (!ds->load(meshTS, varFilesTS, varNames, nParts)) {
150 delete ds;
151 break;
152 }
153 }
154
155 string baseName(esdFile);
156 size_t dot = baseName.rfind('.');
157 if (dot != baseName.npos)
158 baseName.erase(dot, baseName.length()-dot);
159
160 ostringstream siloFile;
161 siloFile << baseName;
162 if (nTimesteps > 1)
163 siloFile << "." << timeStep;
164 siloFile << ".silo";
165
166 ds->saveAsSilo(siloFile.str(), writeMultiMesh);
167
168 // keep mesh from first timestep if it should be reused
169 if (writeMeshOnce && nTimesteps > 1 && timeStep == 0) {
170 meshFromTzero = ds->extractMesh();
171 meshFile = siloFile.str();
172 MeshBlocks::iterator meshIt;
173 for (meshIt = meshFromTzero.begin(); meshIt != meshFromTzero.end();
174 meshIt++)
175 {
176 // Prepend Silo mesh paths with the filename of the mesh to be
177 // used
178 string fullSiloPath = meshFile + string(":");
179 fullSiloPath += (*meshIt)->getSiloPath();
180 (*meshIt)->setSiloPath(fullSiloPath);
181 }
182 }
183
184 delete ds;
185 }
186
187 // clean up
188 MeshBlocks::iterator meshIt;
189 for (meshIt = meshFromTzero.begin(); meshIt != meshFromTzero.end(); meshIt++)
190 delete *meshIt;
191
192 cout << "All done." << endl;
193
194 return 0;
195 }
196

  ViewVC Help
Powered by ViewVC 1.1.26