1 | |
2 | |
3 | |
4 | |
5 | |
6 | |
7 | |
8 | #include "JEventProcessor_ST_Tresolution.h" |
9 | #include "TRIGGER/DTrigger.h" |
10 | using namespace jana; |
11 | |
12 | |
13 | |
14 | #include <JANA/JApplication.h> |
15 | #include <JANA/JFactory.h> |
16 | extern "C"{ |
17 | void InitPlugin(JApplication *app){ |
18 | InitJANAPlugin(app); |
19 | app->AddProcessor(new JEventProcessor_ST_Tresolution()); |
20 | } |
21 | } |
22 | |
23 | |
24 | |
25 | |
26 | |
27 | JEventProcessor_ST_Tresolution::JEventProcessor_ST_Tresolution() |
28 | { |
29 | |
30 | } |
31 | |
32 | |
33 | |
34 | |
35 | JEventProcessor_ST_Tresolution::~JEventProcessor_ST_Tresolution() |
36 | { |
37 | |
38 | } |
39 | |
40 | |
41 | |
42 | |
43 | jerror_t JEventProcessor_ST_Tresolution::init(void) |
44 | { |
45 | |
46 | |
47 | |
48 | |
49 | |
50 | |
51 | |
52 | |
53 | |
54 | |
55 | |
56 | TDirectory *main = gDirectory(TDirectory::CurrentDirectory()); |
57 | gDirectory(TDirectory::CurrentDirectory())->mkdir("ST_Tresolution")->cd(); |
58 | |
59 | h2_CorrectedTime_z = new TH2I*[NCHANNELS]; |
60 | |
61 | int NoBins_time = 200; |
62 | int NoBins_z = 60; |
63 | double time_lower_limit = -10.0; |
64 | double time_upper_limit = 10.0; |
65 | double z_lower_limit = 0.0; |
66 | double z_upper_limit = 60.0; |
67 | for (Int_t i = 0; i < NCHANNELS; i++) |
68 | { |
69 | h2_CorrectedTime_z[i] = new TH2I(Form("h2_CorrectedTime_z_%i", i+1), "Corrected Time vs. Path length along the paddle; Path length along the paddle (cm); Propagation Time (ns)", NoBins_z,z_lower_limit,z_upper_limit, NoBins_time, time_lower_limit, time_upper_limit); |
70 | } |
71 | |
72 | |
73 | main->cd(); |
74 | |
75 | return NOERROR; |
76 | } |
77 | |
78 | |
79 | |
80 | |
81 | jerror_t JEventProcessor_ST_Tresolution::brun(JEventLoop *eventLoop, int32_t runnumber) |
82 | { |
83 | |
84 | |
85 | |
86 | vector<double> locRFPeriodVector; |
87 | eventLoop->GetCalib("PHOTON_BEAM/RF/rf_period", locRFPeriodVector); |
88 | dRFBunchPeriod = locRFPeriodVector[0]; |
89 | |
90 | |
91 | map<string,double> target_params; |
92 | if (eventLoop->GetCalib("/TARGET/target_parms", target_params)) |
| |
93 | jout << "Error loading /TARGET/target_parms/ !" << endl; |
94 | if (target_params.find("TARGET_Z_POSITION") != target_params.end()) |
| |
95 | z_target_center = target_params["TARGET_Z_POSITION"]; |
96 | else |
97 | jerr << "Unable to get TARGET_Z_POSITION from /TARGET/target_parms !" << endl; |
98 | |
99 | |
100 | DApplication* dapp = dynamic_cast<DApplication*>(eventLoop->GetJApplication()); |
| 3 | | Calling 'JEventLoop::GetJApplication' | |
|
| 4 | | Returning from 'JEventLoop::GetJApplication' | |
|
| |
101 | if(!dapp) |
| |
| |
102 | _DBG_std::cerr<<"plugins/Calibration/ST_Tresolution/JEventProcessor_ST_Tresolution.cc" <<":"<<102<<" "<<"Cannot get DApplication from JEventLoop! (are you using a JApplication based program?)"<<endl; |
103 | DGeometry* locGeometry = dapp->GetDGeometry(eventLoop->GetJEvent().GetRunNumber()); |
| 8 | | Called C++ object pointer is null |
|
104 | sc_angle_corr = 1.; |
105 | if(locGeometry->GetStartCounterGeom(sc_pos, sc_norm)) { |
106 | double theta = sc_norm[0][sc_norm[0].size()-2].Theta(); |
107 | sc_angle_corr = 1./cos(M_PI_21.57079632679489661923 - theta); |
108 | } |
109 | |
110 | |
111 | if(eventLoop->GetCalib("START_COUNTER/propagation_time_corr", propagation_time_corr)) |
112 | jout << "Error loading /START_COUNTER/propagation_time_corr !" << endl; |
113 | |
114 | |
115 | if(eventLoop->GetCalib("START_COUNTER/PTC_Boundary", PTC_Boundary)) |
116 | jout << "Error loading /START_COUNTER/PTC_Boundary !" << endl; |
117 | |
118 | |
119 | trackingFOMCut = 0.0027; |
120 | |
121 | return NOERROR; |
122 | } |
123 | |
124 | |
125 | |
126 | |
127 | jerror_t JEventProcessor_ST_Tresolution::evnt(JEventLoop *loop, uint64_t eventnumber) |
128 | { |
129 | |
130 | |
131 | |
132 | |
133 | |
134 | |
135 | |
136 | |
137 | |
138 | |
139 | |
140 | |
141 | |
142 | |
143 | |
144 | |
145 | const DTrigger* locTrigger = NULL__null; |
146 | loop->GetSingle(locTrigger); |
147 | if(locTrigger->Get_L1FrontPanelTriggerBits() != 0) |
148 | return NOERROR; |
149 | |
150 | double speed_light = 29.9792458; |
151 | |
152 | vector<const DSCHit *> scHitVector; |
153 | loop->Get(scHitVector); |
154 | |
155 | |
156 | const DRFTime* thisRFTime = NULL__null; |
157 | vector <const DRFTime*> RFTimeVector; |
158 | auto dRFTimeFactory = static_cast<DRFTime_factory*>(loop->Get(RFTimeVector)); |
159 | if (RFTimeVector.size() != 0) |
160 | thisRFTime = RFTimeVector[0]; |
161 | |
162 | |
163 | vector<const DChargedTrack *> chargedTrackVector; |
164 | loop->Get(chargedTrackVector); |
165 | |
166 | |
167 | const DDetectorMatches* locDetectorMatches = NULL__null; |
168 | loop->GetSingle(locDetectorMatches); |
169 | |
170 | |
171 | const DEventRFBunch *thisRFBunch = NULL__null; |
172 | loop->GetSingle(thisRFBunch); |
173 | |
174 | |
175 | const DParticleID *dParticleID = NULL__null; |
176 | loop->GetSingle(dParticleID); |
177 | |
178 | for (uint32_t i = 0; i < chargedTrackVector.size(); i++) |
179 | { |
180 | |
181 | const DChargedTrack *thisChargedTrack = chargedTrackVector[i]; |
182 | |
183 | const DTrackTimeBased *timeBasedTrack = thisChargedTrack->Get_BestTrackingFOM()->Get_TrackTimeBased(); |
184 | |
185 | |
186 | |
187 | |
188 | if(timeBasedTrack->FOM < trackingFOMCut) continue; |
189 | |
190 | |
191 | shared_ptr<const DSCHitMatchParams> locBestSCHitMatchParams; |
192 | bool foundSC = dParticleID->Get_BestSCMatchParams(timeBasedTrack, locDetectorMatches, locBestSCHitMatchParams); |
193 | if (!foundSC) continue; |
194 | |
195 | |
196 | DVector3 vertex = timeBasedTrack->position(); |
197 | double z_v = vertex.z(); |
198 | double r_v = vertex.Perp(); |
199 | |
200 | bool z_vertex_cut = fabs(z_target_center - z_v) <= 15.0; |
201 | bool r_vertex_cut = r_v < 0.5; |
202 | |
203 | if (!z_vertex_cut) continue; |
204 | if (!r_vertex_cut) continue; |
205 | vector<shared_ptr<const DSCHitMatchParams>> st_params; |
206 | bool st_match = locDetectorMatches->Get_SCMatchParams(timeBasedTrack, st_params); |
207 | |
208 | if (!st_match) continue; |
209 | |
210 | DVector3 IntersectionPoint, IntersectionMomentum; |
211 | shared_ptr<DSCHitMatchParams> locSCHitMatchParams; |
212 | vector<DTrackFitter::Extrapolation_t>extrapolations=timeBasedTrack->extrapolations.at(SYS_START); |
213 | bool sc_match_pid = dParticleID->Cut_MatchDistance(extrapolations, st_params[0]->dSCHit, st_params[0]->dSCHit->t, locSCHitMatchParams, true, &IntersectionPoint, &IntersectionMomentum); |
214 | |
215 | if(!sc_match_pid) continue; |
216 | |
217 | if (thisRFBunch->dNumParticleVotes < 2) continue; |
218 | |
219 | |
220 | |
221 | double locCenteredRFTime = thisRFTime->dTime; |
222 | |
223 | double locCenterToVertexRFTime = (timeBasedTrack->z() - z_target_center)*(1.0/speed_light); |
224 | double locVertexRFTime = locCenteredRFTime + locCenterToVertexRFTime; |
225 | int sc_index= locSCHitMatchParams->dSCHit->sector - 1; |
226 | |
227 | double sc_pos_soss = sc_pos[sc_index][0].z(); |
228 | double sc_pos_eoss = sc_pos[sc_index][1].z(); |
229 | double sc_pos_eobs = sc_pos[sc_index][11].z(); |
230 | double sc_pos_eons = sc_pos[sc_index][12].z(); |
231 | |
232 | double st_time = st_params[0]->dSCHit->t; |
233 | |
234 | double FlightTime = locSCHitMatchParams->dFlightTime; |
235 | |
236 | double st_corr_FlightTime = st_time - FlightTime; |
237 | |
238 | |
239 | double locSCzIntersection = IntersectionPoint.z(); |
240 | |
241 | |
242 | |
243 | |
244 | double incpt_ss = propagation_time_corr[sc_index][0]; |
245 | double slope_ss = propagation_time_corr[sc_index][1]; |
246 | double incpt_bs = propagation_time_corr[sc_index][2]; |
247 | double slope_bs = propagation_time_corr[sc_index][3]; |
248 | double incpt_ns = propagation_time_corr[sc_index][4]; |
249 | double slope_ns = propagation_time_corr[sc_index][5]; |
250 | |
251 | double Bound1 = PTC_Boundary[0][0]; |
252 | double Bound2 = PTC_Boundary[1][0]; |
253 | |
254 | |
255 | |
256 | |
257 | |
258 | |
259 | |
260 | |
261 | double SS_Length = sc_pos_eoss - sc_pos_soss; |
262 | |
263 | |
264 | |
265 | |
266 | |
267 | japp->RootFillLock(this); |
268 | |
269 | |
270 | if (sc_pos_soss < locSCzIntersection && locSCzIntersection <= sc_pos_eoss) |
271 | { |
272 | double path_ss = locSCzIntersection - sc_pos_soss; |
273 | double Corr_Time_ss = st_corr_FlightTime - (incpt_ss + (slope_ss * path_ss)); |
274 | double SC_RFShiftedTime = dRFTimeFactory->Step_TimeToNearInputTime(locVertexRFTime, Corr_Time_ss); |
275 | h2_CorrectedTime_z[sc_index]->Fill(path_ss, Corr_Time_ss -SC_RFShiftedTime); |
276 | } |
277 | |
278 | if(sc_pos_eoss < locSCzIntersection && locSCzIntersection <= sc_pos_eobs) |
279 | { |
280 | |
281 | double path_bs = SS_Length + (locSCzIntersection - sc_pos_eoss)*sc_angle_corr; |
282 | double Corr_Time_bs = st_corr_FlightTime - (incpt_ss + (slope_ss * path_bs)); |
283 | double SC_RFShiftedTime = dRFTimeFactory->Step_TimeToNearInputTime(locVertexRFTime, Corr_Time_bs); |
284 | h2_CorrectedTime_z[sc_index]->Fill(path_bs,Corr_Time_bs - SC_RFShiftedTime); |
285 | } |
286 | |
287 | if(sc_pos_eobs < locSCzIntersection && locSCzIntersection <= sc_pos_eons) |
288 | { |
289 | |
290 | double path_ns = SS_Length + (locSCzIntersection - sc_pos_eoss)*sc_angle_corr; |
291 | if (path_ns <= Bound1) |
292 | { |
293 | double Corr_Time_ns = st_corr_FlightTime - (incpt_ss + (slope_ss * path_ns)); |
294 | double SC_RFShiftedTime = dRFTimeFactory->Step_TimeToNearInputTime(locVertexRFTime, Corr_Time_ns); |
295 | h2_CorrectedTime_z[sc_index]->Fill(path_ns,Corr_Time_ns - SC_RFShiftedTime); |
296 | } |
297 | else if ((Bound1 < path_ns)&&(path_ns <= Bound2)) |
298 | { |
299 | double Corr_Time_ns = st_corr_FlightTime - (incpt_bs + (slope_bs * path_ns)); |
300 | double SC_RFShiftedTime = dRFTimeFactory->Step_TimeToNearInputTime(locVertexRFTime, Corr_Time_ns); |
301 | h2_CorrectedTime_z[sc_index]->Fill(path_ns,Corr_Time_ns - SC_RFShiftedTime); |
302 | } |
303 | else |
304 | { |
305 | double Corr_Time_ns = st_corr_FlightTime - (incpt_ns + (slope_ns * path_ns)); |
306 | double SC_RFShiftedTime = dRFTimeFactory->Step_TimeToNearInputTime(locVertexRFTime, Corr_Time_ns); |
307 | h2_CorrectedTime_z[sc_index]->Fill(path_ns,Corr_Time_ns - SC_RFShiftedTime); |
308 | } |
309 | |
310 | } |
311 | japp->RootFillUnLock(this); |
312 | |
313 | } |
314 | |
315 | |
316 | return NOERROR; |
317 | } |
318 | |
319 | |
320 | |
321 | |
322 | jerror_t JEventProcessor_ST_Tresolution::erun(void) |
323 | { |
324 | |
325 | |
326 | |
327 | return NOERROR; |
328 | } |
329 | |
330 | |
331 | |
332 | |
333 | jerror_t JEventProcessor_ST_Tresolution::fini(void) |
334 | { |
335 | |
336 | return NOERROR; |
337 | } |
338 | |