Bug Summary

File:programs/Simulation/HDGeant/hitFCal.c
Location:line 55, column 10
Description:Value stored to 'zeroHat' during its initialization is never read

Annotated Source Code

1/*
2 * hitFCal - registers hits for forward calorimeter
3 *
4 * This is a part of the hits package for the
5 * HDGeant simulation program for Hall D.
6 *
7 * version 1.0 -Richard Jones July 16, 2001
8 *
9 * changes: Wed Jun 20 13:19:56 EDT 2007 B. Zihlmann
10 * add ipart to the function hitForwardEMcal
11 *
12 * 3/23/2012 B. Schaefer
13 * Removed radiation hard insert functionality
14 */
15
16#include <stdlib.h>
17#include <stdio.h>
18#include <math.h>
19
20
21#include <hddm_s.h>
22#include <geant3.h>
23#include <bintree.h>
24
25#include "calibDB.h"
26extern s_HDDM_t* thisInputEvent;
27
28
29static float ATTEN_LENGTH = 100.;
30static float C_EFFECTIVE = 15.; // cm/ns
31static float WIDTH_OF_BLOCK = 4.; //cm
32static float LENGTH_OF_BLOCK = 45.; //cm
33static float TWO_HIT_RESOL = 75.; // ns
34static int MAX_HITS = 100; // maximum hits per block
35static float THRESH_MEV = 5.;
36static float ACTIVE_RADIUS = 120.;
37static int CENTRAL_ROW = 29;
38static int CENTRAL_COLUMN = 29;
39
40
41binTree_t* forwardEMcalTree = 0;
42static int blockCount = 0;
43static int showerCount = 0;
44static int initialized = 0;
45
46
47/* register hits during tracking (from gustep) */
48
49void hitForwardEMcal (float xin[4], float xout[4],
50 float pin[5], float pout[5], float dEsum,
51 int track, int stack, int history, int ipart)
52{
53 float x[3], t;
54 float xfcal[3];
55 float zeroHat[] = {0,0,0};
Value stored to 'zeroHat' during its initialization is never read
56
57 if (!initialized){
58
59 mystr_t strings[50];
60 float values[50];
61 int nvalues = 50;
62 int status = GetConstants("FCAL/fcal_parms", &nvalues, values, strings);
63
64 if (!status) {
65
66 int ncounter = 0;
67 int i;
68 for ( i=0;i<(int)nvalues;i++){
69 //printf("%d %s \n",i,strings[i].str);
70 if (!strcmp(strings[i].str,"FCAL_ATTEN_LENGTH")) {
71 ATTEN_LENGTH = values[i];
72 ncounter++;
73 }
74
75 if (!strcmp(strings[i].str,"FCAL_C_EFFECTIVE")) {
76 C_EFFECTIVE = values[i];
77 ncounter++;
78 }
79 if (!strcmp(strings[i].str,"FCAL_WIDTH_OF_BLOCK")) {
80 WIDTH_OF_BLOCK = values[i];
81 ncounter++;
82 }
83 if (!strcmp(strings[i].str,"FCAL_LENGTH_OF_BLOCK")) {
84 LENGTH_OF_BLOCK = values[i];
85 ncounter++;
86 }
87 if (!strcmp(strings[i].str,"FCAL_TWO_HIT_RESOL")) {
88 TWO_HIT_RESOL = values[i];
89 ncounter++;
90 }
91 if (!strcmp(strings[i].str,"FCAL_MAX_HITS")) {
92 MAX_HITS = (int)values[i];
93 ncounter++;
94 }
95 if (!strcmp(strings[i].str,"FCAL_THRESH_MEV")) {
96 THRESH_MEV = values[i];
97 ncounter++;
98 }
99 if (!strcmp(strings[i].str,"FCAL_ACTIVE_RADIUS")) {
100 ACTIVE_RADIUS = values[i];
101 ncounter++;
102 }
103 if (!strcmp(strings[i].str,"FCAL_CENTRAL_ROW")) {
104 CENTRAL_ROW = (int)values[i];
105 ncounter++;
106 }
107 if (!strcmp(strings[i].str,"FCAL_CENTRAL_COLUMN")) {
108 CENTRAL_COLUMN = (int)values[i];
109 ncounter++;
110 }
111 }
112 const int nparams=10;
113 if (ncounter==nparams){
114 printf("FCAL: ALL parameters loaded from Data Base\n");
115 } else if (ncounter<nparams){
116 printf("FCAL: NOT ALL necessary parameters found in Data Base %d out of %d\n",ncounter,nparams);
117 } else {
118 printf("FCAL: SOME parameters found more than once in Data Base\n");
119 }
120 }
121 initialized = 1;
122 }
123
124 x[0] = (xin[0] + xout[0])/2;
125 x[1] = (xin[1] + xout[1])/2;
126 x[2] = (xin[2] + xout[2])/2;
127 t = (xin[3] + xout[3])/2 * 1e9;
128 transformCoord(x,"global",xfcal,"FCAL")transformcoord_(x,"global",xfcal,"FCAL",strlen("global"),strlen
("FCAL"));
129
130 /* post the hit to the truth tree */
131
132 if ((history == 0) && (pin[3] > THRESH_MEV/1e3))
133 {
134 s_FcalTruthShowers_t* showers;
135 int mark = (1<<30) + showerCount;
136 void** twig = getTwig(&forwardEMcalTree, mark);
137 if (*twig == 0)
138 {
139 s_ForwardEMcal_t* cal = *twig = make_s_ForwardEMcal();
140 cal->fcalTruthShowers = showers = make_s_FcalTruthShowers(1);
141 int a = thisInputEvent->physicsEvents->in[0].reactions->in[0].vertices->in[0].products->mult;
142 showers->in[0].primary = (stack <= a);
143 showers->in[0].track = track;
144 showers->in[0].t = xin[3]*1e9;
145 showers->in[0].x = xin[0];
146 showers->in[0].y = xin[1];
147 showers->in[0].z = xin[2];
148 showers->in[0].px = pin[0]*pin[4];
149 showers->in[0].py = pin[1]*pin[4];
150 showers->in[0].pz = pin[2]*pin[4];
151 showers->in[0].E = pin[3];
152 showers->in[0].ptype = ipart;
153 showers->mult = 1;
154 showerCount++;
155 }
156 }
157
158 /* post the hit to the hits tree, mark block as hit */
159
160 if (dEsum > 0)
161 {
162 int nhit;
163 s_FcalTruthHits_t* hits;
164 int row = getrow_();
165 int column = getcolumn_();
166
167 float dist = 0.5*LENGTH_OF_BLOCK-xfcal[2];
168 float dEcorr = dEsum * exp(-dist/ATTEN_LENGTH);
169
170
171
172 float tcorr = t + dist/C_EFFECTIVE;
173 int mark = ((row+1)<<16) + (column+1);
174 void** twig = getTwig(&forwardEMcalTree, mark);
175 if (*twig == 0)
176 {
177 s_ForwardEMcal_t* cal = *twig = make_s_ForwardEMcal();
178 s_FcalBlocks_t* blocks = make_s_FcalBlocks(1);
179 blocks->mult = 1;
180 blocks->in[0].row = row;
181 blocks->in[0].column = column;
182 blocks->in[0].fcalTruthHits = hits = make_s_FcalTruthHits(MAX_HITS);
183 cal->fcalBlocks = blocks;
184 blockCount++;
185 }
186 else
187 {
188 s_ForwardEMcal_t* cal = *twig;
189 hits = cal->fcalBlocks->in[0].fcalTruthHits;
190 }
191
192 for (nhit = 0; nhit < hits->mult; nhit++)
193 {
194 if (fabs(hits->in[nhit].t - tcorr) < TWO_HIT_RESOL)
195 {
196 break;
197 }
198 }
199 if (nhit < hits->mult) /* merge with former hit */
200 {
201 hits->in[nhit].t =
202 (hits->in[nhit].t * hits->in[nhit].E + tcorr*dEcorr)
203 / (hits->in[nhit].E += dEcorr);
204 }
205 else if (nhit < MAX_HITS) /* create new hit */
206 {
207 hits->in[nhit].t = tcorr;
208 hits->in[nhit].E = dEcorr;
209 hits->mult++;
210 }
211 else
212 {
213 fprintf(stderrstderr,"HDGeant error in hitforwardEMcal: ");
214 fprintf(stderrstderr,"max hit count %d exceeded, truncating!\n",MAX_HITS);
215 exit(2);
216 }
217 }
218}
219
220/* entry point from fortran */
221
222void hitforwardemcal_(float* xin, float* xout,
223 float* pin, float* pout, float* dEsum,
224 int* track, int* stack, int* history, int* ipart)
225{
226 hitForwardEMcal(xin,xout,pin,pout,*dEsum,*track,*stack,*history, *ipart);
227}
228
229
230/* pick and package the hits for shipping */
231
232s_ForwardEMcal_t* pickForwardEMcal ()
233{
234 s_ForwardEMcal_t* box;
235 s_ForwardEMcal_t* item;
236
237#if TESTING_CAL_CONTAINMENT
238 double Etotal = 0;
239#endif
240 if ((blockCount == 0) && (showerCount == 0))
241 {
242 return HDDM_NULL(void*)&hddm_s_nullTarget;
243 }
244
245 box = make_s_ForwardEMcal();
246 box->fcalBlocks = make_s_FcalBlocks(blockCount);
247 box->fcalTruthShowers = make_s_FcalTruthShowers(showerCount);
248 while (item = (s_ForwardEMcal_t*) pickTwig(&forwardEMcalTree))
249 {
250 s_FcalBlocks_t* blocks = item->fcalBlocks;
251 int block;
252 s_FcalTruthShowers_t* showers = item->fcalTruthShowers;
253 int shower;
254 for (block=0; block < blocks->mult; ++block)
255 {
256 int row = blocks->in[block].row;
257 int column = blocks->in[block].column;
258 float y0 = (row - CENTRAL_ROW)*WIDTH_OF_BLOCK;
259 float x0 = (column - CENTRAL_COLUMN)*WIDTH_OF_BLOCK;
260 float dist = sqrt(x0*x0+y0*y0);
261
262 s_FcalTruthHits_t* hits = blocks->in[block].fcalTruthHits;
263
264 /* compress out the hits outside the active region */
265 if (dist < ACTIVE_RADIUS)
266 {
267 int m = box->fcalBlocks->mult;
268
269 /* compress out the hits below threshold */
270 int i,iok;
271 for (iok=i=0; i < hits->mult; i++)
272 {
273 if (hits->in[i].E >= THRESH_MEV/1e3)
274 {
275#if TESTING_CAL_CONTAINMENT
276 Etotal += hits->in[i].E;
277#endif
278 if (iok < i)
279 {
280 hits->in[iok] = hits->in[i];
281 }
282 ++iok;
283 }
284 }
285 if (iok)
286 {
287 hits->mult = iok;
288 box->fcalBlocks->in[m] = blocks->in[block];
289 box->fcalBlocks->mult++;
290 }
291 else if (hits != HDDM_NULL(void*)&hddm_s_nullTarget)
292 {
293 FREE(hits)free(hits);
294 }
295 }
296 else if (hits != HDDM_NULL(void*)&hddm_s_nullTarget)
297 {
298 FREE(hits)free(hits);
299 }
300 }
301
302 for (shower=0; shower < showers->mult; ++shower)
303 {
304 int m = box->fcalTruthShowers->mult++;
305 box->fcalTruthShowers->in[m] = showers->in[shower];
306 }
307 if (blocks != HDDM_NULL(void*)&hddm_s_nullTarget)
308 {
309 FREE(blocks)free(blocks);
310 }
311 if (showers != HDDM_NULL(void*)&hddm_s_nullTarget)
312 {
313 FREE(showers)free(showers);
314 }
315 FREE(item)free(item);
316 }
317
318 blockCount = showerCount = 0;
319
320 if ((box->fcalBlocks != HDDM_NULL(void*)&hddm_s_nullTarget) &&
321 (box->fcalBlocks->mult == 0))
322 {
323 FREE(box->fcalBlocks)free(box->fcalBlocks);
324 box->fcalBlocks = HDDM_NULL(void*)&hddm_s_nullTarget;
325 }
326 if ((box->fcalTruthShowers != HDDM_NULL(void*)&hddm_s_nullTarget) &&
327 (box->fcalTruthShowers->mult == 0))
328 {
329 FREE(box->fcalTruthShowers)free(box->fcalTruthShowers);
330 box->fcalTruthShowers = HDDM_NULL(void*)&hddm_s_nullTarget;
331 }
332 if ((box->fcalBlocks->mult == 0) &&
333 (box->fcalTruthShowers->mult == 0))
334 {
335 FREE(box)free(box);
336 box = HDDM_NULL(void*)&hddm_s_nullTarget;
337 }
338#if TESTING_CAL_CONTAINMENT
339 printf("FCal energy sum: %f\n",Etotal/0.614);
340#endif
341 return box;
342}