OPAL (Object Oriented Parallel Accelerator Library)  2021.1.99
OPAL
H5PartWrapperForPT.cpp
Go to the documentation of this file.
1 //
2 // Class H5PartWrapperForPT
3 // A class that manages all calls to H5Part for the Parallel-T tracker.
4 //
5 // Copyright (c) 200x-2021, Paul Scherrer Institut, Villigen PSI, Switzerland
6 // All rights reserved
7 //
8 // This file is part of OPAL.
9 //
10 // OPAL is free software: you can redistribute it and/or modify
11 // it under the terms of the GNU General Public License as published by
12 // the Free Software Foundation, either version 3 of the License, or
13 // (at your option) any later version.
14 //
15 // You should have received a copy of the GNU General Public License
16 // along with OPAL. If not, see <https://www.gnu.org/licenses/>.
17 //
19 
20 #include "OPALconfig.h"
23 #include "Algorithms/Vektor.h"
24 #include "Utilities/Options.h"
25 #include "Utilities/Util.h"
26 #include "Physics/Physics.h"
27 
28 #include "h5core/h5_types.h"
29 
30 #include <cmath>
31 #include <sstream>
32 #include <set>
33 
34 H5PartWrapperForPT::H5PartWrapperForPT(const std::string& fileName, h5_int32_t flags):
35  H5PartWrapper(fileName, flags)
36 { }
37 
38 H5PartWrapperForPT::H5PartWrapperForPT(const std::string& fileName, int restartStep, std::string sourceFile, h5_int32_t flags):
39  H5PartWrapper(fileName, restartStep, sourceFile, flags)
40 {
41  if (restartStep == -1) {
42  restartStep = H5GetNumSteps(file_m) - 1 ;
43  OpalData::getInstance()->setRestartStep(restartStep);
44  }
45 }
46 
48 { }
49 
51  h5_int64_t numFileAttributes = H5GetNumFileAttribs(file_m);
52 
53  const h5_size_t lengthAttributeName = 256;
54  char attributeName[lengthAttributeName];
55  h5_int64_t attributeType;
56  h5_size_t numAttributeElements;
57  std::set<std::string> attributeNames;
58 
59  for (h5_int64_t i = 0; i < numFileAttributes; ++ i) {
60  REPORTONERROR(H5GetFileAttribInfo(file_m,
61  i,
62  attributeName,
63  lengthAttributeName,
64  &attributeType,
65  &numAttributeElements));
66 
67  attributeNames.insert(attributeName);
68  }
69 
70  if (attributeNames.find("dump frequency") != attributeNames.end()) {
71  h5_int64_t dumpfreq;
72  READFILEATTRIB(Int64, file_m, "dump frequency", &dumpfreq);
74  }
75 
76  if (attributeNames.find("dPhiGlobal") != attributeNames.end()) {
77  h5_float64_t dphi;
78  READFILEATTRIB(Float64, file_m, "dPhiGlobal", &dphi);
80  }
81 
82  if (attributeNames.find("nAutoPhaseCavities") != attributeNames.end()) {
83  auto opal = OpalData::getInstance();
84  h5_float64_t phase;
85  h5_int64_t numAutoPhaseCavities = 0;
86  if (!H5HasFileAttrib(file_m, "nAutoPhaseCavities") ||
87  H5ReadFileAttribInt64(file_m, "nAutoPhaseCavities", &numAutoPhaseCavities) != H5_SUCCESS) {
88  numAutoPhaseCavities = 0;
89  } else {
90  for (long i = 0; i < numAutoPhaseCavities; ++ i) {
91  std::string elementName = "Cav-" + std::to_string(i + 1) + "-name";
92  std::string elementPhase = "Cav-" + std::to_string(i + 1) + "-value";
93 
94  char name[128];
95  READFILEATTRIB(String, file_m, elementName.c_str(), name);
96  READFILEATTRIB(Float64, file_m, elementPhase.c_str(), &phase);
97 
98  opal->setMaxPhase(name, phase);
99  }
100  }
101  }
102 }
103 
104 void H5PartWrapperForPT::readStep(PartBunchBase<double, 3>* bunch, h5_ssize_t firstParticle, h5_ssize_t lastParticle) {
105  h5_ssize_t numStepsInSource = H5GetNumSteps(file_m);
106  h5_ssize_t readStep = numStepsInSource - 1;
107  REPORTONERROR(H5SetStep(file_m, readStep));
108 
109  readStepHeader(bunch);
110  readStepData(bunch, firstParticle, lastParticle);
111 }
112 
114  double actualT;
115  READSTEPATTRIB(Float64, file_m, "TIME", &actualT);
116  bunch->setT(actualT);
117 
118  double spos;
119  READSTEPATTRIB(Float64, file_m, "SPOS", &spos);
120  bunch->set_sPos(spos);
121 
122  h5_int64_t ltstep;
123  READSTEPATTRIB(Int64, file_m, "LocalTrackStep", &ltstep);
124  bunch->setLocalTrackStep((long long)(ltstep + 1));
125 
126  h5_int64_t gtstep;
127  READSTEPATTRIB(Int64, file_m, "GlobalTrackStep", &gtstep);
128  bunch->setGlobalTrackStep((long long)(gtstep + 1));
129 
130  Vector_t RefPartR;
131  READSTEPATTRIB(Float64, file_m, "RefPartR", (h5_float64_t *)&RefPartR);
132  bunch->RefPartR_m = RefPartR;
133 
134  Vector_t RefPartP;
135  READSTEPATTRIB(Float64, file_m, "RefPartP", (h5_float64_t *)&RefPartP);
136  bunch->RefPartP_m = RefPartP;
137 
138  Vector_t TaitBryant;
139  READSTEPATTRIB(Float64, file_m, "TaitBryantAngles", (h5_float64_t *)&TaitBryant);
140  Quaternion rotTheta(std::cos(0.5 * TaitBryant[0]), 0, std::sin(0.5 * TaitBryant[0]), 0);
141  Quaternion rotPhi(std::cos(0.5 * TaitBryant[1]), std::sin(0.5 * TaitBryant[1]), 0, 0);
142  Quaternion rotPsi(std::cos(0.5 * TaitBryant[2]), 0, 0, std::sin(0.5 * TaitBryant[2]));
143  Quaternion rotation = rotTheta * (rotPhi * rotPsi);
144  bunch->toLabTrafo_m = CoordinateSystemTrafo(-rotation.conjugate().rotate(RefPartR), rotation);
145 }
146 
147 void H5PartWrapperForPT::readStepData(PartBunchBase<double, 3>* bunch, h5_ssize_t firstParticle, h5_ssize_t lastParticle) {
148  h5_ssize_t numParticles = getNumParticles();
149  if (lastParticle >= numParticles || firstParticle > lastParticle) {
150  throw OpalException("H5PartWrapperForPT::readStepData",
151  "the provided particle numbers don't match the number of particles in the file");
152  }
153 
154  REPORTONERROR(H5PartSetView(file_m, firstParticle, lastParticle));
155 
156  numParticles = lastParticle - firstParticle + 1;
157 
158  std::vector<char> buffer(numParticles * sizeof(h5_float64_t));
159  char* buffer_ptr = Util::c_data(buffer);
160  h5_float64_t* f64buffer = reinterpret_cast<h5_float64_t*>(buffer_ptr);
161  h5_int32_t* i32buffer = reinterpret_cast<h5_int32_t*>(buffer_ptr);
162 
163  READDATA(Float64, file_m, "x", f64buffer);
164  for (long int n = 0; n < numParticles; ++ n) {
165  bunch->R[n](0) = f64buffer[n];
166  bunch->Bin[n] = 0;
167  }
168 
169  READDATA(Float64, file_m, "y", f64buffer);
170  for (long int n = 0; n < numParticles; ++ n) {
171  bunch->R[n](1) = f64buffer[n];
172  }
173 
174  READDATA(Float64, file_m, "z", f64buffer);
175  for (long int n = 0; n < numParticles; ++ n) {
176  bunch->R[n](2) = f64buffer[n];
177  }
178 
179  READDATA(Float64, file_m, "px", f64buffer);
180  for (long int n = 0; n < numParticles; ++ n) {
181  bunch->P[n](0) = f64buffer[n];
182  }
183 
184  READDATA(Float64, file_m, "py", f64buffer);
185  for (long int n = 0; n < numParticles; ++ n) {
186  bunch->P[n](1) = f64buffer[n];
187  }
188 
189  READDATA(Float64, file_m, "pz", f64buffer);
190  for (long int n = 0; n < numParticles; ++ n) {
191  bunch->P[n](2) = f64buffer[n];
192  }
193 
194  READDATA(Float64, file_m, "q", f64buffer);
195  for (long int n = 0; n < numParticles; ++ n) {
196  bunch->Q[n] = f64buffer[n];
197  }
198 
199  READDATA(Int32, file_m, "id", i32buffer);
200  for (long int n = 0; n < numParticles; ++ n) {
201  bunch->ID[n] = i32buffer[n];
202  }
203 
204  REPORTONERROR(H5PartSetView(file_m, -1, -1));
205 }
206 
208  std::stringstream OPAL_version;
209  OPAL_version << OPAL_PROJECT_NAME << " " << OPAL_PROJECT_VERSION << " # git rev. " << Util::getGitRevision();
210  WRITESTRINGFILEATTRIB(file_m, "OPAL_version", OPAL_version.str().c_str());
211 
212  WRITESTRINGFILEATTRIB(file_m, "idUnit", "1");
213  WRITESTRINGFILEATTRIB(file_m, "xUnit", "m");
214  WRITESTRINGFILEATTRIB(file_m, "yUnit", "m");
215  WRITESTRINGFILEATTRIB(file_m, "zUnit", "m");
216 
217  WRITESTRINGFILEATTRIB(file_m, "pxUnit", "#beta#gamma");
218  WRITESTRINGFILEATTRIB(file_m, "pyUnit", "#beta#gamma");
219  WRITESTRINGFILEATTRIB(file_m, "pzUnit", "#beta#gamma");
220 
221  WRITESTRINGFILEATTRIB(file_m, "ptypeUnit", "1");
222  WRITESTRINGFILEATTRIB(file_m, "poriginUnit", "1");
223  WRITESTRINGFILEATTRIB(file_m, "qUnit", "C");
224 
225  if (Options::ebDump) {
226  WRITESTRINGFILEATTRIB(file_m, "ExUnit", "MV/m");
227  WRITESTRINGFILEATTRIB(file_m, "EyUnit", "MV/m");
228  WRITESTRINGFILEATTRIB(file_m, "EzUnit", "MV/m");
229 
230  WRITESTRINGFILEATTRIB(file_m, "BxUnit", "T");
231  WRITESTRINGFILEATTRIB(file_m, "ByUnit", "T");
232  WRITESTRINGFILEATTRIB(file_m, "BzUnit", "T");
233  }
234 
235  if (Options::rhoDump) {
236  WRITESTRINGFILEATTRIB(file_m, "rhoUnit", "C/m^3");
237  }
238 
239  WRITESTRINGFILEATTRIB(file_m, "TaitBryantAnglesUnit", "rad");
240 
241  WRITESTRINGFILEATTRIB(file_m, "SPOSUnit", "m");
242  WRITESTRINGFILEATTRIB(file_m, "TIMEUnit", "s");
243  WRITESTRINGFILEATTRIB(file_m, "#gammaUnit", "1");
244  WRITESTRINGFILEATTRIB(file_m, "ENERGYUnit", "MeV");
245  WRITESTRINGFILEATTRIB(file_m, "#varepsilonUnit", "m rad");
246  WRITESTRINGFILEATTRIB(file_m, "#varepsilon-geomUnit", "m rad");
247 
248  WRITESTRINGFILEATTRIB(file_m, "#sigmaUnit", "1");
249  WRITESTRINGFILEATTRIB(file_m, "RMSXUnit", "m");
250  WRITESTRINGFILEATTRIB(file_m, "centroidUnit", "m");
251  WRITESTRINGFILEATTRIB(file_m, "RMSPUnit", "#beta#gamma");
252  WRITESTRINGFILEATTRIB(file_m, "MEANPUnit", "#beta#gamma");
253 
254  WRITESTRINGFILEATTRIB(file_m, "minXUnit", "m");
255  WRITESTRINGFILEATTRIB(file_m, "maxXUnit", "m");
256  WRITESTRINGFILEATTRIB(file_m, "minPUnit", "#beta#gamma");
257  WRITESTRINGFILEATTRIB(file_m, "maxPUnit", "#beta#gamma");
258 
259  WRITESTRINGFILEATTRIB(file_m, "dEUnit", "MeV");
260 
261  WRITESTRINGFILEATTRIB(file_m, "MASSUnit", "GeV");
262  WRITESTRINGFILEATTRIB(file_m, "CHARGEUnit", "C");
263 
264  WRITESTRINGFILEATTRIB(file_m, "E-refUnit", "MV/m");
265  WRITESTRINGFILEATTRIB(file_m, "B-refUnit", "T");
266 
267  WRITESTRINGFILEATTRIB(file_m, "StepUnit", "1");
268  WRITESTRINGFILEATTRIB(file_m, "LocalTrackStepUnit", "1");
269  WRITESTRINGFILEATTRIB(file_m, "GlobalTrackStepUnit", "1");
270  WRITESTRINGFILEATTRIB(file_m, "NumBunchUnit", "1");
271  WRITESTRINGFILEATTRIB(file_m, "NumPartUnit", "1");
272 
273  WRITESTRINGFILEATTRIB(file_m, "RefPartRUnit", "m");
274  WRITESTRINGFILEATTRIB(file_m, "RefPartPUnit", "#beta#gamma");
275  WRITESTRINGFILEATTRIB(file_m, "SteptoLastInjUnit", "1");
276 
277  WRITESTRINGFILEATTRIB(file_m, "dump frequencyUnit", "1")
278  WRITESTRINGFILEATTRIB(file_m, "dPhiGlobalUnit", "rad")
279 
280 
281  h5_int64_t dumpfreq = Options::psDumpFreq;
282  WRITEFILEATTRIB(Int64, file_m, "dump frequency", &dumpfreq, 1);
283 
284 
286  h5_float64_t dphi = OpalData::getInstance()->getGlobalPhaseShift();
287  WRITEFILEATTRIB(Float64, file_m, "dPhiGlobal", &dphi, 1);
288 }
289 
290 void H5PartWrapperForPT::writeStep(PartBunchBase<double, 3>* bunch, const std::map<std::string, double>& additionalStepAttributes) {
291  if (bunch->getTotalNum() == 0) return;
292 
293  open(H5_O_APPENDONLY);
294  bunch->calcBeamParameters();
295  writeStepHeader(bunch, additionalStepAttributes);
296  writeStepData(bunch);
297  close();
298 }
299 
300 void H5PartWrapperForPT::writeStepHeader(PartBunchBase<double, 3>* bunch, const std::map<std::string, double>& additionalStepAttributes) {
301  double actPos = bunch->get_sPos();
302  double t = bunch->getT();
303  Vector_t rmin = bunch->get_origin();
304  Vector_t rmax = bunch->get_maxExtent();
305  Vector_t centroid = bunch->get_centroid();
306 
307  Vector_t maxP(0.0);
308  Vector_t minP(0.0);
309 
310  Vector_t xsigma = bunch->get_rrms();
311  Vector_t psigma = bunch->get_prms();
312  Vector_t vareps = bunch->get_norm_emit();
313  Vector_t geomvareps = bunch->get_emit();
314  Vector_t RefPartR = bunch->RefPartR_m;
315  Vector_t RefPartP = bunch->RefPartP_m;
317  Vector_t pmean = bunch->get_pmean();
318 
319  double meanEnergy = bunch->get_meanKineticEnergy();
320  double energySpread = bunch->getdE();
321  double I_0 = 4.0 * Physics::pi * Physics::epsilon_0 * Physics::c * bunch->getM() / bunch->getQ();
322  double sigma = ((xsigma[0] * xsigma[0]) + (xsigma[1] * xsigma[1])) /
323  (2.0 * bunch->get_gamma() * I_0 * (geomvareps[0] * geomvareps[0] + geomvareps[1] * geomvareps[1]));
324 
325  h5_int64_t localTrackStep = (h5_int64_t)bunch->getLocalTrackStep();
326  h5_int64_t globalTrackStep = (h5_int64_t)bunch->getGlobalTrackStep();
327 
328  double mass = 1.0e-9 * bunch->getM();
329  double charge = bunch->getCharge();
330 
331  h5_int64_t numBunch = 1;
332  h5_int64_t SteptoLastInj = 0;
333 
334  bunch->get_PBounds(minP, maxP);
335 
336  /* ------------------------------------------------------------------------ */
337 
338  REPORTONERROR(H5SetStep(file_m, numSteps_m));
339 
340  char const* OPALFlavour = "opal-t";
341  WRITESTRINGSTEPATTRIB(file_m, "OPAL_flavour", OPALFlavour);
342  WRITESTEPATTRIB(Float64, file_m, "SPOS", &actPos, 1);
343  WRITESTEPATTRIB(Float64, file_m, "RefPartR", (h5_float64_t*)&RefPartR, 3);
344  WRITESTEPATTRIB(Float64, file_m, "centroid", (h5_float64_t*)&centroid, 3);
345  WRITESTEPATTRIB(Float64, file_m, "RMSX", (h5_float64_t*)&xsigma, 3);
346 
347  WRITESTEPATTRIB(Float64, file_m, "RefPartP", (h5_float64_t*)&RefPartP, 3);
348  WRITESTEPATTRIB(Float64, file_m, "MEANP", (h5_float64_t*)&pmean, 3);
349  WRITESTEPATTRIB(Float64, file_m, "RMSP", (h5_float64_t*)&psigma, 3);
350  WRITESTEPATTRIB(Float64, file_m, "TaitBryantAngles", (h5_float64_t*)&TaitBryant, 3);
351 
352  WRITESTEPATTRIB(Float64, file_m, "#varepsilon", (h5_float64_t*)&vareps, 3);
353  WRITESTEPATTRIB(Float64, file_m, "#varepsilon-geom", (h5_float64_t*)&geomvareps, 3);
354 
355  WRITESTEPATTRIB(Float64, file_m, "minX", (h5_float64_t*)&rmin, 3);
356  WRITESTEPATTRIB(Float64, file_m, "maxX", (h5_float64_t*)&rmax, 3);
357 
358  WRITESTEPATTRIB(Float64, file_m, "minP", (h5_float64_t*)&minP, 3);
359  WRITESTEPATTRIB(Float64, file_m, "maxP", (h5_float64_t*)&maxP, 3);
360 
361  WRITESTEPATTRIB(Int64, file_m, "Step", &numSteps_m, 1);
362  WRITESTEPATTRIB(Int64, file_m, "LocalTrackStep", &localTrackStep, 1);
363  WRITESTEPATTRIB(Int64, file_m, "GlobalTrackStep", &globalTrackStep, 1);
364 
365  WRITESTEPATTRIB(Float64, file_m, "#sigma", &sigma, 1);
366 
367  WRITESTEPATTRIB(Float64, file_m, "TIME", &t, 1);
368 
369  WRITESTEPATTRIB(Float64, file_m, "ENERGY", &meanEnergy, 1);
370  WRITESTEPATTRIB(Float64, file_m, "dE", &energySpread, 1);
371 
373  WRITESTEPATTRIB(Float64, file_m, "MASS", &mass, 1);
374 
375  WRITESTEPATTRIB(Float64, file_m, "CHARGE", &charge, 1);
376 
377  WRITESTEPATTRIB(Int64, file_m, "NumBunch", &numBunch, 1);
378 
379  WRITESTEPATTRIB(Int64, file_m, "SteptoLastInj", &SteptoLastInj, 1);
380 
381  try {
382  Vector_t referenceB(additionalStepAttributes.at("B-ref_x"),
383  additionalStepAttributes.at("B-ref_z"),
384  additionalStepAttributes.at("B-ref_y"));
385  Vector_t referenceE(additionalStepAttributes.at("E-ref_x"),
386  additionalStepAttributes.at("E-ref_z"),
387  additionalStepAttributes.at("E-ref_y"));
388 
389  WRITESTEPATTRIB(Float64, file_m, "B-ref", (h5_float64_t*)&referenceB, 3);
390  WRITESTEPATTRIB(Float64, file_m, "E-ref", (h5_float64_t*)&referenceE, 3);
391  } catch (std::out_of_range & m) {
392  ERRORMSG(m.what() << endl);
393 
394  throw OpalException("H5PartWrapperForPC::writeStepHeader",
395  "some additional step attribute not found");
396  }
397 
398  ++ numSteps_m;
399 }
400 
402  size_t numLocalParticles = bunch->getLocalNum();
403 
404  REPORTONERROR(H5PartSetNumParticles(file_m, numLocalParticles));
405 
406  std::vector<char> buffer(numLocalParticles * sizeof(h5_float64_t));
407  char* buffer_ptr = Util::c_data(buffer);
408  h5_float64_t* f64buffer = reinterpret_cast<h5_float64_t*>(buffer_ptr);
409  h5_int64_t* i64buffer = reinterpret_cast<h5_int64_t*>(buffer_ptr);
410  h5_int32_t* i32buffer = reinterpret_cast<h5_int32_t*>(buffer_ptr);
411 
412  for (size_t i = 0; i < numLocalParticles; ++ i)
413  f64buffer[i] = bunch->R[i](0);
414  WRITEDATA(Float64, file_m, "x", f64buffer);
415 
416  for (size_t i = 0; i < numLocalParticles; ++ i)
417  f64buffer[i] = bunch->R[i](1);
418  WRITEDATA(Float64, file_m, "y", f64buffer);
419 
420  for (size_t i = 0; i < numLocalParticles; ++ i)
421  f64buffer[i] = bunch->R[i](2);
422  WRITEDATA(Float64, file_m, "z", f64buffer);
423 
424  for (size_t i = 0; i < numLocalParticles; ++ i)
425  f64buffer[i] = bunch->P[i](0);
426  WRITEDATA(Float64, file_m, "px", f64buffer);
427 
428  for (size_t i = 0; i < numLocalParticles; ++ i)
429  f64buffer[i] = bunch->P[i](1);
430  WRITEDATA(Float64, file_m, "py", f64buffer);
431 
432  for (size_t i = 0; i < numLocalParticles; ++ i)
433  f64buffer[i] = bunch->P[i](2);
434  WRITEDATA(Float64, file_m, "pz", f64buffer);
435 
436  for (size_t i = 0; i < numLocalParticles; ++ i)
437  f64buffer[i] = bunch->Q[i];
438  WRITEDATA(Float64, file_m, "q", f64buffer);
439 
440  for (size_t i = 0; i < numLocalParticles; ++ i)
441  i64buffer[i] = bunch->ID[i];
442  WRITEDATA(Int64, file_m, "id", i64buffer);
443 
444  for (size_t i = 0; i < numLocalParticles; ++ i)
445  i32buffer[i] = (h5_int32_t) bunch->PType[i];
446  WRITEDATA(Int32, file_m, "ptype", i32buffer);
447 
448  for (size_t i = 0; i < numLocalParticles; ++ i)
449  i32buffer[i] = (h5_int32_t) bunch->POrigin[i];
450  WRITEDATA(Int32, file_m, "porigin", i32buffer);
451 
452  if (Options::ebDump) {
453  for (size_t i = 0; i < numLocalParticles; ++ i)
454  f64buffer[i] = bunch->Ef[i](0);
455  WRITEDATA(Float64, file_m, "Ex", f64buffer);
456 
457  for (size_t i = 0; i < numLocalParticles; ++ i)
458  f64buffer[i] = bunch->Ef[i](1);
459  WRITEDATA(Float64, file_m, "Ey", f64buffer);
460 
461  for (size_t i = 0; i < numLocalParticles; ++ i)
462  f64buffer[i] = bunch->Ef[i](2);
463  WRITEDATA(Float64, file_m, "Ez", f64buffer);
464 
465  for (size_t i = 0; i < numLocalParticles; ++ i)
466  f64buffer[i] = bunch->Bf[i](0);
467  WRITEDATA(Float64, file_m, "Bx", f64buffer);
468 
469  for (size_t i = 0; i < numLocalParticles; ++ i)
470  f64buffer[i] = bunch->Bf[i](1);
471  WRITEDATA(Float64, file_m, "By", f64buffer);
472 
473  for (size_t i = 0; i < numLocalParticles; ++ i)
474  f64buffer[i] = bunch->Bf[i](2);
475  WRITEDATA(Float64, file_m, "Bz", f64buffer);
476 
477  }
478 
480  if (Options::rhoDump) {
481  NDIndex<3> idx = bunch->getFieldLayout().getLocalNDIndex();
482  NDIndex<3> elem;
483  h5_err_t herr = H5Block3dSetView(
484  file_m,
485  idx[0].min(), idx[0].max(),
486  idx[1].min(), idx[1].max(),
487  idx[2].min(), idx[2].max());
488  reportOnError(herr, __FILE__, __LINE__);
489 
490  std::unique_ptr<h5_float64_t[]> data(new h5_float64_t[(idx[0].max() + 1) * (idx[1].max() + 1) * (idx[2].max() + 1)]);
491 
492  int ii = 0;
493  // h5block uses the fortran convention of storing data:
494  // INTEGER, DIMENSION(2,3) :: a
495  // => {a(1,1), a(2,1), a(1,2), a(2,2), a(1,3), a(2,3)}
496  for (int i = idx[2].min(); i <= idx[2].max(); ++ i) {
497  for (int j = idx[1].min(); j <= idx[1].max(); ++ j) {
498  for (int k = idx[0].min(); k <= idx[0].max(); ++ k) {
499  data[ii] = bunch->getRho(k, j, i);
500  ++ ii;
501  }
502  }
503  }
504  herr = H5Block3dWriteScalarFieldFloat64(file_m, "rho", data.get());
505  reportOnError(herr, __FILE__, __LINE__);
506 
508  herr = H5Block3dSetFieldOrigin(file_m, "rho",
509  (h5_float64_t)bunch->get_origin()(0),
510  (h5_float64_t)bunch->get_origin()(1),
511  (h5_float64_t)bunch->get_origin()(2));
512  reportOnError(herr, __FILE__, __LINE__);
513 
514  herr = H5Block3dSetFieldSpacing(file_m, "rho",
515  (h5_float64_t)bunch->get_hr()(0),
516  (h5_float64_t)bunch->get_hr()(1),
517  (h5_float64_t)bunch->get_hr()(2));
518  reportOnError(herr, __FILE__, __LINE__);
519 
520  }
521 }
Tps< T > cos(const Tps< T > &x)
Cosine.
Definition: TpsMath.h:129
Tps< T > sin(const Tps< T > &x)
Sine.
Definition: TpsMath.h:111
#define OPAL_PROJECT_VERSION
Definition: OPALconfig.h:5
#define OPAL_PROJECT_NAME
Definition: OPALconfig.h:2
#define WRITESTRINGFILEATTRIB(file, name, value)
Definition: H5PartWrapper.h:24
#define WRITESTEPATTRIB(type, file, name, value, length)
Definition: H5PartWrapper.h:29
#define REPORTONERROR(rc)
Definition: H5PartWrapper.h:22
#define READSTEPATTRIB(type, file, name, value)
Definition: H5PartWrapper.h:27
#define WRITEDATA(type, file, name, value)
Definition: H5PartWrapper.h:32
#define WRITEFILEATTRIB(type, file, name, value, length)
Definition: H5PartWrapper.h:25
#define WRITESTRINGSTEPATTRIB(file, name, value)
Definition: H5PartWrapper.h:28
#define READDATA(type, file, name, value)
Definition: H5PartWrapper.h:31
#define READFILEATTRIB(type, file, name, value)
Definition: H5PartWrapper.h:23
T::PETE_Expr_t::PETE_Return_t max(const PETE_Expr< T > &expr, NDIndex< D > &loc)
Definition: ReductionLoc.h:84
T::PETE_Expr_t::PETE_Return_t min(const PETE_Expr< T > &expr, NDIndex< D > &loc)
Definition: ReductionLoc.h:76
Inform & endl(Inform &inf)
Definition: Inform.cpp:42
#define ERRORMSG(msg)
Definition: IpplInfo.h:350
const std::string name
constexpr double epsilon_0
The permittivity of vacuum in As/Vm.
Definition: Physics.h:57
constexpr double c
The velocity of light in m/s.
Definition: Physics.h:51
constexpr double pi
The value of.
Definition: Physics.h:30
int psDumpFreq
The frequency to dump the phase space, i.e.dump data when steppsDumpFreq==0.
Definition: Options.cpp:39
bool rhoDump
Definition: Options.cpp:63
bool ebDump
Definition: Options.cpp:65
Vector_t getTaitBryantAngles(Quaternion rotation, const std::string &)
Definition: Util.cpp:102
T * c_data(std::vector< T, A > &v)
Definition: Util.h:201
std::string getGitRevision()
Definition: Util.cpp:18
ParticleAttrib< Vector_t > Ef
double get_meanKineticEnergy() const
ParticlePos_t & R
ParticleAttrib< int > Bin
Vector_t RefPartP_m
double getQ() const
Access to reference data.
size_t getLocalNum() const
void set_sPos(double s)
void setLocalTrackStep(long long n)
step in a TRACK command
size_t getTotalNum() const
double get_gamma() const
ParticleAttrib< Vector_t > P
ParticleAttrib< ParticleType > PType
ParticleAttrib< ParticleOrigin > POrigin
Vector_t get_rrms() const
Vector_t RefPartR_m
ParticleAttrib< double > Q
void calcBeamParameters()
virtual double getRho(int x, int y, int z)=0
Vector_t get_origin() const
Vector_t get_prms() const
double getdE() const
virtual FieldLayout_t & getFieldLayout()=0
double getCharge() const
get the total charge per simulation particle
CoordinateSystemTrafo toLabTrafo_m
long long getLocalTrackStep() const
Vector_t get_norm_emit() const
Vector_t get_maxExtent() const
void get_PBounds(Vector_t &min, Vector_t &max) const
void setGlobalTrackStep(long long n)
step in multiple TRACK commands
Vector_t get_centroid() const
Vector_t get_pmean() const
ParticleAttrib< Vector_t > Bf
long long getGlobalTrackStep() const
void setT(double t)
virtual Vector_t get_hr() const
ParticleIndex_t & ID
Vector_t get_emit() const
double get_sPos() const
double getM() const
double getT() const
void setRestartDumpFreq(const int &N)
set the dump frequency as found in restart file
Definition: OpalData.cpp:340
double getGlobalPhaseShift()
units: (sec)
Definition: OpalData.cpp:451
void setGlobalPhaseShift(double shift)
units: (sec)
Definition: OpalData.cpp:446
static OpalData * getInstance()
Definition: OpalData.cpp:195
void setRestartStep(int s)
store the location where to restart
Definition: OpalData.cpp:319
Quaternion getRotation() const
Vector_t rotate(const Vector_t &) const
Definition: Quaternion.cpp:122
Quaternion conjugate() const
Definition: Quaternion.h:105
size_t getNumParticles() const
void open(h5_int32_t flags)
h5_int64_t numSteps_m
Definition: H5PartWrapper.h:76
h5_file_t file_m
Definition: H5PartWrapper.h:73
static void reportOnError(h5_int64_t rc, const char *file, int line)
Definition: H5PartWrapper.h:83
void writeStepData(PartBunchBase< double, 3 > *)
void readStepHeader(PartBunchBase< double, 3 > *)
virtual void writeHeader()
virtual void writeStep(PartBunchBase< double, 3 > *, const std::map< std::string, double > &additionalStepAttributes)
void readStepData(PartBunchBase< double, 3 > *, h5_ssize_t, h5_ssize_t)
H5PartWrapperForPT(const std::string &fileName, h5_int32_t flags=H5_O_WRONLY)
virtual void readHeader()
virtual void readStep(PartBunchBase< double, 3 > *, h5_ssize_t firstParticle, h5_ssize_t lastParticle)
void writeStepHeader(PartBunchBase< double, 3 > *, const std::map< std::string, double > &)
The base class for all OPAL exceptions.
Definition: OpalException.h:28
NDIndex< Dim > getLocalNDIndex()