1 | #include "AMPTOOLS_DATAIO/etapiPlotGenerator.h" |
2 | #include "IUAmpTools/Histogram1D.h" |
3 | #include "IUAmpTools/Kinematics.h" |
4 | #include "TLorentzVector.h" |
5 | #include "TLorentzRotation.h" |
6 | |
7 | etapiPlotGenerator::etapiPlotGenerator( const FitResults& results ) : |
8 | PlotGenerator( results ) |
9 | { |
10 | bookHistogram( khm12, new Histogram1D( 400, 1.0, 5.0, "hm12", "Mass( 1 2 )" ) ); |
11 | bookHistogram( khm13, new Histogram1D( 400, 1.0, 5.0, "hm13", "Mass( 1 3 )" ) ); |
12 | bookHistogram( khm23, new Histogram1D( 100, 1.04, 2.56, "hm23", "M(#eta#pi) [GeV/c^{2}]" ) ); |
13 | bookHistogram( khm1, new Histogram1D( 200, 0.5, 1.5, "hm1", "Mass( 1 )" ) ); |
14 | bookHistogram( khm2, new Histogram1D( 200, 0.0, 1.0, "hm2", "Mass( 2 )" ) ); |
15 | bookHistogram( khm3, new Histogram1D( 200, 0.0, 1.0, "hm3", "Mass( 3 )" ) ); |
16 | bookHistogram( kdltz, new Histogram2D( 80, 0.0, 25.0, 80, 0.0, 9.0, "dltz", "Dalitz Plot" ) ); |
17 | bookHistogram( cosT, new Histogram1D( 220, -1.1, 1.1, "cosT", "CosTheta") ); |
18 | bookHistogram( phiAng, new Histogram1D(200, -3.2, 3.2, "phiAng", "#phi") ); |
19 | bookHistogram( cosT_lab, new Histogram1D( 220, -1.1, 1.1, "cosT_lab", "CosThetaLab") ); |
20 | bookHistogram( phiAng_lab, new Histogram1D(200, -3.2, 3.2, "phiAng_lab", "#phi_{lab}") ); |
21 | bookHistogram( cosT_m23_lab, new Histogram2D(200, 0.6, 2.6, 100, -1.0, 1.0, "cosTLab_m23", "cos(#theta_{lab}) vs. Mass(#eta #pi^{0})" ) ); |
22 | bookHistogram( phi_m23_lab, new Histogram2D(200, 0.6, 2.6, 100, -3.2, 3.2, "PhiLab_m23", "#phi_{lab} vs. Mass(#eta #pi^{0})" ) ); |
23 | bookHistogram( PhiT, new Histogram1D(200, -3.2, 3.2, "PhiT", "#Phi") ); |
24 | bookHistogram( Omega, new Histogram1D(100, -3.2, 3.2, "Omega", "#Omega") ); |
25 | bookHistogram( cosT_m23, new Histogram2D(200, 0.7, 2.7, 100, -1.0, 1.0, "cosT_m23", "cos(#theta) vs. Mass(#eta #pi^{0})" ) ); |
26 | bookHistogram( cosT_phi, new Histogram2D(40, -3.2, 3.2, 100, -1.0, 1.0, "cosT_phi", "cos(#theta) vs. #phi" ) ); |
27 | bookHistogram( cosT_Phi, new Histogram2D(40, -3.2, 3.2, 100, -1.0, 1.0, "cosT_Phi", "cos(#theta) vs. #Phi" ) ); |
28 | } |
29 | |
30 | void etapiPlotGenerator::projectEvent( Kinematics* kin ){ |
31 | |
32 | |
33 | |
34 | |
35 | projectEvent( kin, "" ); |
36 | } |
37 | |
38 | void etapiPlotGenerator::projectEvent( Kinematics* kin, const string& reactionName ){ |
39 | |
40 | |
41 | |
42 | |
43 | |
44 | int nargs = cfgInfo()->amplitudeList( reactionName, "", "" ).at(0)->factors().at(0).size(); |
45 | TVector3 eps; |
46 | TLorentzVector P0; |
47 | |
48 | if(nargs==7) { |
49 | double polAngle = stod(cfgInfo()->amplitudeList( reactionName, "", "" ).at(0)->factors().at(0).at(6)); |
| Value stored to 'polAngle' during its initialization is never read |
50 | polAngle = 0.0; |
51 | eps.SetXYZ(cos(polAngle*TMath::DegToRad()), sin(polAngle*TMath::DegToRad()), 0.0); |
52 | P0 = kin->particle(0); |
53 | } else if(nargs==5) { |
54 | P0.SetPxPyPzE(0., 0., kin->particle(0).E(),kin->particle(0).E()); |
55 | eps.SetXYZ(kin->particle(0).Px(),kin->particle(0).Py(), 0.); |
56 | } else { |
57 | std::cout << "Zlm does not have required number of arguments (4 OR 6, but found " << nargs << ")! Aborting." << std::endl; |
58 | std::cout << cfgInfo()->amplitudeList( reactionName, "", "" ).at(0)->factors().at(0).at(0) << std::endl; |
59 | return; |
60 | } |
61 | |
62 | TLorentzVector P1 = kin->particle(1); |
63 | TLorentzVector P2 = kin->particle(2); |
64 | TLorentzVector P3 = kin->particle(3); |
65 | |
66 | fillHistogram( khm12, (P1+P2).M() ); |
67 | fillHistogram( khm13, (P1+P3).M() ); |
68 | fillHistogram( khm23, (P2+P3).M() ); |
69 | fillHistogram( khm1, (P1).M() ); |
70 | fillHistogram( khm2, (P2).M() ); |
71 | fillHistogram( khm3, (P3).M() ); |
72 | fillHistogram( kdltz, (P1+P2).M2(), (P2+P3).M2() ); |
73 | |
74 | TLorentzVector resonance = P2 + P3; |
75 | TLorentzRotation resRestBoost( -resonance.BoostVector() ); |
76 | |
77 | TLorentzVector beam_res = resRestBoost * P0; |
78 | TLorentzVector recoil_res = resRestBoost * P1; |
79 | TLorentzVector p3_res = resRestBoost * P3; |
80 | |
81 | |
82 | GDouble locCosT_lab = resonance.CosTheta(); |
83 | GDouble locPhi_lab = resonance.Phi(); |
84 | |
85 | fillHistogram( cosT_lab, locCosT_lab ); |
86 | fillHistogram( phiAng_lab, locPhi_lab ); |
87 | fillHistogram( cosT_m23_lab, resonance.M(), locCosT_lab ); |
88 | fillHistogram( phi_m23_lab, resonance.M(), locPhi_lab ); |
89 | |
90 | |
91 | TVector3 z = -1. * recoil_res.Vect().Unit(); |
92 | |
93 | |
94 | |
95 | TVector3 y = (P0.Vect().Unit().Cross(-P1.Vect().Unit())).Unit(); |
96 | |
97 | TVector3 x = y.Cross(z); |
98 | |
99 | TVector3 angles( (p3_res.Vect()).Dot(x), |
100 | (p3_res.Vect()).Dot(y), |
101 | (p3_res.Vect()).Dot(z) ); |
102 | |
103 | Double_t cosTheta = angles.CosTheta(); |
104 | Double_t phi = angles.Phi(); |
105 | |
106 | GDouble Phi = atan2(y.Dot(eps), P0.Vect().Unit().Dot(eps.Cross(y))); |
107 | |
108 | TVector3 eta = P3.Vect(); |
109 | GDouble omega = atan2(y.Dot(eta), P0.Vect().Unit().Dot(eta.Cross(y))); |
110 | |
111 | |
112 | fillHistogram( PhiT, Phi); |
113 | fillHistogram( cosT, cosTheta); |
114 | fillHistogram( phiAng, phi); |
115 | fillHistogram( Omega, omega); |
116 | fillHistogram( cosT_m23, (P2+P3).M(), cosTheta); |
117 | fillHistogram( cosT_phi, phi, cosTheta); |
118 | fillHistogram( cosT_Phi, Phi, cosTheta); |
119 | } |