Hey Eugene,
I looked at the macro (didnt run it) and have to admit I dont see
anything wrong with it...
Have you changed anything recently in this , assuming you have run it
before?
Have you tried piping in events from the flr to your ntuple without cuts
etc?
I suppose there is a problem either with the ntuple construction, the
Struct , SetBranch etc,
or with the cuts. I didnt see anything wrong with the former. I have
been using a similar macro
and the only problems I have ever had were related to messing up the
order of variables
(I think they have to be in the same order in the
struct Recon{
Float_t Ee,pe,thetae,phie,betae;
Float_t Ep,pp,thetap,phip,betap;
Float_t pm,Em,Mm,Q2;
Float_t pe_pp,pe_te,pe_tp;
Float_t pp_pe,pp_tp,pp_te;
Float_t te_pe,te_pp,te_tp;
Float_t tp_pp,tp_pe,tp_te;
Float_t Q2_pe,Q2_pp,Q2_te,Q2_tp,X2Q2,dQ2;
Float_t pe_Q2,pp_Q2,te_Q2,tp_Q2;
} recon; // calculated variables
AND in the Branch:
kc->Branch("recon",&recon,"Ee:pe:thetae:phie:betae:"
"Ep:pp:thetap:phip:betap:pm:Em:Mm:Q2:"
"pe_pp:pe_te:pe_tp:"
"pp_pe:pp_tp:pp_te:"
"te_pe:te_pp:te_tp:"
"tp_pp:tp_pe:tp_te:"
"Q2_pe:Q2_pp:Q2_te:Q2_tp:X2Q2:dQ2:"
"pe_Q2:pp_Q2:te_Q2:tp_Q2");
}
, and in this case they seem ok).
Another problem was related to the length of the Branch (somehow it
wouldnt write more than
a certain number of variables, have you tried using fewer variables in
your recon? There may be
a limit...).
Other than that, the problem could be in the cuts:
if ( LooseCut() ) {
if ( eLeftpRight() ){
...
}
if ( pLefteRight() ){
...
}
}
etc. I am a little suspicious of the double equality thing, which I have used before too but
seemed unreliable somehow...
emom = recon.pe = entry.pwr;
I would try without that too and see what happens.
I cant see anything obvious though.
n
Eugene J. Geis wrote:
>Hi,
>
>I'm pulling out my hair. I'm using the flr files and putting together
>a truncated dataset using simple methods such as new TTree() and
>SetAddress()... blah blah. There is one point in my macro where I
>simply make a new variable for electron momentum and, depending
>on cuts, read in the value, recon.pe = entry.pwr, where entry.pwr
>is linked to the flr entry pwr, etc.
> I'm then doing kinematic calculations and writing in OTHER entries
>in my branches of recon and entry. Somehow, I'm getting a split in
>the data down the middle of my momentum distribution that
>looks like the included graph c1.ps where I've plotted pe:pp (p_elec vs.
>p_prot). This "software mistake" is actually rewritten all the way down
>into the pwr entry that is linked directly to the address of the flr pwr
>entry. When I look at the actual data in the flr files, they have no split like
>this. I'm completely lost... I've rewritten the macro three times and
>cleared out a load of code and I have no freaking clue. If anyone has a
>few minutes to see what the F&%$ I'm doing wrong, or what weird C++/
>ROOT error I might have haphazardly invoked, PLEASE let me know...
>
>eugene
>
>--------------------------------------------------------------------------
>Eugene Geis
>PhD Student, Physics Department, ASU
>Research Affiliate, MIT-Bates Laboratory of Nuclear Science
>eugene.geis@asu.edu
>
>
>const Float_t eBeam = 0.8500; //electron beam energy
>
>const Int_t nSects = 2 ; //number of sectors
>const Int_t nSpins = 6 ; //number of beam-target spin combinations
>const Int_t nX2Bins = 20 ;
>
>const Float_t mE = 0.000511; // mass of electron
>const Float_t mD = 1.875613; //mass of deuteron
>const Float_t mP = 0.938272; //mass of proton
>const Float_t mN = 0.939565; //mass of neutron
>const Float_t Vsc = 30.00000; // speed of light in vacuum
>
>const Float_t deg2rad = (4.0*atan(1.0))/180.0; //[rad] = deg2rad * [deg]
>const Float_t rad2deg = 180.0/(4.0*atan(1.0)); //[deg] = rad2deg * [rad]
>
>//*****************************************************************************
>
>struct Entry { //linked to flr entries
> Float_t ttl, ttr, atl, atr, ntl, ntr, ptl, ptr;
> Float_t pwl, twl, fwl, zwl, qwl, xwl;
> Float_t pwr, twr, fwr, zwr, qwr, xwr;
> Float_t ncl, ncr, acl, acr, tcl, tcr;
> Float_t Bwl, Bwr, trig;
> Float_t nrun;
>} entry; // variables inhereted from reconstruction
>
> struct Recon{
> Float_t Ee,pe,thetae,phie,betae;
> Float_t Ep,pp,thetap,phip,betap;
> Float_t pm,Em,Mm,Q2;
> Float_t pe_pp,pe_te,pe_tp;
> Float_t pp_pe,pp_tp,pp_te;
> Float_t te_pe,te_pp,te_tp;
> Float_t tp_pp,tp_pe,tp_te;
> Float_t Q2_pe,Q2_pp,Q2_te,Q2_tp,X2Q2,dQ2;
> Float_t pe_Q2,pp_Q2,te_Q2,tp_Q2;
> } recon; // calculated variables
>
>TTree *kc = new TTree("kc", "KineCorr NTuple"); // usefull for debuging
>
>Bool_t mc = 0; //true if the data is Monte Carlo
>Bool_t whgen = 0; //true if the data is white-generated
>Bool_t ECUT = 0; //true if I do nice momentum vs. theta cut on electron
>Int_t nEvents=999999999;
>
>Float_t omega,q,phiq,thetaq,Q2_center,tester,tester2;
>Int_t sect,npass;
>Float_t pnx,pny,pnz;
>Float_t qx,qy,qz;
>Float_t pmx,pmy,pmz;
>Float_t theta,phi;
>Float_t pB[4] = {0};
>Float_t pQ[4] = {0};
>Float_t test[6];
>Float_t Q2pe, Q2pp, Q2te, Q2tp, Q2center, d_Q2; //deleted function calling these
>Float_t emom,pmom,eth,pth; // one attempt at debugging.
>//*****************************************************************************
>
>
>TFile *fmc;
>
>
>//*****************************************************************************
>
>void InitAll(){
>
> gStyle->SetOptStat(0);
> gStyle->SetOptTitle(0);
> gStyle->SetOptDate(1);
> gStyle->SetOptFit(1);
>
> int arc=gApplication->Argc();
> char** arv=gApplication->Argv();
>
>// for(Int_t i=0;i<arc;i++){
>// TString arg(arv[i]);
>// if(arg=="-m"){mc=true;cout<<"Monte Carlo File"<<endl;}
>// else if((arg=="-n")&&(i+1<arc)){sscanf(arv[++i],"%i",&nEvents);}
>// else if (arg=="-white"){whgen=true;cout<<"White Generator: "<<whgen<<endl;}
>// else if (arg=="-cut"){ECUT=true;}
>// }
>
> gROOT->ProcessLine(".x flr.C");
>
> TMatrix ch = charge[0]; //array holding various charges
>
>
> if(ch(0,0)> 0){ // beam target
> printf("Charge in (beam-target) polarization states:\n");
> printf(" ( h, Pz, Pzz) charge\n");
> printf(" -------- -----------\n");
> printf(" ( +1, +1, +1) %10.4f\n",ch(0,0));
> printf(" ( -1, +1, +1) %10.4f\n",ch(1,0));
> printf(" ( +1, -1, +1) %10.4f\n",ch(0,1));
> printf(" ( -1, -1, +1) %10.4f\n",ch(1,1));
> printf(" ( +1, 0, -2) %10.4f\n",ch(0,2));
> printf(" ( -1, 0, -2) %10.4f\n",ch(1,2));
> //printf(" ( -0, +0) %10.4f %10.4f\n",ch(0,2), ch(1,2));
> printf("\n");
> }
> else{
> printf("\nERROR: NO CHARGE\n\n");
> exit(1);
> }
>
>
> SetAddress(chain,&entry);
> kc->Branch("entry",&entry,"ttl:ttr:atl:atr:ntl:ntr:ptl:ptr:" // set TTree entries
> "pwl:twl:fwl:zwl:qwl:xwl:"
> "pwr:twr:fwr:zwr:qwr:xwr:"
> "ncl:ncr:acl:acr:tcl:tcr:"
> "Bwl:Bwr:trig:nrun");
> kc->Branch("recon",&recon,"Ee:pe:thetae:phie:betae:"
> "Ep:pp:thetap:phip:betap:pm:Em:Mm:Q2:"
> "pe_pp:pe_te:pe_tp:"
> "pp_pe:pp_tp:pp_te:"
> "te_pe:te_pp:te_tp:"
> "tp_pp:tp_pe:tp_te:"
> "Q2_pe:Q2_pp:Q2_te:Q2_tp:X2Q2:dQ2:"
> "pe_Q2:pp_Q2:te_Q2:tp_Q2");
>
>
Eugene J. Geis wrote:
>Hi,
>
>I'm pulling out my hair. I'm using the flr files and putting together
>a truncated dataset using simple methods such as new TTree() and
>SetAddress()... blah blah. There is one point in my macro where I
>simply make a new variable for electron momentum and, depending
>on cuts, read in the value, recon.pe = entry.pwr, where entry.pwr
>is linked to the flr entry pwr, etc.
> I'm then doing kinematic calculations and writing in OTHER entries
>in my branches of recon and entry. Somehow, I'm getting a split in
>the data down the middle of my momentum distribution that
>looks like the included graph c1.ps where I've plotted pe:pp (p_elec vs.
>p_prot). This "software mistake" is actually rewritten all the way down
>into the pwr entry that is linked directly to the address of the flr pwr
>entry. When I look at the actual data in the flr files, they have no split like
>this. I'm completely lost... I've rewritten the macro three times and
>cleared out a load of code and I have no freaking clue. If anyone has a
>few minutes to see what the F&%$ I'm doing wrong, or what weird C++/
>ROOT error I might have haphazardly invoked, PLEASE let me know...
>
>eugene
>
>--------------------------------------------------------------------------
>Eugene Geis
>PhD Student, Physics Department, ASU
>Research Affiliate, MIT-Bates Laboratory of Nuclear Science
>eugene.geis@asu.edu
>
>
>const Float_t eBeam = 0.8500; //electron beam energy
>
>const Int_t nSects = 2 ; //number of sectors
>const Int_t nSpins = 6 ; //number of beam-target spin combinations
>const Int_t nX2Bins = 20 ;
>
>const Float_t mE = 0.000511; // mass of electron
>const Float_t mD = 1.875613; //mass of deuteron
>const Float_t mP = 0.938272; //mass of proton
>const Float_t mN = 0.939565; //mass of neutron
>const Float_t Vsc = 30.00000; // speed of light in vacuum
>
>const Float_t deg2rad = (4.0*atan(1.0))/180.0; //[rad] = deg2rad * [deg]
>const Float_t rad2deg = 180.0/(4.0*atan(1.0)); //[deg] = rad2deg * [rad]
>
>//*****************************************************************************
>
>struct Entry { //linked to flr entries
> Float_t ttl, ttr, atl, atr, ntl, ntr, ptl, ptr;
> Float_t pwl, twl, fwl, zwl, qwl, xwl;
> Float_t pwr, twr, fwr, zwr, qwr, xwr;
> Float_t ncl, ncr, acl, acr, tcl, tcr;
> Float_t Bwl, Bwr, trig;
> Float_t nrun;
>} entry; // variables inhereted from reconstruction
>
> struct Recon{
> Float_t Ee,pe,thetae,phie,betae;
> Float_t Ep,pp,thetap,phip,betap;
> Float_t pm,Em,Mm,Q2;
> Float_t pe_pp,pe_te,pe_tp;
> Float_t pp_pe,pp_tp,pp_te;
> Float_t te_pe,te_pp,te_tp;
> Float_t tp_pp,tp_pe,tp_te;
> Float_t Q2_pe,Q2_pp,Q2_te,Q2_tp,X2Q2,dQ2;
> Float_t pe_Q2,pp_Q2,te_Q2,tp_Q2;
> } recon; // calculated variables
>
>TTree *kc = new TTree("kc", "KineCorr NTuple"); // usefull for debuging
>
>Bool_t mc = 0; //true if the data is Monte Carlo
>Bool_t whgen = 0; //true if the data is white-generated
>Bool_t ECUT = 0; //true if I do nice momentum vs. theta cut on electron
>Int_t nEvents=999999999;
>
>Float_t omega,q,phiq,thetaq,Q2_center,tester,tester2;
>Int_t sect,npass;
>Float_t pnx,pny,pnz;
>Float_t qx,qy,qz;
>Float_t pmx,pmy,pmz;
>Float_t theta,phi;
>Float_t pB[4] = {0};
>Float_t pQ[4] = {0};
>Float_t test[6];
>Float_t Q2pe, Q2pp, Q2te, Q2tp, Q2center, d_Q2; //deleted function calling these
>Float_t emom,pmom,eth,pth; // one attempt at debugging.
>//*****************************************************************************
>
>
>TFile *fmc;
>
>
>//*****************************************************************************
>
>void InitAll(){
>
> gStyle->SetOptStat(0);
> gStyle->SetOptTitle(0);
> gStyle->SetOptDate(1);
> gStyle->SetOptFit(1);
>
> int arc=gApplication->Argc();
> char** arv=gApplication->Argv();
>
>// for(Int_t i=0;i<arc;i++){
>// TString arg(arv[i]);
>// if(arg=="-m"){mc=true;cout<<"Monte Carlo File"<<endl;}
>// else if((arg=="-n")&&(i+1<arc)){sscanf(arv[++i],"%i",&nEvents);}
>// else if (arg=="-white"){whgen=true;cout<<"White Generator: "<<whgen<<endl;}
>// else if (arg=="-cut"){ECUT=true;}
>// }
>
> gROOT->ProcessLine(".x flr.C");
>
> TMatrix ch = charge[0]; //array holding various charges
>
>
> if(ch(0,0)> 0){ // beam target
> printf("Charge in (beam-target) polarization states:\n");
> printf(" ( h, Pz, Pzz) charge\n");
> printf(" -------- -----------\n");
> printf(" ( +1, +1, +1) %10.4f\n",ch(0,0));
> printf(" ( -1, +1, +1) %10.4f\n",ch(1,0));
> printf(" ( +1, -1, +1) %10.4f\n",ch(0,1));
> printf(" ( -1, -1, +1) %10.4f\n",ch(1,1));
> printf(" ( +1, 0, -2) %10.4f\n",ch(0,2));
> printf(" ( -1, 0, -2) %10.4f\n",ch(1,2));
> //printf(" ( -0, +0) %10.4f %10.4f\n",ch(0,2), ch(1,2));
> printf("\n");
> }
> else{
> printf("\nERROR: NO CHARGE\n\n");
> exit(1);
> }
>
>
> SetAddress(chain,&entry);
> kc->Branch("entry",&entry,"ttl:ttr:atl:atr:ntl:ntr:ptl:ptr:" // set TTree entries
> "pwl:twl:fwl:zwl:qwl:xwl:"
> "pwr:twr:fwr:zwr:qwr:xwr:"
> "ncl:ncr:acl:acr:tcl:tcr:"
> "Bwl:Bwr:trig:nrun");
> kc->Branch("recon",&recon,"Ee:pe:thetae:phie:betae:"
> "Ep:pp:thetap:phip:betap:pm:Em:Mm:Q2:"
> "pe_pp:pe_te:pe_tp:"
> "pp_pe:pp_tp:pp_te:"
> "te_pe:te_pp:te_tp:"
> "tp_pp:tp_pe:tp_te:"
> "Q2_pe:Q2_pp:Q2_te:Q2_tp:X2Q2:dQ2:"
> "pe_Q2:pp_Q2:te_Q2:tp_Q2");
>
>}
>
>int SetAddress(TChain** chain,struct Entry* entry) { // link to flr
> chain[0]->SetBranchStatus("ttl", 1);
> chain[0]->SetBranchStatus("ttr", 1);
> chain[0]->SetBranchStatus("atl", 1);
> chain[0]->SetBranchStatus("atr", 1);
> chain[0]->SetBranchStatus("ntl", 1);
> chain[0]->SetBranchStatus("ntr", 1);
> chain[0]->SetBranchStatus("ptl", 1);
> chain[0]->SetBranchStatus("ptr", 1);
>
> chain[0]->SetBranchStatus("pwl", 1);
> chain[0]->SetBranchStatus("twl", 1);
> chain[0]->SetBranchStatus("fwl", 1);
> chain[0]->SetBranchStatus("zwl", 1);
> chain[0]->SetBranchStatus("qwl", 1);
> chain[0]->SetBranchStatus("xwl", 1);
>
> chain[0]->SetBranchStatus("pwr", 1);
> chain[0]->SetBranchStatus("twr", 1);
> chain[0]->SetBranchStatus("fwr", 1);
> chain[0]->SetBranchStatus("zwr", 1);
> chain[0]->SetBranchStatus("qwr", 1);
> chain[0]->SetBranchStatus("xwr", 1);
>
> chain[0]->SetBranchStatus("ncl", 1);
> chain[0]->SetBranchStatus("ncr", 1);
> chain[0]->SetBranchStatus("acl", 1);
> chain[0]->SetBranchStatus("acr", 1);
> chain[0]->SetBranchStatus("tcl", 1);
> chain[0]->SetBranchStatus("tcr", 1);
>
> chain[0]->SetBranchStatus("Bwl", 1);
> chain[0]->SetBranchStatus("Bwr", 1);
>
> chain[0]->SetBranchStatus("trig", 1);
> chain[0]->SetBranchStatus("nrun", 1);
>
> chain[0]->SetBranchAddress("ttl", &entry->ttl);
> chain[0]->SetBranchAddress("ttr", &entry->ttr);
> chain[0]->SetBranchAddress("atl", &entry->atl);
> chain[0]->SetBranchAddress("atr", &entry->atr);
> chain[0]->SetBranchAddress("ntl", &entry->ntl);
> chain[0]->SetBranchAddress("ntr", &entry->ntr);
> chain[0]->SetBranchAddress("ptl", &entry->ptl);
> chain[0]->SetBranchAddress("ptr", &entry->ptr);
>
> chain[0]->SetBranchAddress("pwl", &entry->pwl);
> chain[0]->SetBranchAddress("twl", &entry->twl);
> chain[0]->SetBranchAddress("fwl", &entry->fwl);
> chain[0]->SetBranchAddress("zwl", &entry->zwl);
> chain[0]->SetBranchAddress("qwl", &entry->qwl);
> chain[0]->SetBranchAddress("xwl", &entry->xwl);
>
> chain[0]->SetBranchAddress("pwr", &entry->pwr);
> chain[0]->SetBranchAddress("twr", &entry->twr);
> chain[0]->SetBranchAddress("fwr", &entry->fwr);
> chain[0]->SetBranchAddress("zwr", &entry->zwr);
> chain[0]->SetBranchAddress("qwr", &entry->qwr);
> chain[0]->SetBranchAddress("xwr", &entry->xwr);
>
> chain[0]->SetBranchAddress("ncl", &entry->ncl);
> chain[0]->SetBranchAddress("ncr", &entry->ncr);
> chain[0]->SetBranchAddress("acl", &entry->acl);
> chain[0]->SetBranchAddress("acr", &entry->acr);
> chain[0]->SetBranchAddress("tcl", &entry->tcl);
> chain[0]->SetBranchAddress("tcr", &entry->tcr);
>
> chain[0]->SetBranchAddress("Bwl", &entry->Bwl);
> chain[0]->SetBranchAddress("Bwr", &entry->Bwr);
>
> chain[0]->SetBranchAddress("trig", &entry->trig);
> chain[0]->SetBranchAddress("nrun", &entry->nrun);
>
>}
>
>//*****************************************************************************
>
>bool LooseCut() {
> return (fabs(entry.zwl) < 20 && fabs(entry.zwr) < 20 && entry.trig==1);
>}
>
>bool eLeftpRight() {
> return (entry.qwl == -1 && entry.ntl >-1 && entry.ntr > -1 && entry.qwr == 1);
>}
>
>bool pLefteRight() {
> return (entry.qwr == -1 && entry.ntr >-1 && entry.ntl > -1 && entry.qwl == 1);
>}
>
>bool MmCut(Float_t tol){
> return (fabs(recon.Mm) < tol);
>}
>
>bool CopCut(Float_t tol){
> return (fabs(entry.fwl-entry.fwr+180.0) < tol);
>}
>
>Bool_t BetaCut(Float_t beta,Float_t min,Float_t max){
> return((beta>min)&&(beta<max));
>}
>
>//*****************************************************************************
>
>void Blast2Q(Float_t *pB, Float_t *pQ, Float_t theta, Float_t phi) { // rotation
>
> pQ[0] = pB[0]*(-cos(theta*deg2rad)*cos(phi*deg2rad))
> + pB[1]*(-cos(theta*deg2rad)*sin(phi*deg2rad))
> + pB[2]*sin(theta*deg2rad);
> pQ[1] = pB[0]*sin(phi*deg2rad) + pB[1]*(-cos(phi*deg2rad));
> pQ[2] = pB[0]*(sin(theta*deg2rad)*cos(phi*deg2rad))
> + pB[1]*(sin(theta*deg2rad)*sin(phi*deg2rad))
> + pB[2]*cos(theta*deg2rad);
> pQ[3] = pB[3];
>}
>
>//**************************************************************************
>void KineCalc() { /// calculate kinematics
> recon.Ee = sqrt(emom*emom+mE*mE);
>
> recon.Ep = sqrt(pmom*pmom+mP*mP);
>
> pB[0] = pnx = pmom*sin(pth*deg2rad)*cos(recon.phip*deg2rad);
> pB[1] = pny = pmom*sin(pth*deg2rad)*sin(recon.phip*deg2rad);
> pB[2] = pnz = pmom*cos(pth*deg2rad);
>
> qx = -emom*sin(eth*deg2rad)*cos(recon.phie*deg2rad);
> qy = -emom*sin(eth*deg2rad)*sin(recon.phie*deg2rad);
> qz = eBeam - emom*cos(eth*deg2rad);
> omega = eBeam-recon.Ee;
>
> q = sqrt(qx*qx+qy*qy+qz*qz);
> phiq = (abs(recon.phie) < 90) ? 180 + recon.phie : recon.phie - 180;
> thetaq = acos(qz/q)*rad2deg;
>
> Blast2Q(pB,pQ,thetaq,phiq);
>
> pmx = -pQ[0];
> pmy = -pQ[1];
> pmz = q - pQ[2];
> recon.pm = sqrt(q*q - 2*q*pQ[2] + pmom*pmom);
> recon.Em = (mP+eBeam-recon.Ee-recon.Ep);
> recon.Mm = sqrt(recon.Em*recon.Em-recon.pm*recon.pm);
> recon.Q2 = q*q-omega*omega;
>
> if (fabs(recon.Mm)<0.1) { // calculate elastic relations
>
> recon.pe_pp = mP+eBeam-sqrt(pmom**2+mP**2);
>
> recon.pe_te = (eBeam*mP+mE**2)/(mP+eBeam-eBeam*cos(eth*deg2rad));
>
> recon.pp_pe = sqrt((mP+eBeam-emom)**2-mP**2);
>
> recon.pp_tp = 2*mP*(eBeam+mP)/((eBeam*cos(pth*deg2rad))*(((eBeam+mP)/
>(eBeam*cos(pth*deg2rad)))**2-1));
>
> recon.pe_tp = mP+eBeam-sqrt(recon.pp_tp**2+mP**2);
>
> recon.pp_te = sqrt((mP+eBeam-recon.pe_te)**2-mP**2);
>
> test[0] = (emom*(eBeam+mP)-eBeam*mP-mE**2)/(eBeam*emom);
> test[1] = ((mP+eBeam-sqrt(pmom**2+mP**2))*(eBeam+mP)-eBeam*mP-mE**2)/
>(eBeam*(mP+eBeam-sqrt(pmom**2+mP**2)));
> test[2] = (recon.pe_tp*(eBeam+mP)-eBeam*mP-mE**2)/(eBeam*recon.pe_tp);
> test[3] = (eBeam+mP)*(sqrt(pmom**2+mP**2)-mP)/(eBeam*pmom);
> test[4] = (eBeam+mP)*(sqrt((sqrt((mP+eBeam-emom)**2-mP**2))**2+mP**2)-mP)/
>(eBeam*(sqrt((mP+eBeam-emom)**2-mP**2)));
> test[5] = (eBeam+mP)*(sqrt(recon.pp_te**2+mP**2)-mP)/(eBeam*recon.pp_te);
>
>
>
> if (fabs(test[0])<=1.0) recon.te_pe = rad2deg*acos(test[0]); // if acos() can be invoked
> else recon.te_pe = -1;
>
> if (fabs(test[1])<=1.0) recon.te_pp = rad2deg*acos(test[1]);
> else recon.te_pp = -1;
>
> if (fabs(test[2])<=1.0) recon.te_tp = rad2deg*acos(test[2]);
> else recon.te_tp = -1;
>
> if (fabs(test[3])<=1.0) recon.tp_pp = rad2deg*acos(test[3]);
> else recon.tp_pp = -1;
>
> if (fabs(test[4])<=1.0) recon.tp_pe = rad2deg*acos(test[4]);
> else recon.tp_pe = -1;
>
> if (fabs(test[5])<=1.0) recon.tp_te = rad2deg*acos(test[5]);
> else recon.tp_te = -1;
>
> if (recon.te_pe>0&&recon.te_pp>-1&&recon.te_tp>-1) { // if reasonable values are found
> recon.Q2_pe = 4*0.85*emom*(sin(recon.te_pe*deg2rad/2.0))**2;
> recon.Q2_pp = 4*0.85*recon.pe_pp*(sin(recon.te_pp*deg2rad/2.0))**2;
> recon.Q2_te = 4*0.85*recon.pe_te*(sin(eth*deg2rad/2.0))**2;
> recon.Q2_tp = 4*0.85*recon.pe_tp*(sin(recon.te_tp*deg2rad/2.0))**2;
> }
> kc->Fill(); // fill ntuple
> npass++;
> }
>}
>
>
>
>void show_corr() {
>
> InitAll(); // Initialize ntuples and data
> npass = 0;
> for(Int_t i=0;((i<chain[0]->GetEntries())&&(i<nEvents));++i){
> if(i%50000==0){printf("\n%9d/%6d/ ",i,npass);}
> else if(i%10000==0) cout<<":"<<flush;
> else if(i%1000==0) cout<<"."<<flush;
>
> sect = -1;
> chain[0]->GetEntry(i);
>
> if ( LooseCut() ) { // vertex
> if ( eLeftpRight() ){ // sector dependent ep-elastic
>
> emom = recon.pe = entry.pwl; // MOMENTUM with global variable attempt at debug
> eth = recon.thetae = entry.twl;
> recon.phie = entry.fwl;
> recon.betae = entry.Bwl;
> pmom = recon.pp = entry.pwr; // MOMENTUM
> pth = recon.thetap = entry.twr;
> recon.phip = entry.fwr;
> recon.betap = entry.Bwr;
>
> sect = 0;
> }
> if ( pLefteRight() ){
> emom = recon.pe = entry.pwr; //MOMENTUM
> eth = recon.thetae = entry.twr;
> recon.phie = entry.fwr;
> recon.betae = entry.Bwr;
> pmom = recon.pp = entry.pwl; //MOMENTUM
> pth = recon.thetap = entry.twl;
> recon.phip = entry.fwl;
> recon.betap = entry.Bwl;
>
> sect = 1;
>
> }
>
> if(sect>-1&&CopCut(5.0)&&recon.betae>0.99){ // IF passed cuts, calculate
> KineCalc(); // calculate kinematics
> }
>
> }
> }
> cout << "I'm calling output file."<<endl;
>
> fout = new TFile("/scratch/bud08/egeis/data/kc-1.root","recreate"); // WRITE FILE.
> cout << "I'm writing kc."<<endl;
> kc->Write();
> cout << "I'm closing the output file."<<endl;
> fout->Close();
>
>}
>
>
>
This archive was generated by hypermail 2.1.2 : Mon Feb 24 2014 - 14:07:32 EST