Commit 3aa578f6 authored by Daniel Winklehner's avatar Daniel Winklehner
Browse files

Fixed minor issues with saving/restarting in local frame in CyclotronTracker.

parent aa2ff1af
......@@ -284,7 +284,7 @@ void ParallelCyclotronTracker::closeFiles() {
*/
void ParallelCyclotronTracker::visitRing(const Ring &ring) {
*gmsg << "Adding Ring" << endl;
*gmsg << "* ----------------------------- Adding Ring ------------------------------ *" << endl;
if (opalRing_m != NULL)
delete opalRing_m;
......@@ -346,7 +346,7 @@ void ParallelCyclotronTracker::visitRing(const Ring &ring) {
*/
void ParallelCyclotronTracker::visitCyclotron(const Cyclotron &cycl) {
*gmsg << "* --------- Cyclotron ------------------------------" << endl;
*gmsg << "* -------------------------- Adding Cyclotron ---------------------------- *" << endl;
Cyclotron *elptr = dynamic_cast<Cyclotron *>(cycl.clone());
myElements.push_back(elptr);
......@@ -396,10 +396,6 @@ void ParallelCyclotronTracker::visitCyclotron(const Cyclotron &cycl) {
// Restart a run:
} else {
// TEMP -DW
*gmsg << "previousH5Local = " << previousH5Local << endl;
*gmsg << "psDumpLocalFrame = " << Options::psDumpLocalFrame << endl;
// If the user wants to save the restarted run in local frame,
// make sure the previous h5 file was local too
if (Options::psDumpLocalFrame) {
......@@ -420,7 +416,10 @@ void ParallelCyclotronTracker::visitCyclotron(const Cyclotron &cycl) {
}
}
referencePtot = bega;
// Adjust some of the reference variables from the h5 file
referencePhi *= PIOVER180;
referencePsi *= PIOVER180;
referencePtot = bega;
if(referenceTheta <= -180.0 || referenceTheta > 180.0) {
......@@ -1006,7 +1005,7 @@ void ParallelCyclotronTracker::buildupFieldList(double BcParameter[], std::strin
(localpair->second).second = elptr;
// always put cyclotron as the first element in the list.
if(ElementType == "OPALRING" || ElementType == "CYCLOTRON") { // CYCLOTRON is TEMP -DW
if(ElementType == "OPALRING" || ElementType == "CYCLOTRON") {
sindex = FieldDimensions.begin();
} else {
sindex = FieldDimensions.end();
......@@ -5153,11 +5152,13 @@ void ParallelCyclotronTracker::initDistInGlobalFrame() {
checkNumPart(std::string("* After repartition: "));
itsBunch->calcBeamParameters_cycl();
*gmsg << endl << "* ********* Bunch After Calculation of Beam Parameterss in local frame: ********** *" << endl;
*gmsg << *itsBunch << endl;
itsBunch->R *= Vector_t(1000.0); // m --> mm
localToGlobal(itsBunch->R, phi, psi, meanR);
localToGlobal(itsBunch->P, phi, psi, Vector_t(0.0));
localToGlobal(itsBunch->P, phi, psi);
// Save initial distribution if not a restart
if(!OpalData::getInstance()->inRestartRun()) {
......@@ -5170,8 +5171,17 @@ void ParallelCyclotronTracker::initDistInGlobalFrame() {
step_m += 1;
}
// Print out the Bunch information at beginning of the run
// Print out the Bunch information at beginning of the run. Because the bunch information
// displays in units of m we have to change back and forth one more time.
// Furthermore it is my opinion that the same units should be used throughout OPAL. -DW
itsBunch->R *= Vector_t(0.001); // mm --> m
itsBunch->calcBeamParameters_cycl();
*gmsg << endl << "* ********* Bunch After Calculation of Beam Parameters in global frame: ********** *" << endl;
*gmsg << *itsBunch << endl;
itsBunch->R *= Vector_t(1000.0); // m --> mm
}
void ParallelCyclotronTracker::singleParticleDump() {
......@@ -5324,7 +5334,7 @@ void ParallelCyclotronTracker::bunchDumpStatData(){
// If we are saving in local frame, bunch and fields at the bunch center have to be rotated
// TODO: Make decision if we maybe want to always save statistics data in local frame? -DW
if(Options::psDumpLocalFrame) {
// -------------------- ----------- Do Transformations -------------------------------------- //
// -------------------- ----------- Do Transformations ---------------------------------- //
// Bunch (local) azimuth at meanR w.r.t. y-axis
phi = calculateAngle(meanP(0), meanP(1)) - 0.5 * pi;
......@@ -5428,7 +5438,7 @@ void ParallelCyclotronTracker::bunchDumpPhaseSpaceData() {
// at ref. R/Th/Z
psi / PIOVER180, // P_mean elevation
// at ref. R/Th/Z
0); // Flag localFrame
false); // Flag localFrame
// Tell user in which mode we are dumping
*gmsg << endl << "* Phase space dump " << lastDumpedStep_m
......@@ -5460,7 +5470,7 @@ void ParallelCyclotronTracker::bunchDumpPhaseSpaceData() {
// at ref. R/Th/Z
psi / PIOVER180, // P_mean elevation
// at ref. R/Th/Z
1); // Flag localFrame
true); // Flag localFrame
localToGlobal(itsBunch->R, phi, psi, meanR * Vector_t(0.001));
localToGlobal(itsBunch->P, phi, psi);
......
......@@ -868,6 +868,7 @@ void Distribution::DoRestartOpalCycl(PartBunch &beam, size_t Np, int restartStep
h5_int64_t localDump = 0;
previousH5Local_m = false;
rc = H5ReadStepAttribInt64(H5file, "LOCAL", &localDump);
if(rc != H5_SUCCESS) {
......@@ -879,8 +880,8 @@ void Distribution::DoRestartOpalCycl(PartBunch &beam, size_t Np, int restartStep
OPAL 1.3.0!");
} else {
if (localDump > 0) previousH5Local_m = true;
if (localDump == 1) previousH5Local_m = true;
}
double actualT;
......
......@@ -930,48 +930,52 @@ void DataSink::writePhaseSpace(PartBunch &beam, Vector_t FDext[], double sposHea
int DataSink::writePhaseSpace_cycl(PartBunch &beam, Vector_t FDext[], double meanEnergy,
double refPr, double refPt, double refPz,
double refR, double refTheta, double refZ,
double azimuth, double elevation, int localFrame) {
double azimuth, double elevation, bool local) {
if (!doHDF5_m) return -1;
//if (beam.getLocalNum() == 0) return -1; //TEMP for testing -DW
if (!doHDF5_m) return -1;
//if (beam.getLocalNum() == 0) return -1; //TEMP for testing -DW
h5_int64_t rc;
IpplTimings::startTimer(H5PartTimer_m);
h5_int64_t rc;
IpplTimings::startTimer(H5PartTimer_m);
beam.calcBeamParameters_cycl();
beam.calcBeamParameters_cycl();
double t = beam.getT();
double t = beam.getT();
Vektor< double, 3 > rmin = beam.get_origin();
Vektor< double, 3 > rmax = beam.get_maxExtend();
Vektor< double, 3 > centroid = beam.get_centroid();
size_t nLoc = beam.getLocalNum();
Vektor<double, 3 > xsigma = beam.get_rrms();
Vektor<double, 3 > psigma = beam.get_prms();
Vektor<double, 3 > geomvareps = beam.get_emit();
Vektor<double, 3 > vareps = beam.get_norm_emit();
Vektor< double, 3 > rmin = beam.get_origin();
Vektor< double, 3 > rmax = beam.get_maxExtend();
Vektor< double, 3 > centroid = beam.get_centroid();
size_t nLoc = beam.getLocalNum();
Vektor<double, 3 > xsigma = beam.get_rrms();
Vektor<double, 3 > psigma = beam.get_prms();
Vektor<double, 3 > geomvareps = beam.get_emit();
Vektor<double, 3 > vareps = beam.get_norm_emit();
Vektor<double, 3 > RefPartR = beam.RefPart_R;
Vektor<double, 3 > RefPartP = beam.RefPart_P;
Vektor<double, 3 > RefPartR = beam.RefPart_R;
Vektor<double, 3 > RefPartP = beam.RefPart_P;
double energySpread = beam.getdE();
double energySpread = beam.getdE();
double sigma = ((xsigma[0] * xsigma[0]) + (xsigma[1] * xsigma[1])) /
(2.0 * beam.get_gamma() * 17.0e3 * ((geomvareps[0] * geomvareps[0]) + (geomvareps[1] * geomvareps[1])));
double sigma = ((xsigma[0] * xsigma[0]) + (xsigma[1] * xsigma[1])) /
(2.0 * beam.get_gamma() * 17.0e3 * ((geomvareps[0] * geomvareps[0]) + (geomvareps[1] * geomvareps[1])));
Vektor< double, 3 > maxP(0.0);
Vektor< double, 3 > minP(0.0);
beam.get_PBounds(minP, maxP);
Vektor< double, 3 > maxP(0.0);
Vektor< double, 3 > minP(0.0);
beam.get_PBounds(minP, maxP);
std::unique_ptr<char[]> varray(new char[(nLoc)*sizeof(double)]);
double *farray = reinterpret_cast<double *>(varray.get());
h5_int64_t *larray = reinterpret_cast<h5_int64_t *>(varray.get());
std::unique_ptr<char[]> varray(new char[(nLoc)*sizeof(double)]);
double *farray = reinterpret_cast<double *>(varray.get());
h5_int64_t *larray = reinterpret_cast<h5_int64_t *>(varray.get());
double pathLength = beam.getLPath();
double pathLength = beam.getLPath();
h5_int64_t localTrackStep = (h5_int64_t)beam.getLocalTrackStep();
h5_int64_t globalTrackStep = (h5_int64_t)beam.getGlobalTrackStep();
h5_int64_t numBunch = (h5_int64_t)beam.getNumBunch();
h5_int64_t SteptoLastInj = (h5_int64_t)beam.getSteptoLastInj();
h5_int64_t localTrackStep = (h5_int64_t)beam.getLocalTrackStep();
h5_int64_t globalTrackStep = (h5_int64_t)beam.getGlobalTrackStep();
h5_int64_t numBunch = (h5_int64_t)beam.getNumBunch();
h5_int64_t SteptoLastInj = (h5_int64_t)beam.getSteptoLastInj();
h5_int64_t localFrame = 0;
if (local) localFrame = 1;
///Get the particle decomposition from all the compute nodes.
std::unique_ptr<size_t[]> locN(new size_t[Ippl::getNodes()]);
......@@ -1125,7 +1129,7 @@ int DataSink::writePhaseSpace_cycl(PartBunch &beam, Vector_t FDext[], double mea
if(rc != H5_SUCCESS)
ERRORMSG("H5 rc= " << rc << " in " << __FILE__ << " @ line " << __LINE__ << endl);
rc = H5WriteStepAttribInt64(H5file_m, "LOCAL", (h5_int64_t *)&localFrame, 1);
rc = H5WriteStepAttribInt64(H5file_m, "LOCAL", &localFrame, 1);
if(rc != H5_SUCCESS)
ERRORMSG("H5 rc= " << rc << " in " << __FILE__ << " @ line " << __LINE__ << endl);
......
......@@ -191,7 +191,7 @@ public:
int writePhaseSpace_cycl(PartBunch &beam, Vector_t FDext[], double E,
double refPr, double refPt, double refPz,
double refR, double refTheta, double refZ,
double azimuth, double elevation, int localFrame);
double azimuth, double elevation, bool local);
/** \brief Dumps Phase Space for Envelope trakcer to H5 file.
*
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment