Commit 721e061e authored by Erik Asbjørn Mikkelsen Jensen's avatar Erik Asbjørn Mikkelsen Jensen
Browse files

deleted expand-nuchart-QBxn-QBxp-QBa.cxx

parent dc08aa52
// TODO: Finish this script (remember to check whether the resulting values correct?),
// same as expand-nuchart-QEC.cxx
/*
* In this ROOT macro the ROOT Tree 'nuchart.root' generated by treemaker.py is expanded by
* adding Q-values for electron capture to the tree via the formula
* QEC = MEP - MED;
* QEC is the Q-value for electron capture,
* MEP is the mass excess of the Parent which is decaying,
* MED is the mass excess of the Daughter which the Parent decays into.
* If either or both of MEP or MED is an estimate (as given by the original ROOT Tree
* 'nuchart.root'), QEC is also considered an estimate.
* If MEP exists, but there is no available MED, the Q-value for the Parent is flagged as
* incalculable.
*/
#define ELMASS 510.99895 // electron mass in keV
#include <TFile.h>
#include <TTree.h>
using namespace std;
int main() {
// Open tree file and retrieve the tree named 'a'
TFile f("nuchart.root", "update");
TTree *a;
f.GetObject("a", a);
// Locate the tabulated values of mass excesses which are already defined in the tree
struct QEC_params { Double_t QEC; Bool_t QEC_est; Bool_t QEC_calc; };
QEC_params QEC = QEC_params();
a->SetBranchAddress("QEC", &QEC);
struct Sp_params { Double_t Sp; Bool_t Sp_est; Bool_t Sp_calc; };
Sp_params Sp = Sp_params();
a->SetBranchAddress("Sp", &Sp);
// Add a new branch to the tree which is to contain QEC values and the above-mentioned boolean values
struct QBp_params { Double_t QBp; Bool_t QBp_est; Bool_t QBp_calc; };
QBp_params qbp = QBp_params();
TBranch *newBranch = a->Branch("QBp", &qbp.QBp, "QBp/D:QBp_est/O:QBp_calc");
a->BuildIndex("A", "Z"); // prerequisite for searching in the tree based on A and Z values in GetEntryNumberWithIndex() below
Double_t QBP = 0.0, SpD = 0.0;
Double_t QBp = 0.0;
Bool_t QBp_est = 0, QBp_calc = 0;
Long64_t indP = -1, indD = -1;
Long64_t minA = a->GetMinimum("A"), maxA = a->GetMaximum("A"), minZ = a->GetMinimum("Z"), maxZ = a->GetMaximum("Z");
for (Long64_t A = minA; A < maxA; A++) {
for (Long64_t Z = minZ; Z < maxZ; Z++) {
indP = a->GetEntryNumberWithIndex(A, Z); // returns '-1' if the entry does not exist
indD = a->GetEntryNumberWithIndex(A, Z - 1);
if (indP != -1 && indD != -1) { // we can get the mass excesses of both the Parent and the Daughter nuclei
QBp_calc = 1;
a->GetEntry(indP);
if (QEC.QEC_est == 1) {
QBp_est = 1;
}
QBP = QEC.QEC - 2*ELMASS;
a->GetEntry(indD);
if (Sp.Sp_est == 1) {
QBp_est = 1;
} else { // if neither of the ME's were estimates (ME_est = 0) then QEC is not an estimate (QEC_est = 0)
QBp_est = 0;
}
SpD = Sp.Sp;
QBp = QBP - SpD;
qbp.QBp = QBp; // values might not be recorded correctly, unless the values are assigned in this order
qbp.QBp_est = QBp_est;
qbp.QBp_calc = QBp_calc;
} else if (indP != -1) { // only the mass excess of the Parent is known; record the fact that we cannot compute QEC
QBp = NAN; // mimicking assignment order and appearance from just above
QBp_est = 0;
QBp_calc = 0;
qbp.QBp = QBp;
qbp.QBp_est = QBp_est;
qbp.QBp_calc = QBp_calc;
} else { // no entry in the tree corresponds to this case
continue;
}
a->GetEntry(indP); // navigate (back) to the index of the Parent nucleus ...
newBranch->Fill(); // ... and record the QEC value and the boolean values there
}
}
a->Write("", TObject::kOverwrite); // save only the new version of the tree
f.Close();
return EXIT_SUCCESS;
}
Supports Markdown
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