OPAL (Object Oriented Parallel Accelerator Library)  2024.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 
18 H5PartWrapperForPC::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 
29 H5PartWrapperForPC::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 
43 H5PartWrapperForPC::~H5PartWrapperForPC()
44 { }
45 
46 void 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 
73 void 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 
85 void 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 
148 void 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 
225 void 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 
311 void 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 
324 void 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 
464 void 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 }
OpalData::getInstance
static OpalData * getInstance()
Definition: OpalData.cpp:196
sqrt
Tps< T > sqrt(const Tps< T > &x)
Square root.
Definition: TpsMath.h:91
Physics::c
constexpr double c
The velocity of light in m/s.
Definition: Physics.h:45
H5PartWrapper::getNumParticles
size_t getNumParticles() const
Definition: H5PartWrapper.cpp:524
WRITEFILEATTRIB
#define WRITEFILEATTRIB(type, file, name, value, length)
Definition: H5PartWrapper.h:25
H5PartWrapperForPC::readHeader
virtual void readHeader()
Definition: H5PartWrapperForPC.cpp:46
PartBunchBase::setLocalTrackStep
void setLocalTrackStep(long long n)
step in a TRACK command
Definition: PartBunchBase.hpp:1395
OPAL_PROJECT_VERSION
#define OPAL_PROJECT_VERSION
Definition: config.h.in:5
PartBunchBase::setGlobalTrackStep
void setGlobalTrackStep(long long n)
step in multiple TRACK commands
Definition: PartBunchBase.hpp:1383
PartBunchBase::get_emit
Vector_t get_emit() const
Definition: PartBunchBase.hpp:1151
PyOpal::PyOpalObjectNS::attributeName
std::map< AttributeType, std::string > attributeName
Definition: PyOpalObject.cpp:22
Options::psDumpFrame
DumpFrame psDumpFrame
flag to decide in which coordinate frame the phase space will be dumped for OPAL-cycl ...
Definition: Options.cpp:45
PartBunchBase::P
ParticleAttrib< Vector_t > P
Definition: PartBunchBase.h:576
Physics::epsilon_0
constexpr double epsilon_0
The permittivity of vacuum in As/Vm.
Definition: Physics.h:51
PartBunchBase::hasBinning
bool hasBinning() const
Definition: PartBunchBase.hpp:248
PartBunchBase::setT
void setT(double t)
Definition: PartBunchBase.hpp:1031
Util::getGitRevision
std::string getGitRevision()
Definition: Util.cpp:33
H5PartWrapperForPC::writeStep
virtual void writeStep(PartBunchBase< double, 3 > *, const std::map< std::string, double > &additionalStepAttributes)
Definition: H5PartWrapperForPC.cpp:311
PartBunchBase::Ef
ParticleAttrib< Vector_t > Ef
Definition: PartBunchBase.h:580
ERRORMSG
#define ERRORMSG(msg)
Definition: IpplInfo.h:350
H5PartWrapperForPC::writeStepHeader
void writeStepHeader(PartBunchBase< double, 3 > *, const std::map< std::string, double > &)
Definition: H5PartWrapperForPC.cpp:324
PartBunchBase::get_centroid
Vector_t get_centroid() const
Definition: PartBunchBase.hpp:1105
PartBunchBase::bunchNum
ParticleAttrib< short > bunchNum
Definition: PartBunchBase.h:590
PartBunchBase::getRho
virtual double getRho(int x, int y, int z)=0
H5PartWrapperForPC::readStep
virtual void readStep(PartBunchBase< double, 3 > *, h5_ssize_t firstParticle, h5_ssize_t lastParticle)
Definition: H5PartWrapperForPC.cpp:73
PartBunchBase::M
ParticleAttrib< double > M
Definition: PartBunchBase.h:578
OPAL_PROJECT_NAME
#define OPAL_PROJECT_NAME
Definition: config.h.in:2
PartBunchBase::getLocalTrackStep
long long getLocalTrackStep() const
Definition: PartBunchBase.hpp:1407
PartBunchBase::get_pmean
Vector_t get_pmean() const
Definition: PartBunchBase.hpp:1135
PartBunchBase::calcBeamParameters
void calcBeamParameters()
Definition: PartBunchBase.hpp:1267
WRITESTRINGSTEPATTRIB
#define WRITESTRINGSTEPATTRIB(file, name, value)
Definition: H5PartWrapper.h:28
PartBunchBase::get_sPos
double get_sPos() const
Definition: PartBunchBase.hpp:1049
PartBunchBase::getdE
double getdE() const
Definition: PartBunchBase.hpp:1700
H5PartWrapperForPC::writeStepData
void writeStepData(PartBunchBase< double, 3 > *)
Definition: H5PartWrapperForPC.cpp:464
PartBunchBase::get_rrms
Vector_t get_rrms() const
Definition: PartBunchBase.hpp:1111
PartBunchBase::getQ
double getQ() const
Access to reference data.
Definition: PartBunchBase.hpp:1630
PartBunchBase::get_hr
virtual Vector_t get_hr() const
Definition: PartBunchBase.hpp:1232
Physics::pi
constexpr double pi
The value of .
Definition: Physics.h:30
Util::c_data
T * c_data(std::vector< T, A > &v)
Definition: Util.h:247
endl
Inform & endl(Inform &inf)
Definition: Inform.cpp:42
REPORTONERROR
#define REPORTONERROR(rc)
Definition: H5PartWrapper.h:22
PartBunchBase::getM
double getM() const
Definition: PartBunchBase.hpp:1636
PartBunchBase::get_prms
Vector_t get_prms() const
Definition: PartBunchBase.hpp:1129
min
T::PETE_Expr_t::PETE_Return_t min(const PETE_Expr< T > &expr, NDIndex< D > &loc)
Definition: ReductionLoc.h:76
Options::ebDump
bool ebDump
Definition: Options.cpp:65
max
T::PETE_Expr_t::PETE_Return_t max(const PETE_Expr< T > &expr, NDIndex< D > &loc)
Definition: ReductionLoc.h:84
PartBunchBase::setSteptoLastInj
void setSteptoLastInj(int n)
Definition: PartBunchBase.hpp:1506
PartBunchBase::Bf
ParticleAttrib< Vector_t > Bf
Definition: PartBunchBase.h:583
H5PartWrapperForPC::meanMomentum_m
h5_float64_t meanMomentum_m
Definition: H5PartWrapperForPC.h:50
PartBunchBase::get_gamma
double get_gamma() const
Definition: PartBunchBase.hpp:1061
PartBunchBase::set_sPos
void set_sPos(double s)
Definition: PartBunchBase.hpp:1055
READDATA
#define READDATA(type, file, name, value)
Definition: H5PartWrapper.h:31
H5PartWrapper::file_m
h5_file_t file_m
Definition: H5PartWrapper.h:73
WRITEDATA
#define WRITEDATA(type, file, name, value)
Definition: H5PartWrapper.h:32
PartBunchBase::getT
double getT() const
Definition: PartBunchBase.hpp:1043
PartBunchBase::getTotalNum
size_t getTotalNum() const
Definition: PartBunchBase.hpp:1906
OpalException
The base class for all OPAL exceptions.
Definition: OpalException.h:28
PartBunchBase::get_origin
Vector_t get_origin() const
Definition: PartBunchBase.hpp:1093
FieldLayout::getLocalNDIndex
NDIndex< Dim > getLocalNDIndex()
Definition: FieldLayout.hpp:1579
PartBunchBase::getLocalNum
size_t getLocalNum() const
Definition: PartBunchBase.hpp:1911
Options::psDumpFreq
int psDumpFreq
The frequency to dump the phase space, i.e.dump data when steppsDumpFreq==0.
Definition: Options.cpp:39
H5PartWrapperForPC::readStepData
void readStepData(PartBunchBase< double, 3 > *, h5_ssize_t, h5_ssize_t)
Definition: H5PartWrapperForPC.cpp:148
OpalData::setRestartDumpFreq
void setRestartDumpFreq(const int &N)
set the dump frequency as found in restart file
Definition: OpalData.cpp:341
PartBunchBase::getGlobalTrackStep
long long getGlobalTrackStep() const
Definition: PartBunchBase.hpp:1389
PartBunchBase::Q
ParticleAttrib< double > Q
Definition: PartBunchBase.h:577
PartBunchBase::getFieldLayout
virtual FieldLayout_t & getFieldLayout()=0
H5PartWrapper::predecessorOPALFlavour_m
std::string predecessorOPALFlavour_m
Definition: H5PartWrapper.h:75
H5PartWrapperForPC::elevation_m
h5_float64_t elevation_m
Definition: H5PartWrapperForPC.h:52
H5PartWrapperForPC::H5PartWrapperForPC
H5PartWrapperForPC(const std::string &fileName, h5_int32_t flags=H5_O_WRONLY)
Definition: H5PartWrapperForPC.cpp:18
PartBunchBase::get_PBounds
void get_PBounds(Vector_t &min, Vector_t &max) const
Definition: PartBunchBase.hpp:1261
READFILEATTRIB
#define READFILEATTRIB(type, file, name, value)
Definition: H5PartWrapper.h:23
READSTEPATTRIB
#define READSTEPATTRIB(type, file, name, value)
Definition: H5PartWrapper.h:27
PartBunchBase::get_rmean
Vector_t get_rmean() const
Definition: PartBunchBase.hpp:1123
PartBunchBase::R
ParticlePos_t & R
Definition: PartBunchBase.h:572
WRITESTRINGFILEATTRIB
#define WRITESTRINGFILEATTRIB(file, name, value)
Definition: H5PartWrapper.h:24
PartBunchBase::setNumBunch
void setNumBunch(short n)
Definition: PartBunchBase.hpp:1413
PartBunchBase::get_norm_emit
Vector_t get_norm_emit() const
Definition: PartBunchBase.hpp:1157
H5PartWrapperForPC::writeHeader
virtual void writeHeader()
Definition: H5PartWrapperForPC.cpp:225
PartBunchBase::ID
ParticleIndex_t & ID
Definition: PartBunchBase.h:573
WRITESTEPATTRIB
#define WRITESTEPATTRIB(type, file, name, value, length)
Definition: H5PartWrapper.h:29
PartBunchBase::RefPartR_m
Vector_t RefPartR_m
Definition: PartBunchBase.h:592
Units::eV2GeV
constexpr double eV2GeV
Definition: Units.h:71
PartBunchBase::getCharge
double getCharge() const
get the total charge per simulation particle
Definition: PartBunchBase.hpp:1317
Options::rhoDump
bool rhoDump
Definition: Options.cpp:63
PartBunchBase::get_meanKineticEnergy
double get_meanKineticEnergy() const
Definition: PartBunchBase.hpp:1067
PartBunchBase::Bin
ParticleAttrib< int > Bin
Definition: PartBunchBase.h:584
PartBunchBase::getNumBunch
short getNumBunch() const
Definition: PartBunchBase.hpp:1421
PartBunchBase::get_maxExtent
Vector_t get_maxExtent() const
Definition: PartBunchBase.hpp:1099
pow
Tps< T > pow(const Tps< T > &x, int y)
Integer power.
Definition: TpsMath.h:76
PartBunchBase::RefPartP_m
Vector_t RefPartP_m
Definition: PartBunchBase.h:593
H5PartWrapper::open
void open(h5_int32_t flags)
Definition: H5PartWrapper.cpp:71
H5PartWrapper::numSteps_m
h5_int64_t numSteps_m
Definition: H5PartWrapper.h:76
H5PartWrapperForPC::readStepHeader
void readStepHeader(PartBunchBase< double, 3 > *)
Definition: H5PartWrapperForPC.cpp:85
Units::eV2MeV
constexpr double eV2MeV
Definition: Units.h:77
PartBunchBase::getSteptoLastInj
int getSteptoLastInj() const
Definition: PartBunchBase.hpp:1512
Generated on Tue Jul 2 2024 15:05:26 for OPAL (Object Oriented Parallel Accelerator Library) by   doxygen 1.8.5