OPAL (Object Oriented Parallel Accelerator Library) 2022.1
OPAL
H5PartWrapperForPC.cpp
Go to the documentation of this file.
1//
2// Copyright & License: See Copyright.readme in src directory
3//
4
5#include "Structure/H5PartWrapperForPC.h"
6
7#include "AbstractObjects/OpalData.h"
8#include "Algorithms/PartBunchBase.h"
9#include "OPALconfig.h"
10#include "Physics/Physics.h"
11#include "Physics/Units.h"
12#include "Utilities/Options.h"
13#include "Utilities/Util.h"
14
15#include <sstream>
16#include <set>
17
18H5PartWrapperForPC::H5PartWrapperForPC(const std::string &fileName, h5_int32_t flags):
19 H5PartWrapper(fileName, flags),
20 previousH5Local_m(false),
21 referenceMomentum_m(0.0),
22 referenceLocation_m(0.0),
23 meanE_m(0.0),
24 meanMomentum_m(0.0),
25 azimuth_m(0.0),
26 elevation_m(0.0)
27{ }
28
29H5PartWrapperForPC::H5PartWrapperForPC(const std::string &fileName,
30 int restartStep,
31 std::string sourceFile,
32 h5_int32_t flags):
33 H5PartWrapper(fileName, restartStep, sourceFile, flags),
34 previousH5Local_m(false),
35 referenceMomentum_m(0.0),
36 referenceLocation_m(0.0),
37 meanE_m(0.0),
38 meanMomentum_m(0.0),
39 azimuth_m(0.0),
40 elevation_m(0.0)
41{ }
42
43H5PartWrapperForPC::~H5PartWrapperForPC()
44{ }
45
46void H5PartWrapperForPC::readHeader() {
47 h5_int64_t numFileAttributes = H5GetNumFileAttribs(file_m);
48
49 const h5_size_t lengthAttributeName = 256;
50 char attributeName[lengthAttributeName];
51 h5_int64_t attributeType;
52 h5_size_t numAttributeElements;
53 std::set<std::string> attributeNames;
54
55 for (h5_int64_t i = 0; i < numFileAttributes; ++ i) {
56 REPORTONERROR(H5GetFileAttribInfo(file_m,
57 i,
58 attributeName,
59 lengthAttributeName,
60 &attributeType,
61 &numAttributeElements));
62
63 attributeNames.insert(attributeName);
64 }
65
66 if (attributeNames.find("dump frequency") != attributeNames.end()) {
67 h5_int64_t dumpfreq;
68 READFILEATTRIB(Int64, file_m, "dump frequency", &dumpfreq);
69 OpalData::getInstance()->setRestartDumpFreq(dumpfreq);
70 }
71}
72
73void H5PartWrapperForPC::readStep(PartBunchBase<double, 3>* bunch,
74 h5_ssize_t firstParticle,
75 h5_ssize_t lastParticle)
76{
77 h5_ssize_t numStepsInSource = H5GetNumSteps(file_m);
78 h5_ssize_t readStep = numStepsInSource - 1;
79 REPORTONERROR(H5SetStep(file_m, readStep));
80
81 readStepHeader(bunch);
82 readStepData(bunch, firstParticle, lastParticle);
83}
84
85void H5PartWrapperForPC::readStepHeader(PartBunchBase<double, 3>* bunch) {
86 h5_float64_t pathLength;
87 READSTEPATTRIB(Float64, file_m, "SPOS", &pathLength);
88 bunch->set_sPos(pathLength);
89
90 h5_int64_t ltstep;
91 READSTEPATTRIB(Int64, file_m, "LocalTrackStep", &ltstep);
92 bunch->setLocalTrackStep((long long)ltstep);
93
94 h5_int64_t gtstep;
95 READSTEPATTRIB(Int64, file_m, "GlobalTrackStep", &gtstep);
96 bunch->setGlobalTrackStep((long long)gtstep);
97
98 READSTEPATTRIB(Float64, file_m, "ENERGY", &meanE_m);
99
100 double actualT;
101 READSTEPATTRIB(Float64, file_m, "TIME", &actualT);
102 bunch->setT(actualT);
103
104 h5_int64_t SteptoLastInj;
105 READSTEPATTRIB(Int64, file_m, "SteptoLastInj", &SteptoLastInj);
106 bunch->setSteptoLastInj((int)SteptoLastInj);
107
108 h5_int64_t numBunch;
109 READSTEPATTRIB(Int64, file_m, "NumBunch", &numBunch);
110 bunch->setNumBunch((short)numBunch);
111
112 if (predecessorOPALFlavour_m == "opal-cycl") {
113 READSTEPATTRIB(Float64, file_m, "REFPR", (h5_float64_t*) &referenceMomentum_m[0]);
114 READSTEPATTRIB(Float64, file_m, "REFPT", (h5_float64_t*) &referenceMomentum_m[1]);
115 READSTEPATTRIB(Float64, file_m, "REFPZ", (h5_float64_t*) &referenceMomentum_m[2]);
116
117 READSTEPATTRIB(Float64, file_m, "REFR", (h5_float64_t*) &referenceLocation_m[0]);
118 READSTEPATTRIB(Float64, file_m, "REFTHETA", (h5_float64_t*) &referenceLocation_m[1]);
119 READSTEPATTRIB(Float64, file_m, "REFZ", (h5_float64_t*) &referenceLocation_m[2]);
120
121 READSTEPATTRIB(Float64, file_m, "AZIMUTH", &azimuth_m);
122 READSTEPATTRIB(Float64, file_m, "ELEVATION", &elevation_m);
123
124 h5_int64_t localDump = 0;
125 h5_int64_t rc = H5ReadStepAttribInt64(file_m, "LOCAL", &localDump);
126 if(rc != H5_SUCCESS) {
127
128 ERRORMSG("H5 rc= " << rc << " in " << __FILE__ << " @ line " << __LINE__ << endl);
129
130 std::string errormsg =
131 "You are trying to restart from a legacy file that doesn't contain"
132 "information on local/global frame. We are working on legacy support,"
133 "but for now you have to use OPAL 1.3.0!";
134
135 throw OpalException("H5PartWrapperForPC::readStepHeader", errormsg);
136 } else {
137 if (localDump == 1) previousH5Local_m = true;
138 }
139 } else {
140 bunch->setT(0.0);
141 bunch->setLocalTrackStep((long long) 0 );
142 }
143
144 double mass = bunch->getM() * Units::eV2MeV;
145 meanMomentum_m = sqrt(std::pow(meanE_m,2.0) + 2 * meanE_m * mass) / mass;
146}
147
148void H5PartWrapperForPC::readStepData(PartBunchBase<double, 3>* bunch,
149 h5_ssize_t firstParticle,
150 h5_ssize_t lastParticle)
151{
152 h5_ssize_t numParticles = getNumParticles();
153 if (lastParticle >= numParticles || firstParticle > lastParticle) {
154 throw OpalException("H5PartWrapperForPC::readStepData",
155 "the provided particle numbers don't match the number of particles in the file");
156 }
157 const unsigned int yIndex = (predecessorOPALFlavour_m == "opal-t"? 2: 1);
158 const unsigned int zIndex = (predecessorOPALFlavour_m == "opal-t"? 1: 2);
159
160 REPORTONERROR(H5PartSetView(file_m, firstParticle, lastParticle));
161
162 numParticles = lastParticle - firstParticle + 1;
163
164 std::vector<char> buffer(numParticles * sizeof(h5_float64_t));
165 char* buffer_ptr = Util::c_data(buffer);
166 h5_float64_t *f64buffer = reinterpret_cast<h5_float64_t*>(buffer_ptr);
167 h5_int64_t *i64buffer = reinterpret_cast<h5_int64_t*>(buffer_ptr);
168
169 READDATA(Float64, file_m, "x", f64buffer);
170 for(long int n = 0; n < numParticles; ++ n) {
171 bunch->R[n](0) = f64buffer[n];
172 bunch->Bin[n] = 0;
173 }
174
175 READDATA(Float64, file_m, "y", f64buffer);
176 for(long int n = 0; n < numParticles; ++ n) {
177 bunch->R[n](yIndex) = f64buffer[n];
178 }
179
180 READDATA(Float64, file_m, "z", f64buffer);
181 for(long int n = 0; n < numParticles; ++ n) {
182 bunch->R[n](zIndex) = f64buffer[n];
183 }
184
185 READDATA(Float64, file_m, "px", f64buffer);
186 for(long int n = 0; n < numParticles; ++ n) {
187 bunch->P[n](0) = f64buffer[n];
188 }
189
190 READDATA(Float64, file_m, "py", f64buffer);
191 for(long int n = 0; n < numParticles; ++ n) {
192 bunch->P[n](yIndex) = f64buffer[n];
193 }
194
195 READDATA(Float64, file_m, "pz", f64buffer);
196 for(long int n = 0; n < numParticles; ++ n) {
197 bunch->P[n](zIndex) = f64buffer[n];
198 }
199
200 READDATA(Float64, file_m, "q", f64buffer);
201 for(long int n = 0; n < numParticles; ++ n) {
202 bunch->Q[n] = f64buffer[n];
203 }
204
205 READDATA(Float64, file_m, "mass", f64buffer);
206 for(long int n = 0; n < numParticles; ++ n) {
207 bunch->M[n] = f64buffer[n];
208 }
209
210 if ( bunch->getNumBunch() > 1 ) {
211 READDATA(Int64, file_m, "bin", i64buffer);
212 for(long int n = 0; n < numParticles; ++ n) {
213 bunch->Bin[n] = i64buffer[n];
214 }
215 }
216
217 READDATA(Int64, file_m, "bunchNumber", i64buffer);
218 for(long int n = 0; n < numParticles; ++ n) {
219 bunch->bunchNum[n] = i64buffer[n];
220 }
221
222 REPORTONERROR(H5PartSetView(file_m, -1, -1));
223}
224
225void H5PartWrapperForPC::writeHeader() {
226 std::stringstream OPAL_version;
227 OPAL_version << OPAL_PROJECT_NAME << " " << OPAL_PROJECT_VERSION << " # git rev. " << Util::getGitRevision();
228 WRITESTRINGFILEATTRIB(file_m, "OPAL_version", OPAL_version.str().c_str());
229
230 // Units for the step data
231 WRITESTRINGFILEATTRIB(file_m, "xUnit", "m");
232 WRITESTRINGFILEATTRIB(file_m, "yUnit", "m");
233 WRITESTRINGFILEATTRIB(file_m, "zUnit", "m");
234 WRITESTRINGFILEATTRIB(file_m, "pxUnit", "#beta#gamma");
235 WRITESTRINGFILEATTRIB(file_m, "pyUnit", "#beta#gamma");
236 WRITESTRINGFILEATTRIB(file_m, "pzUnit", "#beta#gamma");
237 WRITESTRINGFILEATTRIB(file_m, "qUnit", "C");
238 WRITESTRINGFILEATTRIB(file_m, "massUnit", "GeV");
239 WRITESTRINGFILEATTRIB(file_m, "idUnit", "1");
240 WRITESTRINGFILEATTRIB(file_m, "binUnit", "1");
241 WRITESTRINGFILEATTRIB(file_m, "bunchNumberUnit", "1");
242 if (Options::ebDump) {
243 WRITESTRINGFILEATTRIB(file_m, "Ex", "MV/m");
244 WRITESTRINGFILEATTRIB(file_m, "Ey", "MV/m");
245 WRITESTRINGFILEATTRIB(file_m, "Ez", "MV/m");
246 WRITESTRINGFILEATTRIB(file_m, "Bx", "T");
247 WRITESTRINGFILEATTRIB(file_m, "By", "T");
248 WRITESTRINGFILEATTRIB(file_m, "Bz", "T");
249 }
250
251 // Units for the step attributes
252 WRITESTRINGFILEATTRIB(file_m, "TIMEUnit", "s");
253 WRITESTRINGFILEATTRIB(file_m, "SPOSUnit", "m");
254
255 WRITESTRINGFILEATTRIB(file_m, "RefPartRUnit", "m");
256 WRITESTRINGFILEATTRIB(file_m, "centroidUnit", "m");
257 WRITESTRINGFILEATTRIB(file_m, "RMSXUnit", "m");
258 WRITESTRINGFILEATTRIB(file_m, "minXUnit", "m");
259 WRITESTRINGFILEATTRIB(file_m, "maxXUnit", "m");
260
261 WRITESTRINGFILEATTRIB(file_m, "RefPartPUnit", "#beta#gamma");
262 WRITESTRINGFILEATTRIB(file_m, "RMSPUnit", "#beta#gamma");
263 WRITESTRINGFILEATTRIB(file_m, "minPUnit", "#beta#gamma");
264 WRITESTRINGFILEATTRIB(file_m, "maxPUnit", "#beta#gamma");
265
266 WRITESTRINGFILEATTRIB(file_m, "#varepsilonUnit", "m rad");
267 WRITESTRINGFILEATTRIB(file_m, "#varepsilon-geomUnit", "m rad");
268
269 WRITESTRINGFILEATTRIB(file_m, "StepUnit", "1");
270 WRITESTRINGFILEATTRIB(file_m, "LocalTrackStepUnit", "1");
271 WRITESTRINGFILEATTRIB(file_m, "GlobalTrackStepUnit", "1");
272
273 WRITESTRINGFILEATTRIB(file_m, "#sigmaUnit", "1");
274
275 WRITESTRINGFILEATTRIB(file_m, "ENERGYUnit", "MeV");
276 WRITESTRINGFILEATTRIB(file_m, "dEUnit", "MeV");
277
278 WRITESTRINGFILEATTRIB(file_m, "MASSUnit", "GeV");
279 WRITESTRINGFILEATTRIB(file_m, "CHARGEUnit", "C");
280
281 WRITESTRINGFILEATTRIB(file_m, "NumBunchUnit", "1");
282 WRITESTRINGFILEATTRIB(file_m, "SteptoLastInjUnit", "1");
283 WRITESTRINGFILEATTRIB(file_m, "LOCALUnit", "1");
284
285 // additional step attributes
286 WRITESTRINGFILEATTRIB(file_m, "REFPRUnit", "#beta#gamma");
287 WRITESTRINGFILEATTRIB(file_m, "REFPTUnit", "#beta#gamma");
288 WRITESTRINGFILEATTRIB(file_m, "REFPZUnit", "#beta#gamma");
289 WRITESTRINGFILEATTRIB(file_m, "REFRUnit", "mm");
290 WRITESTRINGFILEATTRIB(file_m, "REFTHETAUnit", "deg");
291 WRITESTRINGFILEATTRIB(file_m, "REFZUnit", "mm");
292 WRITESTRINGFILEATTRIB(file_m, "REFAZIMUTHUnit", "deg");
293 WRITESTRINGFILEATTRIB(file_m, "REFELEVATIONUnit", "deg");
294
295 WRITESTRINGFILEATTRIB(file_m, "spos-headUnit", "m");
296 WRITESTRINGFILEATTRIB(file_m, "spos-refUnit", "m");
297 WRITESTRINGFILEATTRIB(file_m, "spos-tailUnit", "m");
298
299 WRITESTRINGFILEATTRIB(file_m, "B-refUnit", "T");
300 WRITESTRINGFILEATTRIB(file_m, "E-refUnit", "MV/m");
301 WRITESTRINGFILEATTRIB(file_m, "B-headUnit", "T");
302 WRITESTRINGFILEATTRIB(file_m, "E-headUnit", "MV/m");
303 WRITESTRINGFILEATTRIB(file_m, "B-tailUnit", "T");
304 WRITESTRINGFILEATTRIB(file_m, "E-tailUnit", "MV/m");
305
307 h5_int64_t dumpfreq = Options::psDumpFreq;
308 WRITEFILEATTRIB(Int64, file_m, "dump frequency", &dumpfreq, 1);
309}
310
311void H5PartWrapperForPC::writeStep(PartBunchBase<double, 3>* bunch,
312 const std::map<std::string, double> &additionalStepAttributes)
313{
314 if (bunch->getTotalNum() == 0) return;
315
316 open(H5_O_APPENDONLY);
317 writeStepHeader(bunch, additionalStepAttributes);
318 writeStepData(bunch);
319 close();
320
321 ++ numSteps_m;
322}
323
324void H5PartWrapperForPC::writeStepHeader(PartBunchBase<double, 3>* bunch,
325 const std::map<std::string, double> &additionalStepAttributes)
326{
327 bunch->calcBeamParameters();
328
329 double t = bunch->getT();
330 double pathLength = bunch->get_sPos();
331 Vector_t rmin = bunch->get_origin();
332 Vector_t rmax = bunch->get_maxExtent();
333 Vector_t centroid = bunch->get_centroid();
334
335 Vector_t meanR = bunch->get_rmean();
336 Vector_t meanP = bunch->get_pmean();
337 Vector_t xsigma = bunch->get_rrms();
338 Vector_t psigma = bunch->get_prms();
339 Vector_t vareps = bunch->get_norm_emit();
340 Vector_t geomvareps = bunch->get_emit();
341
342 Vector_t RefPartR = bunch->RefPartR_m;
343 Vector_t RefPartP = bunch->RefPartP_m;
344
345 double meanEnergy = bunch->get_meanKineticEnergy();
346 double energySpread = bunch->getdE();
347 double I_0 = 4.0 * Physics::pi * Physics::epsilon_0 * Physics::c * bunch->getM() / bunch->getQ();
348 double sigma = ((xsigma[0] * xsigma[0]) + (xsigma[1] * xsigma[1])) /
349 (2.0 * bunch->get_gamma() * I_0 * (geomvareps[0] * geomvareps[0] + geomvareps[1] * geomvareps[1]));
350
351 h5_int64_t localTrackStep = (h5_int64_t)bunch->getLocalTrackStep();
352 h5_int64_t globalTrackStep = (h5_int64_t)bunch->getGlobalTrackStep();
353
354 h5_int64_t numBunch = (h5_int64_t)bunch->getNumBunch();
355 h5_int64_t SteptoLastInj = (h5_int64_t)bunch->getSteptoLastInj();
356
357 double mass = Units::eV2GeV * bunch->getM();
358 double charge = bunch->getCharge();
359
360 h5_int64_t localFrame = ( Options::psDumpFrame != DumpFrame::GLOBAL );
361
362 double sposHead = 0.0;
363 double sposRef = 0.0;
364 double sposTail = 0.0;
365
366 Vector_t maxP(0.0);
367 Vector_t minP(0.0);
368 bunch->get_PBounds(minP, maxP);
369
370 REPORTONERROR(H5SetStep(file_m, numSteps_m));
371
372 char const *OPALFlavour = "opal-cycl";
373 WRITESTRINGSTEPATTRIB(file_m, "OPAL_flavour", OPALFlavour);
374 WRITESTEPATTRIB(Float64, file_m, "TIME", &t, 1);
375 WRITESTEPATTRIB(Float64, file_m, "SPOS", &pathLength, 1);
376
377 WRITESTEPATTRIB(Float64, file_m, "RefPartR", (h5_float64_t *)&RefPartR, 3);
378 WRITESTEPATTRIB(Float64, file_m, "centroid", (h5_float64_t *)&centroid, 3);
379 WRITESTEPATTRIB(Float64, file_m, "RMSX", (h5_float64_t *)&xsigma, 3);
380 WRITESTEPATTRIB(Float64, file_m, "minX", (h5_float64_t *)&rmin, 3);
381 WRITESTEPATTRIB(Float64, file_m, "maxX", (h5_float64_t *)&rmax, 3);
382
383 WRITESTEPATTRIB(Float64, file_m, "RefPartP", (h5_float64_t *)&RefPartP, 3);
384 WRITESTEPATTRIB(Float64, file_m, "RMSP", (h5_float64_t *)&psigma, 3);
385 WRITESTEPATTRIB(Float64, file_m, "minP", (h5_float64_t *)&minP, 3);
386 WRITESTEPATTRIB(Float64, file_m, "maxP", (h5_float64_t *)&maxP, 3);
387
388 WRITESTEPATTRIB(Float64, file_m, "#varepsilon", (h5_float64_t *)&vareps, 3);
389 WRITESTEPATTRIB(Float64, file_m, "#varepsilon-geom", (h5_float64_t *)&geomvareps, 3);
390
391 WRITESTEPATTRIB(Int64, file_m, "Step", &numSteps_m, 1);
392 WRITESTEPATTRIB(Int64, file_m, "LocalTrackStep", &localTrackStep, 1);
393 WRITESTEPATTRIB(Int64, file_m, "GlobalTrackStep", &globalTrackStep, 1);
394
395 WRITESTEPATTRIB(Float64, file_m, "#sigma", &sigma, 1);
396
397 WRITESTEPATTRIB(Float64, file_m, "ENERGY", &meanEnergy, 1);
398 WRITESTEPATTRIB(Float64, file_m, "dE", &energySpread, 1);
399
401 WRITESTEPATTRIB(Float64, file_m, "MASS", &mass, 1);
402 WRITESTEPATTRIB(Float64, file_m, "CHARGE", &charge, 1);
403
404 WRITESTEPATTRIB(Int64, file_m, "NumBunch", &numBunch, 1);
405 WRITESTEPATTRIB(Int64, file_m, "SteptoLastInj", &SteptoLastInj, 1);
406 WRITESTEPATTRIB(Int64, file_m, "LOCAL", &localFrame, 1);
407
408 try {
409 h5_float64_t refpr = additionalStepAttributes.at("REFPR");
410 h5_float64_t refpt = additionalStepAttributes.at("REFPT");
411 h5_float64_t refpz = additionalStepAttributes.at("REFPZ");
412 h5_float64_t refr = additionalStepAttributes.at("REFR");
413 h5_float64_t reft = additionalStepAttributes.at("REFTHETA");
414 h5_float64_t refz = additionalStepAttributes.at("REFZ");
415 h5_float64_t azimuth = additionalStepAttributes.at("AZIMUTH");
416 h5_float64_t elevation = additionalStepAttributes.at("ELEVATION");
417
418 Vector_t referenceB(additionalStepAttributes.at("B-ref_x"),
419 additionalStepAttributes.at("B-ref_z"),
420 additionalStepAttributes.at("B-ref_y"));
421 Vector_t referenceE(additionalStepAttributes.at("E-ref_x"),
422 additionalStepAttributes.at("E-ref_z"),
423 additionalStepAttributes.at("E-ref_y"));
424 Vector_t headB(additionalStepAttributes.at("B-head_x"),
425 additionalStepAttributes.at("B-head_z"),
426 additionalStepAttributes.at("B-head_y"));
427 Vector_t headE(additionalStepAttributes.at("E-head_x"),
428 additionalStepAttributes.at("E-head_z"),
429 additionalStepAttributes.at("E-head_y"));
430 Vector_t tailB(additionalStepAttributes.at("B-tail_x"),
431 additionalStepAttributes.at("B-tail_z"),
432 additionalStepAttributes.at("B-tail_y"));
433 Vector_t tailE(additionalStepAttributes.at("E-tail_x"),
434 additionalStepAttributes.at("E-tail_z"),
435 additionalStepAttributes.at("E-tail_y"));
436
437 WRITESTEPATTRIB(Float64, file_m, "REFPR", &refpr, 1);
438 WRITESTEPATTRIB(Float64, file_m, "REFPT", &refpt, 1);
439 WRITESTEPATTRIB(Float64, file_m, "REFPZ", &refpz, 1);
440 WRITESTEPATTRIB(Float64, file_m, "REFR", &refr, 1);
441 WRITESTEPATTRIB(Float64, file_m, "REFTHETA", &reft, 1);
442 WRITESTEPATTRIB(Float64, file_m, "REFZ", &refz, 1);
443 WRITESTEPATTRIB(Float64, file_m, "AZIMUTH", &azimuth, 1);
444 WRITESTEPATTRIB(Float64, file_m, "ELEVATION", &elevation, 1);
445
446 WRITESTEPATTRIB(Float64, file_m, "spos-head", &sposHead, 1);
447 WRITESTEPATTRIB(Float64, file_m, "spos-ref", &sposRef, 1);
448 WRITESTEPATTRIB(Float64, file_m, "spos-tail", &sposTail, 1);
449
450 WRITESTEPATTRIB(Float64, file_m, "B-ref", (h5_float64_t *)&referenceB, 3);
451 WRITESTEPATTRIB(Float64, file_m, "E-ref", (h5_float64_t *)&referenceE, 3);
452 WRITESTEPATTRIB(Float64, file_m, "B-head", (h5_float64_t *)&headB, 3);
453 WRITESTEPATTRIB(Float64, file_m, "E-head", (h5_float64_t *)&headE, 3);
454 WRITESTEPATTRIB(Float64, file_m, "B-tail", (h5_float64_t *)&tailB, 3);
455 WRITESTEPATTRIB(Float64, file_m, "E-tail", (h5_float64_t *)&tailE, 3);
456 } catch (std::out_of_range & m) {
457 ERRORMSG(m.what() << endl);
458
459 throw OpalException("H5PartWrapperForPC::writeStepHeader",
460 "some additional step attribute not found");
461 }
462}
463
464void H5PartWrapperForPC::writeStepData(PartBunchBase<double, 3>* bunch) {
465 /*
466 find particle with ID==0
467 and save index in zID
468
469 FIXME After issue 287 is resolved this shouldn't be
470 necessary anymore!
471 */
472
473 size_t IDZero = bunch->getLocalNum();
474 bool found = false;
475#ifndef ENABLE_AMR
476 for(size_t k = 0; k < IDZero; ++ k) {
477 if (bunch->ID[k] == 0) {
478 found = true;
479 IDZero = k;
480 }
481 }
482#endif
483
484 const size_t numLocalParticles = (found? bunch->getLocalNum() - 1: bunch->getLocalNum());
485 const size_t skipID = IDZero;
486
487 std::vector<char> buffer(numLocalParticles * sizeof(h5_float64_t));
488 char* buffer_ptr = Util::c_data(buffer);
489 h5_float64_t *f64buffer = reinterpret_cast<h5_float64_t*>(buffer_ptr);
490 h5_int64_t *i64buffer = reinterpret_cast<h5_int64_t*> (buffer_ptr);
491
492
493 REPORTONERROR(H5PartSetNumParticles(file_m, numLocalParticles));
494
495 for(size_t i = 0; i < skipID; ++ i)
496 f64buffer[i] = bunch->R[i](0);
497 for (size_t i = skipID; i < numLocalParticles; ++ i)
498 f64buffer[i] = bunch->R[i + 1](0);
499
500 WRITEDATA(Float64, file_m, "x", f64buffer);
501
502 for(size_t i = 0; i < skipID; ++ i)
503 f64buffer[i] = bunch->R[i](1);
504 for (size_t i = skipID; i < numLocalParticles; ++ i)
505 f64buffer[i] = bunch->R[i + 1](1);
506
507 WRITEDATA(Float64, file_m, "y", f64buffer);
508
509 for(size_t i = 0; i < skipID; ++ i)
510 f64buffer[i] = bunch->R[i](2);
511 for (size_t i = skipID; i < numLocalParticles; ++ i)
512 f64buffer[i] = bunch->R[i + 1](2);
513
514 WRITEDATA(Float64, file_m, "z", f64buffer);
515
516 for(size_t i = 0; i < skipID; ++ i)
517 f64buffer[i] = bunch->P[i](0);
518 for (size_t i = skipID; i < numLocalParticles; ++ i)
519 f64buffer[i] = bunch->P[i + 1](0);
520
521 WRITEDATA(Float64, file_m, "px", f64buffer);
522
523 for(size_t i = 0; i < skipID; ++ i)
524 f64buffer[i] = bunch->P[i](1);
525 for (size_t i = skipID; i < numLocalParticles; ++ i)
526 f64buffer[i] = bunch->P[i + 1](1);
527
528 WRITEDATA(Float64, file_m, "py", f64buffer);
529
530 for(size_t i = 0; i < skipID; ++ i)
531 f64buffer[i] = bunch->P[i](2);
532 for (size_t i = skipID; i < numLocalParticles; ++ i)
533 f64buffer[i] = bunch->P[i + 1](2);
534
535 WRITEDATA(Float64, file_m, "pz", f64buffer);
536
537 for(size_t i = 0; i < skipID; ++ i)
538 f64buffer[i] = bunch->Q[i];
539 for (size_t i = skipID; i < numLocalParticles; ++ i)
540 f64buffer[i] = bunch->Q[i + 1];
541
542 WRITEDATA(Float64, file_m, "q", f64buffer);
543
544 for(size_t i = 0; i < skipID; ++ i)
545 f64buffer[i] = bunch->M[i];
546 for (size_t i = skipID; i < numLocalParticles; ++ i)
547 f64buffer[i] = bunch->M[i + 1];
548
549 WRITEDATA(Float64, file_m, "mass", f64buffer);
550
551 for(size_t i = 0; i < skipID; ++ i)
552 i64buffer[i] = bunch->ID[i];
553 for (size_t i = skipID; i < numLocalParticles; ++ i)
554 i64buffer[i] = bunch->ID[i + 1];
555
556 WRITEDATA(Int64, file_m, "id", i64buffer);
557
558 if ( bunch->hasBinning() ) {
559 for(size_t i = 0; i < skipID; ++ i)
560 i64buffer[i] = bunch->Bin[i];
561 for (size_t i = skipID; i < numLocalParticles; ++ i)
562 i64buffer[i] = bunch->Bin[i + 1];
563
564 WRITEDATA(Int64, file_m, "bin", i64buffer);
565 }
566
567 for (size_t i = 0; i < skipID; ++ i)
568 i64buffer[i] = bunch->bunchNum[i];
569 for (size_t i = skipID; i < numLocalParticles; ++ i)
570 i64buffer[i] = bunch->bunchNum[i + 1];
571 WRITEDATA(Int64, file_m, "bunchNumber", i64buffer);
572
573 if (Options::ebDump) {
574 for(size_t i = 0; i < skipID; ++ i)
575 f64buffer[i] = bunch->Ef[i](0);
576 for (size_t i = skipID; i < numLocalParticles; ++ i)
577 f64buffer[i] = bunch->Ef[i + 1](0);
578
579 WRITEDATA(Float64, file_m, "Ex", f64buffer);
580
581 for(size_t i = 0; i < skipID; ++ i)
582 f64buffer[i] = bunch->Ef[i](1);
583 for (size_t i = skipID; i < numLocalParticles; ++ i)
584 f64buffer[i] = bunch->Ef[i + 1](1);
585
586 WRITEDATA(Float64, file_m, "Ey", f64buffer);
587
588 for(size_t i = 0; i < skipID; ++ i)
589 f64buffer[i] = bunch->Ef[i](2);
590 for (size_t i = skipID; i < numLocalParticles; ++ i)
591 f64buffer[i] = bunch->Ef[i + 1](2);
592
593 WRITEDATA(Float64, file_m, "Ez", f64buffer);
594
595 for(size_t i = 0; i < skipID; ++ i)
596 f64buffer[i] = bunch->Bf[i](0);
597 for (size_t i = skipID; i < numLocalParticles; ++ i)
598 f64buffer[i] = bunch->Bf[i + 1](0);
599
600 WRITEDATA(Float64, file_m, "Bx", f64buffer);
601
602 for(size_t i = 0; i < skipID; ++ i)
603 f64buffer[i] = bunch->Bf[i](1);
604 for (size_t i = skipID; i < numLocalParticles; ++ i)
605 f64buffer[i] = bunch->Bf[i + 1](1);
606
607 WRITEDATA(Float64, file_m, "By", f64buffer);
608
609 for(size_t i = 0; i < skipID; ++ i)
610 f64buffer[i] = bunch->Bf[i](2);
611 for (size_t i = skipID; i < numLocalParticles; ++ i)
612 f64buffer[i] = bunch->Bf[i + 1](2);
613
614 WRITEDATA(Float64, file_m, "Bz", f64buffer);
615
616 }
617
619 if(Options::rhoDump) {
620 NDIndex<3> idx = bunch->getFieldLayout().getLocalNDIndex();
621 NDIndex<3> elem;
622 REPORTONERROR(H5Block3dSetView(file_m,
623 idx[0].min(), idx[0].max(),
624 idx[1].min(), idx[1].max(),
625 idx[2].min(), idx[2].max()));
626
627 std::unique_ptr<h5_float64_t[]> data(new h5_float64_t[(idx[0].max() + 1) * (idx[1].max() + 1) * (idx[2].max() + 1)]);
628
629 int ii = 0;
630 // h5block uses the fortran convention of storing data:
631 // INTEGER, DIMENSION(2,3) :: a
632 // => {a(1,1), a(2,1), a(1,2), a(2,2), a(1,3), a(2,3)}
633 for(int i = idx[2].min(); i <= idx[2].max(); ++ i) {
634 for(int j = idx[1].min(); j <= idx[1].max(); ++ j) {
635 for(int k = idx[0].min(); k <= idx[0].max(); ++ k) {
636 data[ii] = bunch->getRho(k, j, i);
637 ++ ii;
638 }
639 }
640 }
641 REPORTONERROR(H5Block3dWriteScalarFieldFloat64(file_m, "rho", data.get()));
642
644 REPORTONERROR(H5Block3dSetFieldOrigin(file_m, "rho",
645 (h5_float64_t)bunch->get_origin()(0),
646 (h5_float64_t)bunch->get_origin()(1),
647 (h5_float64_t)bunch->get_origin()(2)));
648
649 REPORTONERROR(H5Block3dSetFieldSpacing(file_m, "rho",
650 (h5_float64_t)bunch->get_hr()(0),
651 (h5_float64_t)bunch->get_hr()(1),
652 (h5_float64_t)bunch->get_hr()(2)));
653 }
654}
OPAL_PROJECT_VERSION
#define OPAL_PROJECT_VERSION
Definition: OPALconfig.h:5
OPAL_PROJECT_NAME
#define OPAL_PROJECT_NAME
Definition: OPALconfig.h:2
pow
Tps< T > pow(const Tps< T > &x, int y)
Integer power.
Definition: TpsMath.h:76
sqrt
Tps< T > sqrt(const Tps< T > &x)
Square root.
Definition: TpsMath.h:91
WRITESTRINGFILEATTRIB
#define WRITESTRINGFILEATTRIB(file, name, value)
Definition: H5PartWrapper.h:24
WRITESTEPATTRIB
#define WRITESTEPATTRIB(type, file, name, value, length)
Definition: H5PartWrapper.h:29
REPORTONERROR
#define REPORTONERROR(rc)
Definition: H5PartWrapper.h:22
READSTEPATTRIB
#define READSTEPATTRIB(type, file, name, value)
Definition: H5PartWrapper.h:27
WRITEDATA
#define WRITEDATA(type, file, name, value)
Definition: H5PartWrapper.h:32
WRITEFILEATTRIB
#define WRITEFILEATTRIB(type, file, name, value, length)
Definition: H5PartWrapper.h:25
WRITESTRINGSTEPATTRIB
#define WRITESTRINGSTEPATTRIB(file, name, value)
Definition: H5PartWrapper.h:28
READDATA
#define READDATA(type, file, name, value)
Definition: H5PartWrapper.h:31
READFILEATTRIB
#define READFILEATTRIB(type, file, name, value)
Definition: H5PartWrapper.h:23
max
T::PETE_Expr_t::PETE_Return_t max(const PETE_Expr< T > &expr, NDIndex< D > &loc)
Definition: ReductionLoc.h:84
min
T::PETE_Expr_t::PETE_Return_t min(const PETE_Expr< T > &expr, NDIndex< D > &loc)
Definition: ReductionLoc.h:76
endl
Inform & endl(Inform &inf)
Definition: Inform.cpp:42
ERRORMSG
#define ERRORMSG(msg)
Definition: IpplInfo.h:350
PyOpal::PyElementNS::attributeName
std::map< AttributeType, std::string > attributeName
Definition: PyElement.cpp:7
Physics::epsilon_0
constexpr double epsilon_0
The permittivity of vacuum in As/Vm.
Definition: Physics.h:51
Physics::c
constexpr double c
The velocity of light in m/s.
Definition: Physics.h:45
Physics::pi
constexpr double pi
The value of.
Definition: Physics.h:30
Units::eV2GeV
constexpr double eV2GeV
Definition: Units.h:71
Units::eV2MeV
constexpr double eV2MeV
Definition: Units.h:77
Options::psDumpFreq
int psDumpFreq
The frequency to dump the phase space, i.e.dump data when steppsDumpFreq==0.
Definition: Options.cpp:39
Options::rhoDump
bool rhoDump
Definition: Options.cpp:63
Options::ebDump
bool ebDump
Definition: Options.cpp:65
Options::psDumpFrame
DumpFrame psDumpFrame
flag to decide in which coordinate frame the phase space will be dumped for OPAL-cycl
Definition: Options.cpp:45
Util::c_data
T * c_data(std::vector< T, A > &v)
Definition: Util.h:227
Util::getGitRevision
std::string getGitRevision()
Definition: Util.cpp:32
PartBunchBase::Ef
ParticleAttrib< Vector_t > Ef
Definition: PartBunchBase.h:578
PartBunchBase::get_meanKineticEnergy
double get_meanKineticEnergy() const
Definition: PartBunchBase.hpp:1067
PartBunchBase::R
ParticlePos_t & R
Definition: PartBunchBase.h:570
PartBunchBase::Bin
ParticleAttrib< int > Bin
Definition: PartBunchBase.h:582
PartBunchBase::RefPartP_m
Vector_t RefPartP_m
Definition: PartBunchBase.h:591
PartBunchBase::getSteptoLastInj
int getSteptoLastInj() const
Definition: PartBunchBase.hpp:1512
PartBunchBase::setNumBunch
void setNumBunch(short n)
Definition: PartBunchBase.hpp:1413
PartBunchBase::getQ
double getQ() const
Access to reference data.
Definition: PartBunchBase.hpp:1630
PartBunchBase::getLocalNum
size_t getLocalNum() const
Definition: PartBunchBase.hpp:1905
PartBunchBase::set_sPos
void set_sPos(double s)
Definition: PartBunchBase.hpp:1055
PartBunchBase::setLocalTrackStep
void setLocalTrackStep(long long n)
step in a TRACK command
Definition: PartBunchBase.hpp:1395
PartBunchBase::M
ParticleAttrib< double > M
Definition: PartBunchBase.h:576
PartBunchBase::getTotalNum
size_t getTotalNum() const
Definition: PartBunchBase.hpp:1900
PartBunchBase::get_gamma
double get_gamma() const
Definition: PartBunchBase.hpp:1061
PartBunchBase::P
ParticleAttrib< Vector_t > P
Definition: PartBunchBase.h:574
PartBunchBase::get_rrms
Vector_t get_rrms() const
Definition: PartBunchBase.hpp:1111
PartBunchBase::RefPartR_m
Vector_t RefPartR_m
Definition: PartBunchBase.h:590
PartBunchBase::Q
ParticleAttrib< double > Q
Definition: PartBunchBase.h:575
PartBunchBase::calcBeamParameters
void calcBeamParameters()
Definition: PartBunchBase.hpp:1267
PartBunchBase::setSteptoLastInj
void setSteptoLastInj(int n)
Definition: PartBunchBase.hpp:1506
PartBunchBase::getRho
virtual double getRho(int x, int y, int z)=0
PartBunchBase::get_origin
Vector_t get_origin() const
Definition: PartBunchBase.hpp:1093
PartBunchBase::get_prms
Vector_t get_prms() const
Definition: PartBunchBase.hpp:1129
PartBunchBase::getdE
double getdE() const
Definition: PartBunchBase.hpp:1694
PartBunchBase::getCharge
double getCharge() const
get the total charge per simulation particle
Definition: PartBunchBase.hpp:1317
PartBunchBase::getLocalTrackStep
long long getLocalTrackStep() const
Definition: PartBunchBase.hpp:1407
PartBunchBase::getNumBunch
short getNumBunch() const
Definition: PartBunchBase.hpp:1421
PartBunchBase::get_norm_emit
Vector_t get_norm_emit() const
Definition: PartBunchBase.hpp:1157
PartBunchBase::get_maxExtent
Vector_t get_maxExtent() const
Definition: PartBunchBase.hpp:1099
PartBunchBase::get_PBounds
void get_PBounds(Vector_t &min, Vector_t &max) const
Definition: PartBunchBase.hpp:1261
PartBunchBase::setGlobalTrackStep
void setGlobalTrackStep(long long n)
step in multiple TRACK commands
Definition: PartBunchBase.hpp:1383
PartBunchBase::get_centroid
Vector_t get_centroid() const
Definition: PartBunchBase.hpp:1105
PartBunchBase::bunchNum
ParticleAttrib< short > bunchNum
Definition: PartBunchBase.h:588
PartBunchBase::get_pmean
Vector_t get_pmean() const
Definition: PartBunchBase.hpp:1135
PartBunchBase::hasBinning
bool hasBinning() const
Definition: PartBunchBase.hpp:248
PartBunchBase::Bf
ParticleAttrib< Vector_t > Bf
Definition: PartBunchBase.h:581
PartBunchBase::getGlobalTrackStep
long long getGlobalTrackStep() const
Definition: PartBunchBase.hpp:1389
PartBunchBase::setT
void setT(double t)
Definition: PartBunchBase.hpp:1031
PartBunchBase::get_hr
virtual Vector_t get_hr() const
Definition: PartBunchBase.hpp:1232
PartBunchBase::ID
ParticleIndex_t & ID
Definition: PartBunchBase.h:571
PartBunchBase::get_emit
Vector_t get_emit() const
Definition: PartBunchBase.hpp:1151
PartBunchBase::getFieldLayout
virtual FieldLayout_t & getFieldLayout()=0
PartBunchBase::get_sPos
double get_sPos() const
Definition: PartBunchBase.hpp:1049
PartBunchBase::getM
double getM() const
Definition: PartBunchBase.hpp:1636
PartBunchBase::getT
double getT() const
Definition: PartBunchBase.hpp:1043
PartBunchBase::get_rmean
Vector_t get_rmean() const
Definition: PartBunchBase.hpp:1123
OpalData::setRestartDumpFreq
void setRestartDumpFreq(const int &N)
set the dump frequency as found in restart file
Definition: OpalData.cpp:341
OpalData::getInstance
static OpalData * getInstance()
Definition: OpalData.cpp:196
H5PartWrapper::getNumParticles
size_t getNumParticles() const
Definition: H5PartWrapper.cpp:525
H5PartWrapper::open
void open(h5_int32_t flags)
Definition: H5PartWrapper.cpp:72
H5PartWrapper::numSteps_m
h5_int64_t numSteps_m
Definition: H5PartWrapper.h:76
H5PartWrapper::predecessorOPALFlavour_m
std::string predecessorOPALFlavour_m
Definition: H5PartWrapper.h:75
H5PartWrapper::file_m
h5_file_t file_m
Definition: H5PartWrapper.h:73
H5PartWrapperForPC::writeStep
virtual void writeStep(PartBunchBase< double, 3 > *, const std::map< std::string, double > &additionalStepAttributes)
Definition: H5PartWrapperForPC.cpp:311
H5PartWrapperForPC::readStepData
void readStepData(PartBunchBase< double, 3 > *, h5_ssize_t, h5_ssize_t)
Definition: H5PartWrapperForPC.cpp:148
H5PartWrapperForPC::elevation_m
h5_float64_t elevation_m
Definition: H5PartWrapperForPC.h:52
H5PartWrapperForPC::readStepHeader
void readStepHeader(PartBunchBase< double, 3 > *)
Definition: H5PartWrapperForPC.cpp:85
H5PartWrapperForPC::writeHeader
virtual void writeHeader()
Definition: H5PartWrapperForPC.cpp:225
H5PartWrapperForPC::H5PartWrapperForPC
H5PartWrapperForPC(const std::string &fileName, h5_int32_t flags=H5_O_WRONLY)
Definition: H5PartWrapperForPC.cpp:18
H5PartWrapperForPC::writeStepData
void writeStepData(PartBunchBase< double, 3 > *)
Definition: H5PartWrapperForPC.cpp:464
H5PartWrapperForPC::meanMomentum_m
h5_float64_t meanMomentum_m
Definition: H5PartWrapperForPC.h:50
H5PartWrapperForPC::readHeader
virtual void readHeader()
Definition: H5PartWrapperForPC.cpp:46
H5PartWrapperForPC::writeStepHeader
void writeStepHeader(PartBunchBase< double, 3 > *, const std::map< std::string, double > &)
Definition: H5PartWrapperForPC.cpp:324
H5PartWrapperForPC::readStep
virtual void readStep(PartBunchBase< double, 3 > *, h5_ssize_t firstParticle, h5_ssize_t lastParticle)
Definition: H5PartWrapperForPC.cpp:73
OpalException
The base class for all OPAL exceptions.
Definition: OpalException.h:28
FieldLayout::getLocalNDIndex
NDIndex< Dim > getLocalNDIndex()
Definition: FieldLayout.hpp:1579
Generated on Thu Oct 20 2022 17:40:12 for OPAL (Object Oriented Parallel Accelerator Library) by doxygen 1.9.3