OPAL (Object Oriented Parallel Accelerator Library)  2.2.0
OPAL
StatWriter.cpp
Go to the documentation of this file.
1 #include "StatWriter.h"
2 
3 #include "Utilities/Timer.h"
4 
5 #include <sstream>
6 
7 StatWriter::StatWriter(const std::string& fname, bool restart)
8  : StatBaseWriter(fname, restart)
9 { }
10 
11 
12 void StatWriter::fillHeader(const losses_t &losses) {
13 
14  if (this->hasColumns()) {
15  return;
16  }
17 
18  columns_m.addColumn("t", "double", "ns", "Time");
19  columns_m.addColumn("s", "double", "m", "Path length");
20  columns_m.addColumn("numParticles", "long", "1", "Number of Macro Particles");
21  columns_m.addColumn("charge", "double", "1", "Bunch Charge");
22  columns_m.addColumn("energy", "double", "MeV", "Mean Bunch Energy");
23 
24  columns_m.addColumn("rms_x", "double", "m", "RMS Beamsize in x");
25  columns_m.addColumn("rms_y", "double", "m", "RMS Beamsize in y");
26  columns_m.addColumn("rms_s", "double", "m", "RMS Beamsize in s");
27 
28  columns_m.addColumn("rms_px", "double", "1", "RMS Normalized Momenta in x");
29  columns_m.addColumn("rms_py", "double", "1", "RMS Normalized Momenta in y");
30  columns_m.addColumn("rms_ps", "double", "1", "RMS Normalized Momenta in s");
31 
32  columns_m.addColumn("emit_x", "double", "m", "Normalized Emittance x");
33  columns_m.addColumn("emit_y", "double", "m", "Normalized Emittance y");
34  columns_m.addColumn("emit_s", "double", "m", "Normalized Emittance s");
35 
36  columns_m.addColumn("mean_x", "double", "m", "Mean Beam Position in x");
37  columns_m.addColumn("mean_y", "double", "m", "Mean Beam Position in y");
38  columns_m.addColumn("mean_s", "double", "m", "Mean Beam Position in s");
39 
40  columns_m.addColumn("ref_x", "double", "m", "x coordinate of reference particle in lab cs");
41  columns_m.addColumn("ref_y", "double", "m", "y coordinate of reference particle in lab cs");
42  columns_m.addColumn("ref_z", "double", "m", "z coordinate of reference particle in lab cs");
43 
44  columns_m.addColumn("ref_px", "double", "1", "x momentum of reference particle in lab cs");
45  columns_m.addColumn("ref_py", "double", "1", "y momentum of reference particle in lab cs");
46  columns_m.addColumn("ref_pz", "double", "1", "z momentum of reference particle in lab cs");
47 
48  columns_m.addColumn("max_x", "double", "m", "Max Beamsize in x");
49  columns_m.addColumn("max_y", "double", "m", "Max Beamsize in y");
50  columns_m.addColumn("max_s", "double", "m", "Max Beamsize in s");
51 
52  columns_m.addColumn("xpx", "double", "1", "Correlation xpx");
53  columns_m.addColumn("ypy", "double", "1", "Correlation ypy");
54  columns_m.addColumn("zpz", "double", "1", "Correlation zpz");
55 
56  columns_m.addColumn("Dx", "double", "m", "Dispersion in x");
57  columns_m.addColumn("DDx", "double", "1", "Derivative of dispersion in x");
58  columns_m.addColumn("Dy", "double", "m", "Dispersion in y");
59  columns_m.addColumn("DDy", "double", "1", "Derivative of dispersion in y");
60 
61  columns_m.addColumn("Bx_ref", "double", "T", "Bx-Field component of ref particle");
62  columns_m.addColumn("By_ref", "double", "T", "By-Field component of ref particle");
63  columns_m.addColumn("Bz_ref", "double", "T", "Bz-Field component of ref particle");
64 
65  columns_m.addColumn("Ex_ref", "double", "MV/m", "Ex-Field component of ref particle");
66  columns_m.addColumn("Ey_ref", "double", "MV/m", "Ey-Field component of ref particle");
67  columns_m.addColumn("Ez_ref", "double", "MV/m", "Ez-Field component of ref particle");
68 
69  columns_m.addColumn("dE", "double", "MeV", "energy spread of the beam");
70  columns_m.addColumn("dt", "double", "ns", "time step size");
71  columns_m.addColumn("partsOutside", "double", "1", "outside n*sigma of the beam");
72 
73  if (OpalData::getInstance()->isInOPALCyclMode() &&
74  Ippl::getNodes() == 1) {
75  columns_m.addColumn("R0_x", "double", "m", "R0 Particle position in x");
76  columns_m.addColumn("R0_y", "double", "m", "R0 Particle position in y");
77  columns_m.addColumn("R0_s", "double", "m", "R0 Particle position in z");
78 
79  columns_m.addColumn("P0_x", "double", "1", "R0 Particle momentum in x");
80  columns_m.addColumn("P0_y", "double", "1", "R0 Particle momentum in y");
81  columns_m.addColumn("P0_s", "double", "1", "R0 Particle momentum in z");
82  }
83 
84  if (OpalData::getInstance()->isInOPALCyclMode()) {
85  columns_m.addColumn("halo_x", "double", "1", "Halo in x");
86  columns_m.addColumn("halo_y", "double", "1", "Halo in y");
87  columns_m.addColumn("halo_z", "double", "1", "Halo in z");
88 
89  columns_m.addColumn("azimuth", "double", "deg",
90  "Azimuth in global coordinates");
91  }
92 
93  for (size_t i = 0; i < losses.size(); ++ i) {
94  columns_m.addColumn(losses[i].first, "long", "1",
95  "Number of lost particles in element");
96  }
97 
98  if ( mode_m == std::ios::app )
99  return;
100 
101  OPALTimer::Timer simtimer;
102  std::string dateStr(simtimer.date());
103  std::string timeStr(simtimer.time());
104 
105  std::stringstream ss;
106  ss << "Statistics data '"
107  << OpalData::getInstance()->getInputFn()
108  << "' " << dateStr << " " << timeStr;
109 
110  this->addDescription(ss.str(), "stat parameters");
111 
112  this->addDefaultParameters();
113 
114 
115  this->addInfo("ascii", 1);
116 }
117 
118 
119 void StatWriter::write(PartBunchBase<double, 3> *beam, Vector_t FDext[],
120  const losses_t &losses, const double& azimuth)
121 {
123  beam->calcBeamParameters();
124 
125  double Ekin = beam->get_meanKineticEnergy();
126 
127  size_t npOutside = 0;
128  if (Options::beamHaloBoundary > 0)
129  npOutside = beam->calcNumPartsOutside(Options::beamHaloBoundary*beam->get_rrms());
130 
131  double pathLength = beam->get_sPos();
132 
135 
136  double Q = beam->getCharge();
137 
138  if (Ippl::myNode() != 0) {
139  return;
140  }
141 
142  fillHeader(losses);
143 
144  this->open();
145 
146  this->writeHeader();
147 
148  columns_m.addColumnValue("t", beam->getT() * 1e9); // 1
149  columns_m.addColumnValue("s", pathLength); // 2
150  columns_m.addColumnValue("numParticles", beam->getTotalNum()); // 3
151  columns_m.addColumnValue("charge", Q); // 4
152  columns_m.addColumnValue("energy", Ekin); // 5
153 
154  columns_m.addColumnValue("rms_x", beam->get_rrms()(0)); // 6
155  columns_m.addColumnValue("rms_y", beam->get_rrms()(1)); // 7
156  columns_m.addColumnValue("rms_s", beam->get_rrms()(2)); // 8
157 
158  columns_m.addColumnValue("rms_px", beam->get_prms()(0)); // 9
159  columns_m.addColumnValue("rms_py", beam->get_prms()(1)); // 10
160  columns_m.addColumnValue("rms_ps", beam->get_prms()(2)); // 11
161 
162  columns_m.addColumnValue("emit_x", beam->get_norm_emit()(0)); // 12
163  columns_m.addColumnValue("emit_y", beam->get_norm_emit()(1)); // 13
164  columns_m.addColumnValue("emit_s", beam->get_norm_emit()(2)); // 14
165 
166  columns_m.addColumnValue("mean_x", beam->get_rmean()(0) ); // 15
167  columns_m.addColumnValue("mean_y", beam->get_rmean()(1) ); // 16
168  columns_m.addColumnValue("mean_s", beam->get_rmean()(2) ); // 17
169 
170  columns_m.addColumnValue("ref_x", beam->RefPartR_m(0)); // 18
171  columns_m.addColumnValue("ref_y", beam->RefPartR_m(1)); // 19
172  columns_m.addColumnValue("ref_z", beam->RefPartR_m(2)); // 20
173 
174  columns_m.addColumnValue("ref_px", beam->RefPartP_m(0)); // 21
175  columns_m.addColumnValue("ref_py", beam->RefPartP_m(1)); // 22
176  columns_m.addColumnValue("ref_pz", beam->RefPartP_m(2)); // 23
177 
178  columns_m.addColumnValue("max_x", beam->get_maxExtent()(0)); // 24
179  columns_m.addColumnValue("max_y", beam->get_maxExtent()(1)); // 25
180  columns_m.addColumnValue("max_s", beam->get_maxExtent()(2)); // 26
181 
182  // Write out Courant Snyder parameters.
183  columns_m.addColumnValue("xpx", beam->get_rprms()(0)); // 27
184  columns_m.addColumnValue("ypy", beam->get_rprms()(1)); // 28
185  columns_m.addColumnValue("zpz", beam->get_rprms()(2)); // 29
186 
187  // Write out dispersion.
188  columns_m.addColumnValue("Dx", beam->get_Dx()); // 30
189  columns_m.addColumnValue("DDx", beam->get_DDx()); // 31
190  columns_m.addColumnValue("Dy", beam->get_Dy()); // 32
191  columns_m.addColumnValue("DDy", beam->get_DDy()); // 33
192 
193  // Write head/reference particle/tail field information.
194  columns_m.addColumnValue("Bx_ref", FDext[0](0)); // 34 B-ref x
195  columns_m.addColumnValue("By_ref", FDext[0](1)); // 35 B-ref y
196  columns_m.addColumnValue("Bz_ref", FDext[0](2)); // 36 B-ref z
197 
198  columns_m.addColumnValue("Ex_ref", FDext[1](0)); // 37 E-ref x
199  columns_m.addColumnValue("Ey_ref", FDext[1](1)); // 38 E-ref y
200  columns_m.addColumnValue("Ez_ref", FDext[1](2)); // 39 E-ref z
201 
202  columns_m.addColumnValue("dE", beam->getdE()); // 40 dE energy spread
203  columns_m.addColumnValue("dt", beam->getdT() * 1e9); // 41 dt time step size
204  columns_m.addColumnValue("partsOutside", npOutside); // 42 number of particles outside n*sigma
205 
206  if (OpalData::getInstance()->isInOPALCyclMode()) {
207  if (Ippl::getNodes() == 1) {
208  if (beam->getLocalNum() > 0) {
209  columns_m.addColumnValue("R0_x", beam->R[0](0));
210  columns_m.addColumnValue("R0_y", beam->R[0](1));
211  columns_m.addColumnValue("R0_s", beam->R[0](2));
212  columns_m.addColumnValue("P0_x", beam->P[0](0));
213  columns_m.addColumnValue("P0_y", beam->P[0](1));
214  columns_m.addColumnValue("P0_s", beam->P[0](2));
215  } else {
216  columns_m.addColumnValue("R0_x", 0.0);
217  columns_m.addColumnValue("R0_y", 0.0);
218  columns_m.addColumnValue("R0_s", 0.0);
219  columns_m.addColumnValue("P0_x", 0.0);
220  columns_m.addColumnValue("P0_y", 0.0);
221  columns_m.addColumnValue("P0_s", 0.0);
222  }
223  }
224  Vector_t halo = beam->get_halo();
225  columns_m.addColumnValue("halo_x", halo(0));
226  columns_m.addColumnValue("halo_y", halo(1));
227  columns_m.addColumnValue("halo_z", halo(2));
228 
229  columns_m.addColumnValue("azimuth", azimuth);
230  }
231 
232  for(size_t i = 0; i < losses.size(); ++ i) {
233  long unsigned int loss = losses[i].second;
234  columns_m.addColumnValue(losses[i].first, loss);
235  }
236 
237  this->writeRow();
238 
239  this->close();
240 }
241 
242 
243 void StatWriter::write(EnvelopeBunch &beam, Vector_t FDext[],
244  double sposHead, double sposRef, double sposTail)
245 {
246  //FIXME https://gitlab.psi.ch/OPAL/src/issues/245
247 
249  beam.calcBeamParameters();
250  beam.gatherLoadBalanceStatistics();
251 
252  double en = beam.get_meanKineticEnergy() * 1e-6;
253  long unsigned int totalnum = beam.getTotalNum();
254  double Q = totalnum * beam.getChargePerParticle();
255 
256 
257  if (Ippl::myNode() != 0) {
258  return;
259  }
260 
261  fillHeader();
262 
263  this->open();
264 
265  this->writeHeader();
266 
267  columns_m.addColumnValue("t", beam.getT() * 1e9); // 1
268  columns_m.addColumnValue("s", sposRef); // 2
269  columns_m.addColumnValue("numParticles", totalnum); // 3
270  columns_m.addColumnValue("charge", Q); // 4
271  columns_m.addColumnValue("energy", en); // 5
272 
273  columns_m.addColumnValue("rms_x", beam.get_rrms()(0)); // 6
274  columns_m.addColumnValue("rms_y", beam.get_rrms()(1)); // 7
275  columns_m.addColumnValue("rms_s", beam.get_rrms()(2)); // 8
276 
277  columns_m.addColumnValue("rms_px", beam.get_prms()(0)); // 9
278  columns_m.addColumnValue("rms_py", beam.get_prms()(1)); // 10
279  columns_m.addColumnValue("rms_ps", beam.get_prms()(2)); // 11
280 
281  columns_m.addColumnValue("emit_x", beam.get_norm_emit()(0)); // 12
282  columns_m.addColumnValue("emit_y", beam.get_norm_emit()(1)); // 13
283  columns_m.addColumnValue("emit_s", beam.get_norm_emit()(2)); // 14
284 
285  columns_m.addColumnValue("mean_x", beam.get_rmean()(0) ); // 15
286  columns_m.addColumnValue("mean_y", beam.get_rmean()(1) ); // 16
287  columns_m.addColumnValue("mean_s", beam.get_rmean()(2) ); // 17
288 
289  columns_m.addColumnValue("ref_x", 0.0); // 18
290  columns_m.addColumnValue("ref_y", 0.0); // 19
291  columns_m.addColumnValue("ref_z", 0.0); // 20
292 
293  columns_m.addColumnValue("ref_px", 0.0); // 21
294  columns_m.addColumnValue("ref_py", 0.0); // 22
295  columns_m.addColumnValue("ref_pz", 0.0); // 23
296 
297  columns_m.addColumnValue("max_x", beam.get_maxExtent()(0)); // 24
298  columns_m.addColumnValue("max_y", beam.get_maxExtent()(1)); // 25
299  columns_m.addColumnValue("max_s", beam.get_maxExtent()(2)); // 26
300 
301  columns_m.addColumnValue("xpx", 0.0); // 27
302  columns_m.addColumnValue("ypy", 0.0); // 28
303  columns_m.addColumnValue("zpz", 0.0); // 29
304 
305  // Write out dispersion.
306  columns_m.addColumnValue("Dx", beam.get_Dx()); // 30
307  columns_m.addColumnValue("DDx", beam.get_DDx()); // 31
308  columns_m.addColumnValue("Dy", beam.get_Dy()); // 32
309  columns_m.addColumnValue("DDy", beam.get_DDy()); // 33
310 
311  // Write head/reference particle/tail field information.
312  columns_m.addColumnValue("Bx_ref", FDext[0](0)); // 34 B-ref x
313  columns_m.addColumnValue("By_ref", FDext[0](1)); // 35 B-ref y
314  columns_m.addColumnValue("Bz_ref", FDext[0](2)); // 36 B-ref z
315 
316  columns_m.addColumnValue("Ex_ref", FDext[1](0)); // 37 E-ref x
317  columns_m.addColumnValue("Ey_ref", FDext[1](1)); // 38 E-ref y
318  columns_m.addColumnValue("Ez_ref", FDext[1](2)); // 39 E-ref z
319 
320  columns_m.addColumnValue("dE", beam.get_dEdt()); // 40 dE energy spread
321  columns_m.addColumnValue("dt", 0.0); // 41 dt time step size
322  columns_m.addColumnValue("partsOutside", 0.0); // 42 number of particles outside n*sigma
323 
324  this->writeRow();
325 
326  this->close();
327 }
PartBunchBase::get_rprms
Vector_t get_rprms() const
Definition: PartBunchBase.hpp:1148
IpplInfo::getNodes
static int getNodes()
Definition: IpplInfo.cpp:773
PartBunchBase::get_Dx
double get_Dx() const
Definition: PartBunchBase.hpp:1206
PartBunchBase::P
ParticleAttrib< Vector_t > P
Definition: PartBunchBase.h:484
PartBunchBase::getT
double getT() const
Definition: PartBunchBase.hpp:1094
EnvelopeBunch::get_dEdt
double get_dEdt()
returns the energy spread
Definition: EnvelopeBunch.h:172
SDDSWriter::hasColumns
bool hasColumns() const
Definition: SDDSWriter.h:190
EnvelopeBunch::get_meanKineticEnergy
double get_meanKineticEnergy()
returns the mean energy
Definition: EnvelopeBunch.h:170
EnvelopeBunch::get_rrms
Vector_t get_rrms()
returns RMS x,y,z
Definition: EnvelopeBunch.h:192
PartBunchBase::get_DDy
double get_DDy() const
Definition: PartBunchBase.hpp:1224
Physics::e
constexpr double e
The value of .
Definition: Physics.h:40
EnvelopeBunch::getT
double getT()
returns the current time of the bunch
Definition: EnvelopeBunch.h:168
PartBunchBase::RefPartR_m
Vector_t RefPartR_m
Definition: PartBunchBase.h:502
SDDSWriter::writeHeader
void writeHeader()
Write SDDS header.
Definition: SDDSWriter.cpp:130
EnvelopeBunch
core of the envelope tracker based on Rene Bakkers BET implementation
Definition: EnvelopeBunch.h:60
SDDSWriter::writeRow
void writeRow()
Definition: SDDSWriter.h:175
EnvelopeBunch::get_norm_emit
Vector_t get_norm_emit()
returns vector with normalized emittance
Definition: EnvelopeBunch.h:180
PartBunchBase::gatherLoadBalanceStatistics
void gatherLoadBalanceStatistics()
Definition: PartBunchBase.hpp:1236
PartBunchBase::get_rrms
Vector_t get_rrms() const
Definition: PartBunchBase.hpp:1142
PartBunchBase::getdT
double getdT() const
Definition: PartBunchBase.hpp:1076
PartBunchBase::get_meanKineticEnergy
double get_meanKineticEnergy() const
Definition: PartBunchBase.hpp:1118
IpplInfo::myNode
static int myNode()
Definition: IpplInfo.cpp:794
SDDSWriter::close
void close()
Definition: SDDSWriter.cpp:123
SDDSWriter::addDescription
void addDescription(const std::string &text, const std::string &content)
Definition: SDDSWriter.h:147
Options::beamHaloBoundary
double beamHaloBoundary
Definition: Options.cpp:24
PartBunchBase::getTotalNum
size_t getTotalNum() const
Definition: PartBunchBase.hpp:2225
SDDSWriter::columns_m
SDDSColumnSet columns_m
Definition: SDDSWriter.h:109
PartBunchBase::get_prms
Vector_t get_prms() const
Definition: PartBunchBase.hpp:1160
OPALTimer::Timer::date
std::string date() const
Return date.
Definition: Timer.cpp:35
PartBunchBase::get_rmean
Vector_t get_rmean() const
Definition: PartBunchBase.hpp:1154
OpalData::getInstance
static OpalData * getInstance()
Definition: OpalData.cpp:209
PartBunchBase::calcBeamParameters
void calcBeamParameters()
Definition: PartBunchBase.hpp:1259
PartBunchBase::get_maxExtent
Vector_t get_maxExtent() const
Definition: PartBunchBase.hpp:1130
PartBunchBase::get_halo
Vector_t get_halo() const
Definition: PartBunchBase.hpp:1194
PartBunchBase::get_Dy
double get_Dy() const
Definition: PartBunchBase.hpp:1212
OpalData::getInputFn
std::string getInputFn()
get opals input filename
Definition: OpalData.cpp:713
PartBunchBase::calcNumPartsOutside
size_t calcNumPartsOutside(Vector_t x)
returns the number of particles outside of a box defined by x
Definition: PartBunchBase.hpp:463
PartBunchBase::RefPartP_m
Vector_t RefPartP_m
Definition: PartBunchBase.h:503
SDDSColumnSet::addColumnValue
void addColumnValue(const std::string &name, const T &val)
Definition: SDDSColumnSet.h:46
StatWriter::write
void write(PartBunchBase< double, 3 > *beam, Vector_t FDext[], const losses_t &losses=losses_t(), const double &azimuth=-1)
Write statistical data.
Definition: StatWriter.cpp:119
StatWriter::fillHeader
void fillHeader(const losses_t &losses=losses_t())
Definition: StatWriter.cpp:12
StatWriter::StatWriter
StatWriter(const std::string &fname, bool restart)
Definition: StatWriter.cpp:7
PartBunchBase::get_norm_emit
Vector_t get_norm_emit() const
Definition: PartBunchBase.hpp:1188
PartBunchBase::get_DDx
double get_DDx() const
Definition: PartBunchBase.hpp:1218
PartBunchBase::getLocalNum
size_t getLocalNum() const
Definition: PartBunchBase.hpp:2230
PartBunchBase::getCharge
double getCharge() const
get the total charge per simulation particle
Definition: PartBunchBase.hpp:1456
SDDSWriter::mode_m
std::ios_base::openmode mode_m
First write to the statistics output file.
Definition: SDDSWriter.h:107
StatWriter::losses_t
std::vector< std::pair< std::string, unsigned int > > losses_t
Definition: StatWriter.h:11
OPALTimer::Timer::time
std::string time() const
Return time.
Definition: Timer.cpp:42
EnvelopeBunch::getChargePerParticle
double getChargePerParticle()
returns charge per slice
Definition: EnvelopeBunch.h:190
SDDSColumnSet::addColumn
void addColumn(const std::string &name, const std::string &type, const std::string &unit, const std::string &desc, std::ios_base::fmtflags flags=std::ios_base::scientific, unsigned short precision=15)
Definition: SDDSColumnSet.cpp:3
SDDSWriter::open
void open()
Definition: SDDSWriter.cpp:113
PartBunchBase::R
ParticlePos_t & R
Definition: PartBunchBase.h:479
SDDSWriter::addDefaultParameters
void addDefaultParameters()
Definition: SDDSWriter.cpp:190
EnvelopeBunch::getTotalNum
int getTotalNum()
returns the total number of slices
Definition: EnvelopeBunch.h:166
EnvelopeBunch::calcBeamParameters
void calcBeamParameters()
calculates envelope statistics
Definition: EnvelopeBunch.cpp:87
SDDSWriter::addInfo
void addInfo(const std::string &mode, const size_t &no_row_counts)
Definition: SDDSWriter.h:168
EnvelopeBunch::get_prms
Vector_t get_prms()
returns RMSP x,y,z
Definition: EnvelopeBunch.h:194
Generated on Sat Jan 4 2020 14:48:07 for OPAL (Object Oriented Parallel Accelerator Library) by   doxygen 1.8.5