Bug Summary

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

Annotated Source Code

1//
2// GlueXSensitiveDetectorCERE - class implementation
3//
4// author: richard.t.jones at uconn.edu
5// version: october 29, 2016
6
7#include "GlueXSensitiveDetectorCERE.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// Geant3-style particle index for optical photon
25#define OPTICAL_PHOTON50 50
26
27// Cutoff on the number of allowed hits per section
28int GlueXSensitiveDetectorCERE::MAX_HITS = 100;
29
30// Minimum hit time difference for two hits on the same tube
31double GlueXSensitiveDetectorCERE::TWO_HIT_TIME_RESOL = 50*ns;
32
33// Minimum photoelectron count for a hit
34double GlueXSensitiveDetectorCERE::THRESH_PE = 2.;
35
36int GlueXSensitiveDetectorCERE::instanceCount = 0;
37G4Mutex GlueXSensitiveDetectorCERE::fMutex = G4MUTEX_INITIALIZER{ { 0, 0, 0, 0, 0, 0, 0, { 0, 0 } } };
38
39GlueXSensitiveDetectorCERE::GlueXSensitiveDetectorCERE(const G4String& name)
40 : G4VSensitiveDetector(name),
41 fTubeHitsMap(0), fPointsMap(0)
42{
43 collectionName.insert("CERETubeHitsCollection");
44 collectionName.insert("CEREPointsCollection");
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 CERE, you must delete all old
50 // objects of this class and create new ones.
51
52 G4AutoLock tuberier(&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) << "CERE: ALL parameters loaded from ccdb" << G4endlstd::endl;
66 }
67 }
68}
69
70GlueXSensitiveDetectorCERE::GlueXSensitiveDetectorCERE(
71 const GlueXSensitiveDetectorCERE &src)
72 : G4VSensitiveDetector(src),
73 fTubeHitsMap(src.fTubeHitsMap), fPointsMap(src.fPointsMap)
74{
75 G4AutoLock tuberier(&fMutex);
76 ++instanceCount;
77}
78
79GlueXSensitiveDetectorCERE &GlueXSensitiveDetectorCERE::operator=(const
80 GlueXSensitiveDetectorCERE &src)
81{
82 G4AutoLock tuberier(&fMutex);
83 *(G4VSensitiveDetector*)this = src;
84 fTubeHitsMap = src.fTubeHitsMap;
85 fPointsMap = src.fPointsMap;
86 return *this;
87}
88
89GlueXSensitiveDetectorCERE::~GlueXSensitiveDetectorCERE()
90{
91 G4AutoLock tuberier(&fMutex);
92 --instanceCount;
93}
94
95void GlueXSensitiveDetectorCERE::Initialize(G4HCofThisEvent* hce)
96{
97 fTubeHitsMap = new
98 GlueXHitsMapCEREtube(SensitiveDetectorName, collectionName[0]);
99 fPointsMap = new
100 GlueXHitsMapCEREpoint(SensitiveDetectorName, collectionName[1]);
101 G4SDManager *sdm = G4SDManager::GetSDMpointer();
102 hce->AddHitsCollection(sdm->GetCollectionID(collectionName[0]), fTubeHitsMap);
103 hce->AddHitsCollection(sdm->GetCollectionID(collectionName[1]), fPointsMap);
104}
105
106G4bool GlueXSensitiveDetectorCERE::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 // Post the hit to the points list in the
138 // order of appearance in the event simulation.
139
140 G4Track *track = step->GetTrack();
141 G4int trackID = track->GetTrackID();
142 int pdgtype = track->GetDynamicParticle()->GetPDGcode();
143 int g3type = GlueXPrimaryGeneratorAction::ConvertPdgToGeant3(pdgtype);
144 GlueXUserTrackInformation *trackinfo = (GlueXUserTrackInformation*)
145 track->GetUserInformation();
146 int itrack = trackinfo->GetGlueXTrackID();
147 if (touch->GetVolume()->GetName() == "CERW") {
148 if (trackinfo->GetGlueXHistory() == 0 && itrack > 0 && xin.dot(pin) > 0) {
149 G4int key = fPointsMap->entries();
150 GlueXHitCEREpoint* lastPoint = (*fPointsMap)[key - 1];
151 if (lastPoint == 0 || lastPoint->track_ != trackID ||
152 fabs(lastPoint->t_ns - t/ns) > 0.1 ||
153 fabs(lastPoint->x_cm - x[0]/cm) > 2. ||
154 fabs(lastPoint->y_cm - x[1]/cm) > 2. ||
155 fabs(lastPoint->z_cm - x[2]/cm) > 2.)
156 {
157 GlueXHitCEREpoint newPoint;
158 newPoint.ptype_G3 = g3type;
159 newPoint.track_ = trackID;
160 newPoint.trackID_ = itrack;
161 newPoint.primary_ = (track->GetParentID() == 0);
162 newPoint.t_ns = t/ns;
163 newPoint.x_cm = x[0]/cm;
164 newPoint.y_cm = x[1]/cm;
165 newPoint.z_cm = x[2]/cm;
166 newPoint.px_GeV = pin[0]/GeV;
167 newPoint.py_GeV = pin[1]/GeV;
168 newPoint.pz_GeV = pin[2]/GeV;
169 newPoint.E_GeV = Ein/GeV;
170 fPointsMap->add(key, newPoint);
171 }
172 }
173 // This sensitive detector is unique in that different volumes are used
174 // to generate truth points and hits, which means that a single entry
175 // to ProcessEvents can make one or the other, but not both.
176 return true;
177 }
178
179 // Post the hit to the hits map, ordered by plane,tube,end index
180
181 if (dEsum > 0) {
182 int sector = GetIdent("sector", touch);
183 int key = GlueXHitCEREtube::GetKey(sector);
184 GlueXHitCEREtube *counter = (*fTubeHitsMap)[key];
185 if (counter == 0) {
186 GlueXHitCEREtube newcounter(sector);
187 fTubeHitsMap->add(key, newcounter);
188 counter = (*fTubeHitsMap)[key];
189 }
190
191 // Add the hit to the hits vector, maintaining strict time ordering
192
193 int merge_hit = 0;
194 std::vector<GlueXHitCEREtube::hitinfo_t>::iterator hiter;
195 for (hiter = counter->hits.begin(); hiter != counter->hits.end(); ++hiter) {
196 if (fabs(hiter->t_ns*ns - t) < TWO_HIT_TIME_RESOL) {
197 merge_hit = 1;
198 break;
199 }
200 else if (hiter->t_ns*ns > t) {
201 break;
202 }
203 if (merge_hit) {
204 // Use the time from the earlier hit but add the charge
205 hiter->pe_ += 1;
206 if (hiter->t_ns*ns > t) {
207 hiter->t_ns = t/ns;
208 }
209 }
210 else {
211 // create new hit
212 hiter = counter->hits.insert(hiter, GlueXHitCEREtube::hitinfo_t());
213 hiter->pe_ = 1;
214 hiter->t_ns = t/ns;
215 }
216 }
217 }
218 return true;
219}
220
221void GlueXSensitiveDetectorCERE::EndOfEvent(G4HCofThisEvent*)
222{
223 std::map<int,GlueXHitCEREtube*> *tubes = fTubeHitsMap->GetMap();
224 std::map<int,GlueXHitCEREpoint*> *points = fPointsMap->GetMap();
225 if (tubes->size() == 0 && points->size() == 0)
226 return;
227 std::map<int,GlueXHitCEREtube*>::iterator siter;
228 std::map<int,GlueXHitCEREpoint*>::iterator piter;
229
230 if (verboseLevel > 1) {
231 G4cout(*G4cout_p) << G4endlstd::endl
232 << "--------> Hits Collection: in this event there are "
233 << tubes->size() << " tubes with hits in the Cerenkov: "
234 << G4endlstd::endl;
235 for (siter = tubes->begin(); siter != tubes->end(); ++siter)
236 siter->second->Print();
237
238 G4cout(*G4cout_p) << G4endlstd::endl
239 << "--------> Hits Collection: in this event there are "
240 << points->size() << " truth points in the Cerenkov: "
241 << G4endlstd::endl;
242 for (piter = points->begin(); piter != points->end(); ++piter)
243 piter->second->Print();
244 }
245
246 // pack hits into ouptut hddm record
247
248 G4EventManager* mgr = G4EventManager::GetEventManager();
249 G4VUserEventInformation* info = mgr->GetUserInformation();
250 hddm_s::HDDM *record = ((GlueXUserEventInformation*)info)->getOutputRecord();
251 if (record == 0) {
252 G4cerr(*G4cerr_p) << "GlueXSensitiveDetectorCERE::EndOfEvent error - "
253 << "hits seen but no output hddm record to save them into, "
254 << "cannot continue!" << G4endlstd::endl;
255 exit(1);
256 }
257
258 if (record->getPhysicsEvents().size() == 0)
259 record->addPhysicsEvents();
260 if (record->getHitViews().size() == 0)
261 record->getPhysicsEvent().addHitViews();
262 hddm_s::HitView &hitview = record->getPhysicsEvent().getHitView();
263 if (hitview.getCerenkovs().size() == 0)
264 hitview.addCerenkovs();
265 hddm_s::Cerenkov &cerenkov = hitview.getCerenkov();
266
267 // Collect and output the cereTruthHits
268 for (siter = tubes->begin(); siter != tubes->end(); ++siter) {
269 std::vector<GlueXHitCEREtube::hitinfo_t> &hits = siter->second->hits;
270 // apply a pulse height threshold cut
271 for (unsigned int ih=0; ih < hits.size(); ++ih) {
272 if (hits[ih].pe_ <= THRESH_PE) {
273 hits.erase(hits.begin() + ih);
274 --ih;
275 }
276 }
277 if (hits.size() > 0) {
278 hddm_s::CereSectionList tube = cerenkov.addCereSections(1);
279 tube(0).setSector(siter->second->sector_);
280 int hitscount = hits.size();
281 if (hitscount > MAX_HITS) {
282 hitscount = MAX_HITS;
283 G4cerr(*G4cerr_p) << "GlueXSensitiveDetectorCERE::EndOfEvent warning: "
284 << "max tube hit count " << MAX_HITS << " exceeded, "
285 << hits.size() - hitscount << " hits discarded."
286 << G4endlstd::endl;
287 }
288 for (int ih=0; ih < hitscount; ++ih) {
289 hddm_s::CereTruthHitList thit = tube(0).addCereTruthHits(1);
290 thit(0).setPe(hits[ih].pe_);
291 thit(0).setT(hits[ih].t_ns);
292 }
293 }
294 }
295
296 // Collect and output the cereTruthPoints
297 for (piter = points->begin(); piter != points->end(); ++piter) {
298 hddm_s::CereTruthPointList point = cerenkov.addCereTruthPoints(1);
299 point(0).setPrimary(piter->second->primary_);
300 point(0).setPtype(piter->second->ptype_G3);
301 point(0).setPx(piter->second->px_GeV);
302 point(0).setPy(piter->second->py_GeV);
303 point(0).setPz(piter->second->pz_GeV);
304 point(0).setE(piter->second->E_GeV);
305 point(0).setX(piter->second->x_cm);
306 point(0).setY(piter->second->y_cm);
307 point(0).setZ(piter->second->z_cm);
308 point(0).setT(piter->second->t_ns);
309 point(0).setTrack(piter->second->track_);
310 hddm_s::TrackIDList tid = point(0).addTrackIDs();
311 tid(0).setItrack(piter->second->trackID_);
312 }
313}
314
315int GlueXSensitiveDetectorCERE::GetIdent(std::string div,
316 const G4VTouchable *touch)
317{
318 const HddsG4Builder* bldr = GlueXDetectorConstruction::GetBuilder();
319 std::map<std::string, std::vector<int> >::const_iterator iter;
320 std::map<std::string, std::vector<int> > *identifiers;
321 int max_depth = touch->GetHistoryDepth();
322 for (int depth = 0; depth < max_depth; ++depth) {
323 G4VPhysicalVolume *pvol = touch->GetVolume(depth);
324 G4LogicalVolume *lvol = pvol->GetLogicalVolume();
325 int volId = fVolumeTable[lvol];
326 if (volId == 0) {
327 volId = bldr->getVolumeId(lvol);
328 fVolumeTable[lvol] = volId;
329 }
330 identifiers = &Refsys::fIdentifierTable[volId];
331 if ((iter = identifiers->find(div)) != identifiers->end()) {
332 int copyNum = touch->GetCopyNumber(depth);
333 copyNum += (dynamic_cast<G4PVPlacement*>(pvol))? -1 : 0;
334 return iter->second[copyNum];
335 }
336 }
337 return -1;
338}