Bug Summary

File:plugins/Calibration/ST_Tresolution/JEventProcessor_ST_Tresolution.cc
Location:line 103, column 28
Description:Called C++ object pointer is null

Annotated Source Code

1// $Id$
2//
3// File: JEventProcessor_ST_Tresolution.cc
4// Created: Fri Jan 8 09:07:34 EST 2016
5// Creator: mkamel (on Linux ifarm1401 2.6.32-431.el6.x86_64 x86_64)
6//
7
8#include "JEventProcessor_ST_Tresolution.h"
9#include "TRIGGER/DTrigger.h"
10using namespace jana;
11
12
13// Routine used to create our JEventProcessor
14#include <JANA/JApplication.h>
15#include <JANA/JFactory.h>
16extern "C"{
17void InitPlugin(JApplication *app){
18 InitJANAPlugin(app);
19 app->AddProcessor(new JEventProcessor_ST_Tresolution());
20}
21} // "C"
22
23
24//------------------
25// JEventProcessor_ST_Tresolution (Constructor)
26//------------------
27JEventProcessor_ST_Tresolution::JEventProcessor_ST_Tresolution()
28{
29
30}
31
32//------------------
33// ~JEventProcessor_ST_Tresolution (Destructor)
34//------------------
35JEventProcessor_ST_Tresolution::~JEventProcessor_ST_Tresolution()
36{
37
38}
39
40//------------------
41// init
42//------------------
43jerror_t JEventProcessor_ST_Tresolution::init(void)
44{
45 // This is called once at program startup. If you are creating
46 // and filling historgrams in this plugin, you should lock the
47 // ROOT mutex like this:
48 //
49 // japp->RootWriteLock();
50 // ... fill historgrams or trees ...
51 // japp->RootUnLock();
52 //
53 // **************** define histograms *************************
54
55 //Create root folder and cd to it, store main dir
56 TDirectory *main = gDirectory(ROOT::Internal::TDirectoryAtomicAdapter(TDirectory::CurrentDirectory
()));
57 gDirectory(ROOT::Internal::TDirectoryAtomicAdapter(TDirectory::CurrentDirectory
()))->mkdir("ST_Tresolution")->cd();
58
59 h2_CorrectedTime_z = new TH2I*[NCHANNELS];
60 // All my Calculations in 2015 were using the binning below
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 // cd back to main directory
73 main->cd();
74
75 return NOERROR;
76}
77
78//------------------
79// brun
80//------------------
81jerror_t JEventProcessor_ST_Tresolution::brun(JEventLoop *eventLoop, int32_t runnumber)
82{
83 // This is called whenever the run number changes
84
85 //RF Period
86 vector<double> locRFPeriodVector;
87 eventLoop->GetCalib("PHOTON_BEAM/RF/rf_period", locRFPeriodVector);
88 dRFBunchPeriod = locRFPeriodVector[0];
89
90 // Obtain the target center along z;
91 map<string,double> target_params;
92 if (eventLoop->GetCalib("/TARGET/target_parms", target_params))
1
Taking true branch
→
93 jout << "Error loading /TARGET/target_parms/ !" << endl;
94 if (target_params.find("TARGET_Z_POSITION") != target_params.end())
2
←
Taking false branch
→
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 // Obtain the Start Counter geometry
100 DApplication* dapp = dynamic_cast<DApplication*>(eventLoop->GetJApplication());
3
←
Calling 'JEventLoop::GetJApplication'
→
4
←
Returning from 'JEventLoop::GetJApplication'
→
5
←
'dapp' initialized here
→
101 if(!dapp)
6
←
Assuming 'dapp' is null
→
7
←
Taking true branch
→
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 // Propagation Time constant
111 if(eventLoop->GetCalib("START_COUNTER/propagation_time_corr", propagation_time_corr))
112 jout << "Error loading /START_COUNTER/propagation_time_corr !" << endl;
113
114 // Propagation Time fit Boundaries
115 if(eventLoop->GetCalib("START_COUNTER/PTC_Boundary", PTC_Boundary))
116 jout << "Error loading /START_COUNTER/PTC_Boundary !" << endl;
117
118 // set some parameters
119 trackingFOMCut = 0.0027; // 3 sigma cut
120
121 return NOERROR;
122}
123
124//------------------
125// evnt
126//------------------
127jerror_t JEventProcessor_ST_Tresolution::evnt(JEventLoop *loop, uint64_t eventnumber)
128{
129 // This is called for every event. Use of common resources like writing
130 // to a file or filling a histogram should be mutex protected. Using
131 // loop->Get(...) to get reconstructed objects (and thereby activating the
132 // reconstruction algorithm) should be done outside of any mutex lock
133 // since multiple threads may call this method at the same time.
134 // Here's an example:
135 //
136 // vector<const MyDataClass*> mydataclasses;
137 // loop->Get(mydataclasses);
138 //
139 // japp->RootWriteLock();
140 // ... fill historgrams or trees ...
141 // japp->RootUnLock();
142
143
144 // select events with physics events, i.e., not LED and other front panel triggers
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 // SC hits
152 vector<const DSCHit *> scHitVector;
153 loop->Get(scHitVector);
154
155 // RF time object (and factory)
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 // Grab charged tracks
163 vector<const DChargedTrack *> chargedTrackVector;
164 loop->Get(chargedTrackVector);
165
166 // Grab the associated detector matches object
167 const DDetectorMatches* locDetectorMatches = NULL__null;
168 loop->GetSingle(locDetectorMatches);
169
170 // Grab the associated RF bunch object
171 const DEventRFBunch *thisRFBunch = NULL__null;
172 loop->GetSingle(thisRFBunch);
173
174 // Grab DParticleID object
175 const DParticleID *dParticleID = NULL__null;
176 loop->GetSingle(dParticleID);
177
178 for (uint32_t i = 0; i < chargedTrackVector.size(); i++)
179 {
180 // Grab the charged track and declare time based track object
181 const DChargedTrack *thisChargedTrack = chargedTrackVector[i];
182 // Grab associated time based track object by selecting charged track with best FOM
183 const DTrackTimeBased *timeBasedTrack = thisChargedTrack->Get_BestTrackingFOM()->Get_TrackTimeBased();
184
185 // Implement quality cuts for the time based tracks
186 //trackingFOMCut = 0.0027; // 3 sigma cut
187 //trackingFOMCut = 0.0001; // 5 sigma cut
188 if(timeBasedTrack->FOM < trackingFOMCut) continue;
189
190 // Grab the ST hit match params object and cut on only tracks matched to the ST
191 shared_ptr<const DSCHitMatchParams> locBestSCHitMatchParams;
192 bool foundSC = dParticleID->Get_BestSCMatchParams(timeBasedTrack, locDetectorMatches, locBestSCHitMatchParams);
193 if (!foundSC) continue;
194
195 // Define vertex vector and cut on target/scattering chamber geometry
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 // Apply vertex cut
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 // If st_match = true, there is a match between this track and the ST
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 // Cut on the number of particle votes to find the best RF time
217 if (thisRFBunch->dNumParticleVotes < 2) continue;
218 // Calculate the TOF estimate of the target time
219
220 // Calculate the RF estimate of the target time
221 double locCenteredRFTime = thisRFTime->dTime;
222 // RF time at center of target
223 double locCenterToVertexRFTime = (timeBasedTrack->z() - z_target_center)*(1.0/speed_light); // Time correction for photon from target center to vertex of track
224 double locVertexRFTime = locCenteredRFTime + locCenterToVertexRFTime;
225 int sc_index= locSCHitMatchParams->dSCHit->sector - 1;
226 // Start Counter geometry in hall coordinates
227 double sc_pos_soss = sc_pos[sc_index][0].z(); // Start of straight section
228 double sc_pos_eoss = sc_pos[sc_index][1].z(); // End of straight section
229 double sc_pos_eobs = sc_pos[sc_index][11].z(); // End of bend section
230 double sc_pos_eons = sc_pos[sc_index][12].z(); // End of nose section
231 //Get the ST time walk corrected time
232 double st_time = st_params[0]->dSCHit->t;
233 // Get the Flight time
234 double FlightTime = locSCHitMatchParams->dFlightTime;
235 //St time corrected for the flight time
236 double st_corr_FlightTime = st_time - FlightTime;
237 // SC_RFShiftedTime = dRFTimeFactory->Step_TimeToNearInputTime(locVertexRFTime, st_corr_FlightTime);
238 // Z intersection of charged track and SC
239 double locSCzIntersection = IntersectionPoint.z();
240 ////////////////////////////////////////////////////////////////////
241 /// Fill the sc time histograms corrected for walk and propagation//
242 ////////////////////////////////////////////////////////////////////
243 // Read the constants from CCDB
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 //Read fit boundary from CCDB
251 double Bound1 = PTC_Boundary[0][0];
252 double Bound2 = PTC_Boundary[1][0];
253 //cout << "Bound1 = " << Bound1 << endl;
254 //cout << "Bound2 = " << Bound2 << endl;
255 ///////////////////////////////////////
256 //Calculate the path along the paddle
257 /////////////////////////////////////
258 //Define some parameters
259 //double Radius = 12.0;
260 //double theta = 18.5 * pi/180.0;
261 double SS_Length = sc_pos_eoss - sc_pos_soss;// same for along z or along the paddle
262 //double BS_Length = Radius * theta ; // along the paddle
263 //double NS_Length = (sc_pos_eons - sc_pos_eobs)/cos(theta);// along the paddle
264
265 // FILL HISTOGRAMS
266 // Since we are filling histograms local to this plugin, it will not interfere with other ROOT operations: can use plugin-wide ROOT fill lock
267 japp->RootFillLock(this); //ACQUIRE ROOT FILL LOCK
268
269 // Straight Sections
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 // Bend Sections
278 if(sc_pos_eoss < locSCzIntersection && locSCzIntersection <= sc_pos_eobs)
279 {
280 //double path_bs = SS_Length + Radius * asin((locSCzIntersection - sc_pos_eoss)/Radius);
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 // Nose Sections
287 if(sc_pos_eobs < locSCzIntersection && locSCzIntersection <= sc_pos_eons)
288 {
289 //double path_ns = SS_Length + BS_Length +((locSCzIntersection - sc_pos_eobs)/cos(theta));
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); //RELEASE ROOT FILL LOCK
312
313 } // sc charged tracks
314
315
316 return NOERROR;
317}
318
319//------------------
320// erun
321//------------------
322jerror_t JEventProcessor_ST_Tresolution::erun(void)
323{
324 // This is called whenever the run number changes, before it is
325 // changed to give you a chance to clean up before processing
326 // events from the next run number.
327 return NOERROR;
328}
329
330//------------------
331// fini
332//------------------
333jerror_t JEventProcessor_ST_Tresolution::fini(void)
334{
335 // Called before program exit after event processing is finished.
336 return NOERROR;
337}
338