ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitHzz4l/Angles/src/ZGAngles.cc
Revision: 1.3
Committed: Mon Dec 17 12:35:22 2012 UTC (12 years, 5 months ago) by dkralph
Content type: text/plain
Branch: MAIN
CVS Tags: compiled, HEAD
Changes since 1.2: +0 -0 lines
Log Message:
many updates

File Contents

# Content
1 // ======================================================================================================
2 //
3 // angles following conventions in Gainer, et. al.
4 // arXiv:1112.1405v2 [hep-ph], 6 May 2012
5 // see also : Gainer, et. al.
6 // arVix:1108.2274v1 [hep-ph], 10 Aug 2011
7 //
8 // ======================================================================================================
9
10 #include "ZGAngles.h"
11 #include <iostream>
12 using namespace std;
13
14 #ifndef STANDALONE
15 #include "CommonDefs.h"
16
17 // --------------------------------------------------------------------------------------------------
18 //
19 // fill angle struct from ZZ-type GenInfo objects
20 //
21 void
22 fillZGAngles( GenInfoStruct &geninfo, ZGAngles &a, bool debug )
23 // --------------------------------------------------------------------------------------------------
24 {
25 ZGLabVectors l;
26
27 l.veckq = geninfo.x_1*TLorentzVector( 0,0,3500, 3500);
28 l.veckqbar = geninfo.x_2*TLorentzVector( 0,0,-3500, 3500);
29
30 if( geninfo.id_1_a < 0 ) {
31 l.veclp.SetPtEtaPhiM( geninfo.pt_1_a, geninfo.eta_1_a, geninfo.phi_1_a,
32 (abs(geninfo.id_1_a) == 11 ? ELECTRON_MASS : MUON_MASS) );
33 l.veclm.SetPtEtaPhiM( geninfo.pt_2_a, geninfo.eta_2_a, geninfo.phi_2_a,
34 (abs(geninfo.id_2_a) == 11 ? ELECTRON_MASS : MUON_MASS) );
35 } else {
36 l.veclp.SetPtEtaPhiM( geninfo.pt_2_a, geninfo.eta_2_a, geninfo.phi_2_a,
37 (abs(geninfo.id_2_a) == 11 ? ELECTRON_MASS : MUON_MASS) );
38 l.veclm.SetPtEtaPhiM( geninfo.pt_1_a, geninfo.eta_1_a, geninfo.phi_1_a,
39 (abs(geninfo.id_1_a) == 11 ? ELECTRON_MASS : MUON_MASS) );
40 }
41
42
43 // l.vecz = (l.veclp+l.veclm);
44 l.vecg.SetPtEtaPhiM( geninfo.vpt_b, geninfo.veta_b, geninfo.vphi_b, geninfo.vmass_b);
45 l.vecz.SetPtEtaPhiM( geninfo.vpt_a, geninfo.veta_a, geninfo.vphi_a, geninfo.vmass_a);
46 l.veczg = (l.vecz+l.vecg);
47
48 getZGAngles( l, a, debug );
49 }
50
51 // --------------------------------------------------------------------------------------------------
52 //
53 // fill angle struct from ZZ-type EventData
54 //
55 void
56 fillZGAngles( EventData &data, ZGAngles &a )
57 // --------------------------------------------------------------------------------------------------
58 {
59 ZGLabVectors l;
60
61 if( data.Z1leptons[0].charge > 0 ) {
62 l.veclp = data.Z1leptons[0].vec;
63 l.veclm = data.Z1leptons[1].vec;
64 l.vecg = data.Z2leptons[0].vec;
65 } else {
66 l.veclp = data.Z1leptons[1].vec;
67 l.veclm = data.Z1leptons[0].vec;
68 l.vecg = data.Z2leptons[0].vec;
69 }
70
71 l.vecz = (l.veclp+l.veclm);
72 l.veczg = (l.vecz+l.vecg);
73
74 getZGAngles( l, a );
75 }
76
77 #endif
78
79
80
81 // --------------------------------------------------------------------------------------------------
82 //
83 // set angles from a struct of lab-frame 4vecs
84 //
85 void
86 getZGAngles( ZGLabVectors &l, ZGAngles &a, bool debug )
87 // --------------------------------------------------------------------------------------------------
88 {
89 if( debug ) {
90 cout << endl << endl;
91 cout << "veczg(x:y:z) :: " << l.veczg.X() <<":"<< l.veczg.Y() <<":"<< l.veczg.Z() << endl;
92 cout << "\tpt: " << l.veczg.Pt() << "\teta: " << l.veczg.Eta() << "\tphi:" << l.veczg.Phi() << endl;
93 cout << "vecz(x:y:z) :: " << l.vecz.X() <<":"<< l.vecz.Y() <<":"<< l.vecz.Z() << endl;
94 cout << "\tpt: " << l.vecz.Pt() << "\teta: " << l.vecz.Eta() << "\tphi:" << l.vecz.Phi() << endl;
95 cout << "veclp(x:y:z) :: " << l.veclp.X() <<":"<< l.veclp.Y() <<":"<< l.veclp.Z() << endl;
96 cout << "\tpt: " << l.veclp.Pt() << "\teta: " << l.veclp.Eta() << "\tphi:" << l.veclp.Phi() << endl;
97 cout << "veclm(x:y:z) :: " << l.veclm.X() <<":"<< l.veclm.Y() <<":"<< l.veclm.Z() << endl;
98 cout << "\tpt: " << l.veclm.Pt() << "\teta: " << l.veclm.Eta() << "\tphi:" << l.veclm.Phi() << endl;
99 cout << "vecg(x:y:z) :: " << l.vecg.X() <<":"<< l.vecg.Y() <<":"<< l.vecg.Z() << endl;
100 cout << "\tpt: " << l.vecg.Pt() << "\teta: " << l.vecg.Eta() << "\tphi:" << l.vecg.Phi() << endl;
101 }
102
103
104
105 //
106 // define the frames & coordinate systems
107 // --------------------------------------------------------------------------------------------------
108
109 TVector3 Xframe = l.veczg.BoostVector();
110 TVector3 Z1frame = l.vecz.BoostVector();
111
112 if( debug ) {
113 cout << "Xboost(x:y:z) :: " << Xframe.X() <<":"<< Xframe.Y() <<":"<< Xframe.Z() << endl;
114 cout << "Zboost(x:y:z) :: " << Z1frame.X() <<":"<< Z1frame.Y() <<":"<< Z1frame.Z() << endl;
115 }
116
117 // "partons"
118 TLorentzVector kq( 0, 0, (l.veczg.E()+l.veczg.Pz())/2, (l.veczg.E()+l.veczg.Pz())/2 );
119 TLorentzVector kqbar( 0, 0, (l.veczg.Pz()-l.veczg.E())/2, (l.veczg.E()-l.veczg.Pz())/2 );
120 TLorentzVector veckq_in_Xframe(kq);
121 TLorentzVector veckqbar_in_Xframe(kqbar);
122 veckq_in_Xframe.Boost(-1*Xframe);
123 veckqbar_in_Xframe.Boost(-1*Xframe);
124 if( debug ) {
125 cout << "kq_X(x:y:z) :: " << veckq_in_Xframe.X() <<":"<< veckq_in_Xframe.Y() <<":"<< veckq_in_Xframe.Z() << endl;
126 cout << "kqb_X(x:y:z) :: " << veckqbar_in_Xframe.X() <<":"<< veckqbar_in_Xframe.Y() <<":"<< veckqbar_in_Xframe.Z() << endl;
127 }
128
129 // Z vectors
130 if( debug ) { cout << "boosting Z,g ..." << endl;}
131 TLorentzVector vecz_in_Xframe (l.vecz);
132 TLorentzVector vecg_in_Xframe (l.vecg);
133 TLorentzVector vecz_in_Z1frame (l.vecz);
134 vecz_in_Xframe.Boost(-1*Xframe);
135 vecg_in_Xframe.Boost(-1*Xframe);
136 vecz_in_Z1frame.Boost(-1*Z1frame);
137
138 if( debug ) { cout << "rotating ..." << endl;}
139 // coord system in the CM frame
140 TVector3 uz_in_Xframe = vecz_in_Xframe.Vect().Unit();
141 TVector3 uy_in_Xframe = (veckq_in_Xframe.Vect().Unit().Cross(uz_in_Xframe.Unit())).Unit();
142 TVector3 ux_in_Xframe = (uy_in_Xframe.Unit().Cross(uz_in_Xframe.Unit())).Unit();
143 TRotation rotation;
144 rotation = rotation.RotateAxes( ux_in_Xframe,uy_in_Xframe,uz_in_Xframe ).Inverse();
145
146 //
147 // for going to the Z frames from the CM frame,
148 // boost after transform
149 //
150 if( debug ) { cout << "xform z ..." << endl; }
151 TLorentzVector vecz_in_Xframe_newcoords(vecz_in_Xframe);
152 vecz_in_Xframe_newcoords.Transform(rotation);
153 TVector3 Z1frame_from_Xframe_newcoords = vecz_in_Xframe_newcoords.BoostVector();
154 if( debug ) {
155 cout << "vecz_in_Xframe_newcoords (x:y:z) :: "
156 << vecz_in_Xframe_newcoords.X() <<":"
157 << vecz_in_Xframe_newcoords.Y() <<":"
158 << vecz_in_Xframe_newcoords.Z() << endl;
159
160 cout << "Z1frame_from_Xframe_newcoords (x:y:z) :: "
161 << Z1frame_from_Xframe_newcoords.X() <<":"
162 << Z1frame_from_Xframe_newcoords.Y() <<":"
163 << Z1frame_from_Xframe_newcoords.Z() << endl;
164 }
165
166
167 //
168 // theta(lm), phi(lm) in Z1 frame
169 // --------------------------------------------------------------------------------------------------
170 TLorentzVector veclm_in_Z1frame(l.veclm);
171 TLorentzVector veclp_in_Z1frame(l.veclp);
172
173 // first boost to CM, then redefine coords
174 if( debug ) { cout << "xform lm..." << endl; }
175 veclm_in_Z1frame.Boost(-1*Xframe);
176 veclm_in_Z1frame.Transform(rotation);
177 veclp_in_Z1frame.Boost(-1*Xframe);
178 veclp_in_Z1frame.Transform(rotation);
179
180 // then boost to Z1
181 if ( debug ) { cout << "boost lm ..." << endl;}
182 veclm_in_Z1frame.Boost(-1*Z1frame_from_Xframe_newcoords);
183 veclp_in_Z1frame.Boost(-1*Z1frame_from_Xframe_newcoords);
184
185 // now get angles
186 a.phi = veclm_in_Z1frame.Phi();
187 a.costheta_lm = veclm_in_Z1frame.CosTheta();
188 a.costheta_lp = veclp_in_Z1frame.CosTheta();
189
190 //
191 // or, just go directly to Z frame, skip CM frame
192 //
193 /*
194 // direct
195 TLorentzVector veclm_in_Z1frame_direct(l.veclm);
196 veclm_in_Z1frame_direct.Boost(-1*Z1frame);
197 TLorentzVector veclp_in_Z1frame_direct(l.veclp);
198 veclp_in_Z1frame_direct.Boost(-1*Z1frame);
199 a.costheta_lm = veclm_in_Z1frame_direct.Vect().Unit().Dot((-1*l.vecz.Vect()).Unit());
200 a.costheta_lp = veclp_in_Z1frame_direct.Vect().Unit().Dot((-1*l.vecz.Vect()).Unit());
201 */
202
203 if( debug ) {
204 TVector3 unitz(0,0,1);
205 cout << "compare :: ct: " << a.costheta_lm
206 << "\tv.ct: " << veclm_in_Z1frame.Vect().CosTheta()
207 << "\tc(t): " << TMath::Cos(veclm_in_Z1frame.Theta())
208 << "\tc(v.t): " << TMath::Cos(veclm_in_Z1frame.Theta())
209 << "\tdot: " << veclm_in_Z1frame.Vect().Unit().Dot(unitz)
210 << endl;
211 cout << "costheta lp vs lm : " << veclp_in_Z1frame.CosTheta() << " vs " << veclm_in_Z1frame.CosTheta() << endl;
212 }
213
214
215 //
216 // Theta(Z1,Zg-lab) in X frame
217 // --------------------------------------------------------------------------------------------------
218 TLorentzVector veczg_in_Xframe(l.veczg);
219 // veczg_in_Xframe.Boost(-1*Xframe);
220 veczg_in_Xframe.Transform(rotation);
221
222 if( debug ) {
223 cout << "veczg_X (x:y:z) :: " << veczg_in_Xframe.X()
224 <<":"<< veczg_in_Xframe.Y()
225 <<":"<< veczg_in_Xframe.Z()
226 << endl;
227 cout << "vecz_X (x:y:z) :: " << vecz_in_Xframe.X()
228 <<":"<< vecz_in_Xframe.Y()
229 <<":"<< vecz_in_Xframe.Z()
230 << endl;
231 cout << "vecg_X (x:y:z) :: " << vecg_in_Xframe.X()
232 <<":"<< vecg_in_Xframe.Y()
233 <<":"<< vecg_in_Xframe.Z()
234 << endl;
235
236 // check
237 TLorentzVector ztmp(l.vecz);
238 ztmp.Boost(-1*Xframe);
239 ztmp.Transform(rotation);
240
241 TLorentzVector gtmp(l.vecg);
242 gtmp.Boost(-1*Xframe);
243 gtmp.Transform(rotation);
244
245 cout << "ztmp.CosTheta(): " << ztmp.CosTheta() << endl;
246 cout << "gtmp.CosTheta(): " << gtmp.CosTheta() << endl;
247 }
248
249
250 TLorentzVector veczg_in_Xframe_newcoords(l.veczg);
251 // veczg_in_Xframe_newcoords.Boost(-1*Xframe);
252 veczg_in_Xframe_newcoords.Transform(rotation);
253 a.cosTheta = (-1*veczg_in_Xframe_newcoords.Vect()).CosTheta();
254
255 if( debug ) { cout << "phi: " << a.phi << "\tct: " << a.costheta_lm << "\tcT:" << a.cosTheta << endl;}
256
257 a.ptg = l.vecg.Pt();
258 a.etag = l.vecg.Eta();
259 a.ptl1 = l.veclp.Pt();
260 a.etal1 = l.veclp.Eta();
261 a.ptl2 = l.veclm.Pt();
262 a.etal2 = l.veclm.Eta();
263 a.ptz = (l.veclp+l.veclm).Pt();
264 a.etaz = (l.veclp+l.veclm).Eta();
265 a.mzg = l.veczg.M();
266 a.mz = l.vecz.M();
267 };