Bug Summary

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

Annotated Source Code

1//
2// GlueXSensitiveDetectorPSC - class implementation
3//
4// author: richard.t.jones at uconn.edu
5// version: october 26, 2016
6
7#include "GlueXSensitiveDetectorPSC.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 GlueXSensitiveDetectorPSC::MAX_HITS = 100;
26
27// Minimum hit time difference for two hits on the same paddle
28double GlueXSensitiveDetectorPSC::TWO_HIT_TIME_RESOL = 25*ns;
29
30// Minimum energy deposition for a hit
31double GlueXSensitiveDetectorPSC::THRESH_MEV = 0.010;
32
33// Geometric parameters
34int GlueXSensitiveDetectorPSC::NUM_MODULES_PER_ARM = 8;
35
36int GlueXSensitiveDetectorPSC::instanceCount = 0;
37G4Mutex GlueXSensitiveDetectorPSC::fMutex = G4MUTEX_INITIALIZER{ { 0, 0, 0, 0, 0, 0, 0, { 0, 0 } } };
38
39GlueXSensitiveDetectorPSC::GlueXSensitiveDetectorPSC(const G4String& name)
40 : G4VSensitiveDetector(name),
41 fCounterHitsMap(0), fPointsMap(0)
42{
43 collectionName.insert("PSCPaddleHitsCollection");
44 collectionName.insert("PSCPointsCollection");
45
46 // The rest of this only needs to happen once, the first time an object
47 // of this type is instantiated for this configuration of geometry and
48 // fields. If the geometry or fields change in such a way as to modify
49 // the drift-time properties of hits in the PSC, you must delete all old
50 // objects of this class and create new ones.
51
52 G4AutoLock barrier(&fMutex);
53 if (instanceCount++ == 0) {
54 int runno = HddmOutput::getRunNo();
55 extern jana::JApplication *japp;
56 if (japp == 0) {
57 G4cerr(*G4cerr_p) << "Error in GlueXSensitiveDetector constructor - "
58 << "jana global DApplication object not set, "
59 << "cannot continue." << G4endlstd::endl;
60 exit(-1);
61 }
62 jana::JCalibration *jcalib = japp->GetJCalibration(runno);
Value stored to 'jcalib' during its initialization is never read
63 if (japp == 0) { // dummy
64 jcalib = 0;
65 G4cout(*G4cout_p) << "PSC: ALL parameters loaded from ccdb" << G4endlstd::endl;
66 }
67 }
68}
69
70GlueXSensitiveDetectorPSC::GlueXSensitiveDetectorPSC(
71 const GlueXSensitiveDetectorPSC &src)
72 : G4VSensitiveDetector(src),
73 fCounterHitsMap(src.fCounterHitsMap), fPointsMap(src.fPointsMap)
74{
75 G4AutoLock barrier(&fMutex);
76 ++instanceCount;
77}
78
79GlueXSensitiveDetectorPSC &GlueXSensitiveDetectorPSC::operator=(const
80 GlueXSensitiveDetectorPSC &src)
81{
82 G4AutoLock barrier(&fMutex);
83 *(G4VSensitiveDetector*)this = src;
84 fCounterHitsMap = src.fCounterHitsMap;
85 fPointsMap = src.fPointsMap;
86 return *this;
87}
88
89GlueXSensitiveDetectorPSC::~GlueXSensitiveDetectorPSC()
90{
91 G4AutoLock barrier(&fMutex);
92 --instanceCount;
93}
94
95void GlueXSensitiveDetectorPSC::Initialize(G4HCofThisEvent* hce)
96{
97 fCounterHitsMap = new
98 GlueXHitsMapPSCpaddle(SensitiveDetectorName, collectionName[0]);
99 fPointsMap = new
100 GlueXHitsMapPSCpoint(SensitiveDetectorName, collectionName[1]);
101 G4SDManager *sdm = G4SDManager::GetSDMpointer();
102 hce->AddHitsCollection(sdm->GetCollectionID(collectionName[0]), fCounterHitsMap);
103 hce->AddHitsCollection(sdm->GetCollectionID(collectionName[1]), fPointsMap);
104}
105
106G4bool GlueXSensitiveDetectorPSC::ProcessHits(G4Step* step,
107 G4TouchableHistory* ROhist)
108{
109 double dEsum = step->GetTotalEnergyDeposit();
110 if (dEsum == 0)
111 return false;
112
113 const G4ThreeVector &pin = step->GetPreStepPoint()->GetMomentum();
114 const G4ThreeVector &xin = step->GetPreStepPoint()->GetPosition();
115 const G4ThreeVector &xout = step->GetPostStepPoint()->GetPosition();
116 double Ein = step->GetPreStepPoint()->GetTotalEnergy();
117 double tin = step->GetPreStepPoint()->GetGlobalTime();
118 double tout = step->GetPostStepPoint()->GetGlobalTime();
119 G4ThreeVector x = (xin + xout) / 2;
120 G4ThreeVector dx = xout - xin;
121 double t = (tin + tout) / 2;
122
123 const G4VTouchable* touch = step->GetPreStepPoint()->GetTouchable();
124 const G4AffineTransform &local_from_global = touch->GetHistory()
125 ->GetTopTransform();
126 G4ThreeVector xlocal = local_from_global.TransformPoint(x);
127
128 // For particles that range out inside the active volume, the
129 // "out" time may sometimes be set to something enormously high.
130 // This screws up the hit. Check for this case here by looking
131 // at tout and making sure it is less than 1 second. If it's
132 // not, then just use tin for "t".
133
134 if (tout > 1.0*s)
135 t = tin;
136
137 double dr = dx.mag();
138 double dEdx = (dr > 1e-3*cm)? dEsum/dr : 0;
139
140 // Post the hit to the points list in the
141 // order of appearance in the event simulation.
142
143 int module = GetIdent("module", touch);
144 int arm = (module - 1) / NUM_MODULES_PER_ARM;
145 module = (module - 1) % NUM_MODULES_PER_ARM + 1;
146 G4Track *track = step->GetTrack();
147 G4int trackID = track->GetTrackID();
148 int pdgtype = track->GetDynamicParticle()->GetPDGcode();
149 int g3type = GlueXPrimaryGeneratorAction::ConvertPdgToGeant3(pdgtype);
150 GlueXUserTrackInformation *trackinfo = (GlueXUserTrackInformation*)
151 track->GetUserInformation();
152 int itrack = trackinfo->GetGlueXTrackID();
153 if (trackinfo->GetGlueXHistory() == 0 && itrack > 0) {
154 G4int key = fPointsMap->entries();
155 GlueXHitPSCpoint* lastPoint = (*fPointsMap)[key - 1];
156 if (lastPoint == 0 || lastPoint->track_ != trackID ||
157 fabs(lastPoint->t_ns - t/ns) > 0.2 ||
158 fabs(lastPoint->x_cm - x[0]/cm) > 5.0 ||
159 fabs(lastPoint->y_cm - x[1]/cm) > 5.0 ||
160 fabs(lastPoint->z_cm - x[2]/cm) > 5.0)
161 {
162 GlueXHitPSCpoint newPoint;
163 newPoint.arm_ = arm;
164 newPoint.module_ = module;
165 newPoint.ptype_G3 = g3type;
166 newPoint.track_ = trackID;
167 newPoint.trackID_ = itrack;
168 newPoint.primary_ = (track->GetParentID() == 0);
169 newPoint.t_ns = t/ns;
170 newPoint.x_cm = x[0]/cm;
171 newPoint.y_cm = x[1]/cm;
172 newPoint.z_cm = x[2]/cm;
173 newPoint.px_GeV = pin[0]/GeV;
174 newPoint.py_GeV = pin[1]/GeV;
175 newPoint.pz_GeV = pin[2]/GeV;
176 newPoint.E_GeV = Ein/GeV;
177 newPoint.dEdx_GeV_cm = dEdx/(GeV/cm);
178 fPointsMap->add(key, newPoint);
179 }
180 }
181
182 // Post the hit to the hits map, ordered by arm,module index
183
184 if (dEsum > 0) {
185 int key = GlueXHitPSCpaddle::GetKey(arm, module);
186 GlueXHitPSCpaddle *paddle = (*fCounterHitsMap)[key];
187 if (paddle == 0) {
188 GlueXHitPSCpaddle newpaddle(arm, module);
189 fCounterHitsMap->add(key, newpaddle);
190 paddle = (*fCounterHitsMap)[key];
191 }
192
193 // Add the hit to the hits vector, maintaining strict time ordering
194
195 int merge_hit = 0;
196 std::vector<GlueXHitPSCpaddle::hitinfo_t>::iterator hiter;
197 for (hiter = paddle->hits.begin(); hiter != paddle->hits.end(); ++hiter) {
198 if (fabs(hiter->t_ns*ns - t) < TWO_HIT_TIME_RESOL) {
199 merge_hit = 1;
200 break;
201 }
202 else if (hiter->t_ns*ns > t) {
203 break;
204 }
205 }
206 if (merge_hit) {
207 // Add the charge, do energy-weighted time averaging
208 hiter->t_ns = (hiter->t_ns * hiter->dE_GeV + t/ns * dEsum/GeV) /
209 (hiter->dE_GeV + dEsum/GeV);
210 hiter->dE_GeV += dEsum/GeV;
211 }
212 else {
213 // create new hit
214 hiter = paddle->hits.insert(hiter, GlueXHitPSCpaddle::hitinfo_t());
215 hiter->dE_GeV = dEsum/GeV;
216 hiter->t_ns = t/ns;
217 hiter->itrack_ = itrack;
218 hiter->ptype_G3 = g3type;
219 }
220 }
221 return true;
222}
223
224void GlueXSensitiveDetectorPSC::EndOfEvent(G4HCofThisEvent*)
225{
226 std::map<int,GlueXHitPSCpaddle*> *paddles = fCounterHitsMap->GetMap();
227 std::map<int,GlueXHitPSCpoint*> *points = fPointsMap->GetMap();
228 if (paddles->size() == 0 && points->size() == 0)
229 return;
230 std::map<int,GlueXHitPSCpaddle*>::iterator siter;
231 std::map<int,GlueXHitPSCpoint*>::iterator piter;
232
233 if (verboseLevel > 1) {
234 G4cout(*G4cout_p) << G4endlstd::endl
235 << "--------> Hits Collection: in this event there are "
236 << paddles->size() << " paddles with hits in the PSC: "
237 << G4endlstd::endl;
238 for (siter = paddles->begin(); siter != paddles->end(); ++siter)
239 siter->second->Print();
240
241 G4cout(*G4cout_p) << G4endlstd::endl
242 << "--------> Hits Collection: in this event there are "
243 << points->size() << " truth points in the PSC: "
244 << G4endlstd::endl;
245 for (piter = points->begin(); piter != points->end(); ++piter)
246 piter->second->Print();
247 }
248
249 // pack hits into ouptut hddm record
250
251 G4EventManager* mgr = G4EventManager::GetEventManager();
252 G4VUserEventInformation* info = mgr->GetUserInformation();
253 hddm_s::HDDM *record = ((GlueXUserEventInformation*)info)->getOutputRecord();
254 if (record == 0) {
255 G4cerr(*G4cerr_p) << "GlueXSensitiveDetectorPSC::EndOfEvent error - "
256 << "hits seen but no output hddm record to save them into, "
257 << "cannot continue!" << G4endlstd::endl;
258 exit(1);
259 }
260
261 if (record->getPhysicsEvents().size() == 0)
262 record->addPhysicsEvents();
263 if (record->getHitViews().size() == 0)
264 record->getPhysicsEvent().addHitViews();
265 hddm_s::HitView &hitview = record->getPhysicsEvent().getHitView();
266 if (hitview.getPairSpectrometerCoarses().size() == 0)
267 hitview.addPairSpectrometerCoarses();
268 hddm_s::PairSpectrometerCoarse &psc = hitview.getPairSpectrometerCoarse();
269
270 // Collect and output the PscTruthHits
271 for (siter = paddles->begin(); siter != paddles->end(); ++siter) {
272 std::vector<GlueXHitPSCpaddle::hitinfo_t> &hits = siter->second->hits;
273 // apply a pulse height threshold cut
274 for (unsigned int ih=0; ih < hits.size(); ++ih) {
275 if (hits[ih].dE_GeV*1e3 <= THRESH_MEV) {
276 hits.erase(hits.begin() + ih);
277 --ih;
278 }
279 }
280 if (hits.size() > 0) {
281 hddm_s::PscPaddleList paddle = psc.addPscPaddles(1);
282 paddle(0).setArm(siter->second->arm_);
283 paddle(0).setModule(siter->second->module_);
284 int hitscount = hits.size();
285 if (hitscount > MAX_HITS) {
286 hitscount = MAX_HITS;
287 G4cerr(*G4cerr_p) << "GlueXSensitiveDetectorPSC::EndOfEvent warning: "
288 << "max hit count " << MAX_HITS << " exceeded, "
289 << hits.size() - hitscount << " hits discarded."
290 << G4endlstd::endl;
291 }
292 for (int ih=0; ih < hitscount; ++ih) {
293 hddm_s::PscTruthHitList thit = paddle(0).addPscTruthHits(1);
294 thit(0).setDE(hits[ih].dE_GeV);
295 thit(0).setT(hits[ih].t_ns);
296 thit(0).setItrack(hits[ih].itrack_);
297 thit(0).setPtype(hits[ih].ptype_G3);
298 }
299 }
300 }
301
302 // Collect and output the pscTruthPoints
303 for (piter = points->begin(); piter != points->end(); ++piter) {
304 hddm_s::PscTruthPointList point = psc.addPscTruthPoints(1);
305 point(0).setE(piter->second->E_GeV);
306 point(0).setDEdx(piter->second->dEdx_GeV_cm);
307 point(0).setPrimary(piter->second->primary_);
308 point(0).setPtype(piter->second->ptype_G3);
309 point(0).setPx(piter->second->px_GeV);
310 point(0).setPy(piter->second->py_GeV);
311 point(0).setPz(piter->second->pz_GeV);
312 point(0).setX(piter->second->x_cm);
313 point(0).setY(piter->second->y_cm);
314 point(0).setZ(piter->second->z_cm);
315 point(0).setT(piter->second->t_ns);
316 point(0).setArm(piter->second->arm_);
317 point(0).setModule(piter->second->module_);
318 point(0).setTrack(piter->second->track_);
319 hddm_s::TrackIDList tid = point(0).addTrackIDs();
320 tid(0).setItrack(piter->second->trackID_);
321 }
322}
323
324int GlueXSensitiveDetectorPSC::GetIdent(std::string div,
325 const G4VTouchable *touch)
326{
327 const HddsG4Builder* bldr = GlueXDetectorConstruction::GetBuilder();
328 std::map<std::string, std::vector<int> >::const_iterator iter;
329 std::map<std::string, std::vector<int> > *identifiers;
330 int max_depth = touch->GetHistoryDepth();
331 for (int depth = 0; depth < max_depth; ++depth) {
332 G4VPhysicalVolume *pvol = touch->GetVolume(depth);
333 G4LogicalVolume *lvol = pvol->GetLogicalVolume();
334 int volId = fVolumeTable[lvol];
335 if (volId == 0) {
336 volId = bldr->getVolumeId(lvol);
337 fVolumeTable[lvol] = volId;
338 }
339 identifiers = &Refsys::fIdentifierTable[volId];
340 if ((iter = identifiers->find(div)) != identifiers->end()) {
341 int copyNum = touch->GetCopyNumber(depth);
342 copyNum += (dynamic_cast<G4PVPlacement*>(pvol))? -1 : 0;
343 return iter->second[copyNum];
344 }
345 }
346 return -1;
347}