Bug Summary

File:scratch/gluex/scan-build-work/hdgeant4/src/GlueXSensitiveDetectorUPV.cc
Location:line 63, column 27
Description:Value stored to 'jcalib' during its initialization is never read

Annotated Source Code

1//
2// GlueXSensitiveDetectorUPV - class implementation
3//
4// author: richard.t.jones at uconn.edu
5// version: october 26, 2016
6
7#include "GlueXSensitiveDetectorUPV.hh"
8#include "GlueXDetectorConstruction.hh"
9#include "GlueXPrimaryGeneratorAction.hh"
10#include "GlueXUserEventInformation.hh"
11#include "GlueXUserTrackInformation.hh"
12#include "HddmOutput.hh"
13
14#include "G4VPhysicalVolume.hh"
15#include "G4PVPlacement.hh"
16#include "G4EventManager.hh"
17#include "G4HCofThisEvent.hh"
18#include "G4Step.hh"
19#include "G4SDManager.hh"
20#include "G4ios.hh"
21
22#include <JANA/JApplication.h>
23
24// Cutoff on the number of allowed hits per counter
25int GlueXSensitiveDetectorUPV::MAX_HITS = 100;
26
27// Light propagation parameters in tof bars
28double GlueXSensitiveDetectorUPV::ATTENUATION_LENGTH = 150.*cm;
29double GlueXSensitiveDetectorUPV::C_EFFECTIVE = 19*cm/ns;
30
31// Minimum hit time difference for two hits on the same bar
32double GlueXSensitiveDetectorUPV::TWO_HIT_TIME_RESOL = 50*ns;
33
34// Minimum energy deposition for a hit
35double GlueXSensitiveDetectorUPV::THRESH_MEV = 5.;
36
37int GlueXSensitiveDetectorUPV::instanceCount = 0;
38G4Mutex GlueXSensitiveDetectorUPV::fMutex = G4MUTEX_INITIALIZER{ { 0, 0, 0, 0, 0, 0, 0, { 0, 0 } } };
39
40GlueXSensitiveDetectorUPV::GlueXSensitiveDetectorUPV(const G4String& name)
41 : G4VSensitiveDetector(name),
42 fBarHitsMap(0), fPointsMap(0)
43{
44 collectionName.insert("UPVBarHitsCollection");
45 collectionName.insert("UPVPointsCollection");
46
47 // The rest of this only needs to happen once, the first time an object
48 // of this type is instantiated for this configuration of geometry and
49 // fields. If the geometry or fields change in such a way as to modify
50 // the drift-time properties of hits in the UPV, you must delete all old
51 // objects of this class and create new ones.
52
53 G4AutoLock barrier(&fMutex);
54 if (instanceCount++ == 0) {
55 int runno = HddmOutput::getRunNo();
56 extern jana::JApplication *japp;
57 if (japp == 0) {
58 G4cerr(*G4cerr_p) << "Error in GlueXSensitiveDetector constructor - "
59 << "jana global DApplication object not set, "
60 << "cannot continue." << G4endlstd::endl;
61 exit(-1);
62 }
63 jana::JCalibration *jcalib = japp->GetJCalibration(runno);
Value stored to 'jcalib' during its initialization is never read
64 if (japp == 0) { // dummy
65 jcalib = 0;
66 G4cout(*G4cout_p) << "UPV: ALL parameters loaded from ccdb" << G4endlstd::endl;
67 }
68 }
69}
70
71GlueXSensitiveDetectorUPV::GlueXSensitiveDetectorUPV(
72 const GlueXSensitiveDetectorUPV &src)
73 : G4VSensitiveDetector(src),
74 fBarHitsMap(src.fBarHitsMap), fPointsMap(src.fPointsMap)
75{
76 G4AutoLock barrier(&fMutex);
77 ++instanceCount;
78}
79
80GlueXSensitiveDetectorUPV &GlueXSensitiveDetectorUPV::operator=(const
81 GlueXSensitiveDetectorUPV &src)
82{
83 G4AutoLock barrier(&fMutex);
84 *(G4VSensitiveDetector*)this = src;
85 fBarHitsMap = src.fBarHitsMap;
86 fPointsMap = src.fPointsMap;
87 return *this;
88}
89
90GlueXSensitiveDetectorUPV::~GlueXSensitiveDetectorUPV()
91{
92 G4AutoLock barrier(&fMutex);
93 --instanceCount;
94}
95
96void GlueXSensitiveDetectorUPV::Initialize(G4HCofThisEvent* hce)
97{
98 fBarHitsMap = new
99 GlueXHitsMapUPVbar(SensitiveDetectorName, collectionName[0]);
100 fPointsMap = new
101 GlueXHitsMapUPVpoint(SensitiveDetectorName, collectionName[1]);
102 G4SDManager *sdm = G4SDManager::GetSDMpointer();
103 hce->AddHitsCollection(sdm->GetCollectionID(collectionName[0]), fBarHitsMap);
104 hce->AddHitsCollection(sdm->GetCollectionID(collectionName[1]), fPointsMap);
105}
106
107G4bool GlueXSensitiveDetectorUPV::ProcessHits(G4Step* step,
108 G4TouchableHistory* ROhist)
109{
110 double dEsum = step->GetTotalEnergyDeposit();
111 if (dEsum == 0)
112 return false;
113
114 const G4ThreeVector &pin = step->GetPreStepPoint()->GetMomentum();
115 const G4ThreeVector &xin = step->GetPreStepPoint()->GetPosition();
116 const G4ThreeVector &xout = step->GetPostStepPoint()->GetPosition();
117 double Ein = step->GetPreStepPoint()->GetTotalEnergy();
118 double tin = step->GetPreStepPoint()->GetGlobalTime();
119 double tout = step->GetPostStepPoint()->GetGlobalTime();
120 G4ThreeVector x = (xin + xout) / 2;
121 G4ThreeVector dx = xout - xin;
122 double t = (tin + tout) / 2;
123
124 const G4VTouchable* touch = step->GetPreStepPoint()->GetTouchable();
125 const G4AffineTransform &local_from_global = touch->GetHistory()
126 ->GetTopTransform();
127 G4ThreeVector xlocal = local_from_global.TransformPoint(x);
128
129 // For particles that range out inside the active volume, the
130 // "out" time may sometimes be set to something enormously high.
131 // This screws up the hit. Check for this case here by looking
132 // at tout and making sure it is less than 1 second. If it's
133 // not, then just use tin for "t".
134
135 if (tout > 1.0*s)
136 t = tin;
137
138 // Post the hit to the points list in the
139 // order of appearance in the event simulation.
140
141 int layer = GetIdent("layer", touch);
142 int row = GetIdent("row", touch);
143 G4Track *track = step->GetTrack();
144 G4int trackID = track->GetTrackID();
145 GlueXUserTrackInformation *trackinfo = (GlueXUserTrackInformation*)
146 track->GetUserInformation();
147 if (trackinfo->GetGlueXHistory() == 0 && Ein/MeV > THRESH_MEV) {
148 G4int key = fPointsMap->entries();
149 GlueXHitUPVpoint* lastPoint = (*fPointsMap)[key - 1];
150 if (lastPoint == 0 || lastPoint->track_ != trackID ||
151 fabs(lastPoint->t_ns - t/ns) > 0.1 ||
152 fabs(lastPoint->x_cm - x[0]/cm) > 2. ||
153 fabs(lastPoint->y_cm - x[1]/cm) > 2. ||
154 fabs(lastPoint->z_cm - x[2]/cm) > 2.)
155 {
156 int pdgtype = track->GetDynamicParticle()->GetPDGcode();
157 int g3type = GlueXPrimaryGeneratorAction::ConvertPdgToGeant3(pdgtype);
158 GlueXHitUPVpoint newPoint;
159 newPoint.ptype_G3 = g3type;
160 newPoint.track_ = trackID;
161 newPoint.trackID_ = trackinfo->GetGlueXTrackID();
162 newPoint.primary_ = (track->GetParentID() == 0);
163 newPoint.t_ns = t/ns;
164 newPoint.x_cm = x[0]/cm;
165 newPoint.y_cm = x[1]/cm;
166 newPoint.z_cm = x[2]/cm;
167 newPoint.px_GeV = pin[0]/GeV;
168 newPoint.py_GeV = pin[1]/GeV;
169 newPoint.pz_GeV = pin[2]/GeV;
170 newPoint.E_GeV = Ein/GeV;
171 fPointsMap->add(key, newPoint);
172 }
173 }
174
175 // Post the hit to the hits map, ordered by layer,row,end index
176
177 if (dEsum > 0) {
178 int key = GlueXHitUPVbar::GetKey(layer, row);
179 GlueXHitUPVbar *counter = (*fBarHitsMap)[key];
180 if (counter == 0) {
181 GlueXHitUPVbar newcounter(layer, row);
182 fBarHitsMap->add(key, newcounter);
183 counter = (*fBarHitsMap)[key];
184 }
185
186 double dxleft = xlocal[0];
187 double dxright = -xlocal[0];
188
189 // Calculate time at the PMT "normalized" to the center,
190 // so a hit in the center will have time "t" at both PMTs.
191 double tleft = t + dxleft / C_EFFECTIVE;
192 double tright = t + dxright / C_EFFECTIVE;
193
194 // calculate energy seen by PM for this track step using attenuation factor
195 double dEleft = dEsum * exp(-dxleft / ATTENUATION_LENGTH);
196 double dEright = dEsum * exp(-dxright / ATTENUATION_LENGTH);
197
198 // Add the hit to the hits vector, maintaining strict time ordering
199
200 if (dEleft > 0) {
201 // add the hit on end=0 (north/left end of the bar)
202 int merge_hit = 0;
203 std::vector<GlueXHitUPVbar::hitinfo_t>::iterator hiter;
204 for (hiter = counter->hits.begin(); hiter != counter->hits.end(); ++hiter) {
205 if (hiter->end_ == 0) {
206 if (fabs(hiter->t_ns*ns - tleft) < TWO_HIT_TIME_RESOL) {
207 merge_hit = 1;
208 break;
209 }
210 else if (hiter->t_ns*ns > tleft) {
211 break;
212 }
213 }
214 }
215 if (merge_hit) {
216 // Use the time from the earlier hit but add the charge
217 hiter->E_GeV += dEleft/GeV;
218 if (hiter->t_ns*ns > tleft) {
219 hiter->t_ns = tleft/ns;
220 }
221 }
222 else {
223 // create new hit
224 hiter = counter->hits.insert(hiter, GlueXHitUPVbar::hitinfo_t());
225 hiter->end_ = 0;
226 hiter->E_GeV = dEleft/GeV;
227 hiter->t_ns = tleft/ns;
228 }
229 }
230
231 if (dEright/MeV > 0) {
232 // add the hit on end=1 (south/bottom end of the bar)
233 int merge_hit = 0;
234 std::vector<GlueXHitUPVbar::hitinfo_t>::iterator hiter;
235 for (hiter = counter->hits.begin(); hiter != counter->hits.end(); ++hiter) {
236 if (hiter->end_ == 1) {
237 if (fabs(hiter->t_ns*ns - tright) < TWO_HIT_TIME_RESOL) {
238 merge_hit = 1;
239 break;
240 }
241 else if (hiter->t_ns*ns > tright) {
242 break;
243 }
244 }
245 }
246 if (merge_hit) {
247 // Use the time from the earlier hit but add the charge
248 hiter->E_GeV += dEright/GeV;
249 if (hiter->t_ns*ns > tright) {
250 hiter->t_ns = tright/ns;
251 }
252 }
253 else {
254 // create new hit
255 hiter = counter->hits.insert(hiter, GlueXHitUPVbar::hitinfo_t());
256 hiter->end_ = 1;
257 hiter->E_GeV = dEright/GeV;
258 hiter->t_ns = tright/ns;
259 }
260 }
261 }
262 return true;
263}
264
265void GlueXSensitiveDetectorUPV::EndOfEvent(G4HCofThisEvent*)
266{
267 std::map<int,GlueXHitUPVbar*> *bars = fBarHitsMap->GetMap();
268 std::map<int,GlueXHitUPVpoint*> *points = fPointsMap->GetMap();
269 if (bars->size() == 0 && points->size() == 0)
270 return;
271 std::map<int,GlueXHitUPVbar*>::iterator siter;
272 std::map<int,GlueXHitUPVpoint*>::iterator piter;
273
274 if (verboseLevel > 1) {
275 G4cout(*G4cout_p) << G4endlstd::endl
276 << "--------> Hits Collection: in this event there are "
277 << bars->size() << " bars with hits in the UPV: "
278 << G4endlstd::endl;
279 for (siter = bars->begin(); siter != bars->end(); ++siter)
280 siter->second->Print();
281
282 G4cout(*G4cout_p) << G4endlstd::endl
283 << "--------> Hits Collection: in this event there are "
284 << points->size() << " truth points in the UPV: "
285 << G4endlstd::endl;
286 for (piter = points->begin(); piter != points->end(); ++piter)
287 piter->second->Print();
288 }
289
290 // pack hits into ouptut hddm record
291
292 G4EventManager* mgr = G4EventManager::GetEventManager();
293 G4VUserEventInformation* info = mgr->GetUserInformation();
294 hddm_s::HDDM *record = ((GlueXUserEventInformation*)info)->getOutputRecord();
295 if (record == 0) {
296 G4cerr(*G4cerr_p) << "GlueXSensitiveDetectorUPV::EndOfEvent error - "
297 << "hits seen but no output hddm record to save them into, "
298 << "cannot continue!" << G4endlstd::endl;
299 exit(1);
300 }
301
302 if (record->getPhysicsEvents().size() == 0)
303 record->addPhysicsEvents();
304 if (record->getHitViews().size() == 0)
305 record->getPhysicsEvent().addHitViews();
306 hddm_s::HitView &hitview = record->getPhysicsEvent().getHitView();
307 if (hitview.getUpstreamEMvetos().size() == 0)
308 hitview.addUpstreamEMvetos();
309 hddm_s::UpstreamEMveto &upv = hitview.getUpstreamEMveto();
310
311 // Collect and output the tofTruthHits
312 for (siter = bars->begin(); siter != bars->end(); ++siter) {
313 std::vector<GlueXHitUPVbar::hitinfo_t> &hits = siter->second->hits;
314 // apply a pulse height threshold cut
315 for (unsigned int ih=0; ih < hits.size(); ++ih) {
316 if (hits[ih].E_GeV*1e3 <= THRESH_MEV) {
317 hits.erase(hits.begin() + ih);
318 --ih;
319 }
320 }
321 if (hits.size() > 0) {
322 hddm_s::UpvPaddleList counter = upv.addUpvPaddles(1);
323 counter(0).setLayer(siter->second->layer_);
324 counter(0).setRow(siter->second->row_);
325 // first the end=0 hits
326 for (int ih=0; ih < (int)hits.size(); ++ih) {
327 if (hits[ih].end_ == 0) {
328 hddm_s::UpvTruthHitList thit = counter(0).addUpvTruthHits(1);
329 thit(0).setEnd(hits[ih].end_);
330 thit(0).setE(hits[ih].E_GeV);
331 thit(0).setT(hits[ih].t_ns);
332 }
333 }
334 // followed by the end=1 hits
335 for (int ih=0; ih < (int)hits.size(); ++ih) {
336 if (hits[ih].end_ == 1) {
337 hddm_s::UpvTruthHitList thit = counter(0).addUpvTruthHits(1);
338 thit(0).setEnd(hits[ih].end_);
339 thit(0).setE(hits[ih].E_GeV);
340 thit(0).setT(hits[ih].t_ns);
341 }
342 }
343 int hitscount = counter(0).getUpvTruthHits().size();
344 if (hitscount > 2 * MAX_HITS) {
345 counter(0).deleteUpvTruthHits(-1, 2 * MAX_HITS);
346 G4cerr(*G4cerr_p) << "GlueXSensitiveDetectorUPV::ENdOfEvent warning: "
347 << "max hit count " << 2 * MAX_HITS << " exceeded, "
348 << hitscount - 2 * MAX_HITS << " hits discarded."
349 << G4endlstd::endl;
350 }
351 }
352 }
353
354 // Collect and output the barTruthShowers
355 for (piter = points->begin(); piter != points->end(); ++piter) {
356 hddm_s::UpvTruthShowerList point = upv.addUpvTruthShowers(1);
357 point(0).setPrimary(piter->second->primary_);
358 point(0).setPtype(piter->second->ptype_G3);
359 point(0).setPx(piter->second->px_GeV);
360 point(0).setPy(piter->second->py_GeV);
361 point(0).setPz(piter->second->pz_GeV);
362 point(0).setE(piter->second->E_GeV);
363 point(0).setX(piter->second->x_cm);
364 point(0).setY(piter->second->y_cm);
365 point(0).setZ(piter->second->z_cm);
366 point(0).setT(piter->second->t_ns);
367 point(0).setTrack(piter->second->track_);
368 hddm_s::TrackIDList tid = point(0).addTrackIDs();
369 tid(0).setItrack(piter->second->trackID_);
370 }
371}
372
373int GlueXSensitiveDetectorUPV::GetIdent(std::string div,
374 const G4VTouchable *touch)
375{
376 const HddsG4Builder* bldr = GlueXDetectorConstruction::GetBuilder();
377 std::map<std::string, std::vector<int> >::const_iterator iter;
378 std::map<std::string, std::vector<int> > *identifiers;
379 int max_depth = touch->GetHistoryDepth();
380 for (int depth = 0; depth < max_depth; ++depth) {
381 G4VPhysicalVolume *pvol = touch->GetVolume(depth);
382 G4LogicalVolume *lvol = pvol->GetLogicalVolume();
383 int volId = fVolumeTable[lvol];
384 if (volId == 0) {
385 volId = bldr->getVolumeId(lvol);
386 fVolumeTable[lvol] = volId;
387 }
388 identifiers = &Refsys::fIdentifierTable[volId];
389 if ((iter = identifiers->find(div)) != identifiers->end()) {
390 int copyNum = touch->GetCopyNumber(depth);
391 copyNum += (dynamic_cast<G4PVPlacement*>(pvol))? -1 : 0;
392 return iter->second[copyNum];
393 }
394 }
395 return -1;
396}