File: | programs/Simulation/HDGeant/hitFTOF.c |
Location: | line 274, column 7 |
Description: | Value stored to 'noMCHits' is never read |
1 | /* |
2 | * hitFTOF - registers hits for forward Time-Of-Flight |
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: -B. Zihlmann June 19. 2007 |
10 | * add hit position to north and south hit structure |
11 | * set THRESH_MEV to zero to NOT concatenate hits. |
12 | * |
13 | * changes: Wed Jun 20 13:19:56 EDT 2007 B. Zihlmann |
14 | * add ipart to the function hitforwardTOF |
15 | * |
16 | * Programmer's Notes: |
17 | * ------------------- |
18 | * 1) In applying the attenuation to light propagating down to both ends |
19 | * of the counters, there has to be some point where the attenuation |
20 | * factor is 1. I chose it to be the midplane, so that in the middle |
21 | * of the counter both ends see the unattenuated dE values. Closer to |
22 | * either end, that end has a larger dE value and the opposite end a |
23 | * lower dE value than the actual deposition. |
24 | * 2) In applying the propagation delay to light propagating down to the |
25 | * ends of the counters, there has to be some point where the timing |
26 | * offset is 0. I chose it to be the midplane, so that for hits in |
27 | * the middle of the counter the t values measure time-of-flight from |
28 | * the t=0 of the event. For hits closer to one end, that end sees |
29 | * a t value smaller than its true time-of-flight, and the other end |
30 | * sees a value correspondingly larger. The average is the true tof. |
31 | */ |
32 | |
33 | #include <stdlib.h> |
34 | #include <stdio.h> |
35 | #include <math.h> |
36 | |
37 | #include <hddm_s.h> |
38 | #include <geant3.h> |
39 | #include <bintree.h> |
40 | |
41 | #include "calibDB.h" |
42 | |
43 | // plastic scintillator specific constants |
44 | static float ATTEN_LENGTH = 150 ; |
45 | static float C_EFFECTIVE = 15.0 ; |
46 | static float BAR_LENGTH = 252.0; // length of the bar |
47 | |
48 | // kinematic constants |
49 | static float TWO_HIT_RESOL = 25. ;// separation time between two different hits |
50 | |
51 | static float THRESH_MEV = 0. ; // do not through away any hits, one can do that later |
52 | |
53 | // maximum particle tracks per counter |
54 | static int MAX_HITS = 25; // was 100 changed to 25 |
55 | |
56 | // maximum MC hits per paddle |
57 | static int MAX_PAD_HITS = 25 ; |
58 | |
59 | // top level pointer of FTOF hit tree |
60 | binTree_t* forwardTOFTree = 0; |
61 | |
62 | static int counterCount = 0; |
63 | static int pointCount = 0; |
64 | static int initialized = 0; |
65 | |
66 | |
67 | /* register hits during tracking (from gustep) */ |
68 | // track is ITRA from GEANT |
69 | // stack is ISTAK from GEANT |
70 | // history is ISTORY from GEANT User flag for current track history (reset to 0 in GLTRAC) |
71 | |
72 | void hitForwardTOF (float xin[4], float xout[4], |
73 | float pin[5], float pout[5], float dEsum, |
74 | int track, int stack, int history, int ipart) { |
75 | float x[3], t; |
76 | float dx[3], dr; |
77 | float dEdx; |
78 | float xlocal[3]; |
79 | float xftof[3]; |
80 | float zeroHat[] = {0,0,0}; |
81 | |
82 | if (!initialized) { |
83 | |
84 | |
85 | mystr_t strings[50]; |
86 | float values[50]; |
87 | int nvalues = 50; |
88 | int status = GetConstants("TOF/tof_parms", &nvalues, values, strings); |
89 | |
90 | if (!status) { |
91 | |
92 | int ncounter = 0; |
93 | int i; |
94 | for ( i=0;i<(int)nvalues;i++){ |
95 | //printf("%d %s \n",i,strings[i].str); |
96 | if (!strcmp(strings[i].str,"TOF_ATTEN_LENGTH")) { |
97 | ATTEN_LENGTH = values[i]; |
98 | ncounter++; |
99 | } |
100 | if (!strcmp(strings[i].str,"TOF_C_EFFECTIVE")) { |
101 | C_EFFECTIVE = values[i]; |
102 | ncounter++; |
103 | } |
104 | if (!strcmp(strings[i].str,"TOF_PADDLE_LENGTH")) { |
105 | BAR_LENGTH = values[i]; |
106 | ncounter++; |
107 | } |
108 | if (!strcmp(strings[i].str,"TOF_TWO_HIT_RESOL")) { |
109 | TWO_HIT_RESOL = values[i]; |
110 | ncounter++; |
111 | } |
112 | if (!strcmp(strings[i].str,"TOF_THRESH_MEV")) { |
113 | THRESH_MEV = values[i]; |
114 | ncounter++; |
115 | } |
116 | if (!strcmp(strings[i].str,"TOF_MAX_HITS")){ |
117 | MAX_HITS = (int)values[i]; |
118 | ncounter++; |
119 | } |
120 | if (!strcmp(strings[i].str,"TOF_MAX_PAD_HITS")) { |
121 | MAX_PAD_HITS = (int)values[i]; |
122 | ncounter++; |
123 | } |
124 | } |
125 | if (ncounter==7){ |
126 | printf("TOF: ALL parameters loaded from Data Base\n"); |
127 | } else if (ncounter<7){ |
128 | printf("TOF: NOT ALL necessary parameters found in Data Base %d out of 7\n",ncounter); |
129 | } else { |
130 | printf("TOF: SOME parameters found more than once in Data Base\n"); |
131 | } |
132 | |
133 | } |
134 | initialized = 1; |
135 | |
136 | } |
137 | |
138 | |
139 | |
140 | // getplane is coded in |
141 | // src/programs/Simulation/HDGeant/hddsGeant3.F |
142 | // this file is automatically generated from the geometry file |
143 | // written in xml format |
144 | // NOTE: there are three files hddsGeant3.F with the same name in |
145 | // the source code tree namely |
146 | // 1) src/programs/Utilities/geantbfield2root/hddsGeant3.F |
147 | // 2) src/programs/Simulation/HDGeant/hddsGeant3.F |
148 | // 3) src/programs/Simulation/hdds/hddsGeant3.F |
149 | // |
150 | // while 2) and 3) are identical 1) is a part of 2) and 3) |
151 | int plane = getplane_(); |
152 | |
153 | // calculate mean location of track and mean time in [ns] units |
154 | // the units of xin xout and x are in [cm] |
155 | x[0] = (xin[0] + xout[0])/2; |
156 | x[1] = (xin[1] + xout[1])/2; |
157 | x[2] = (xin[2] + xout[2])/2; |
158 | t = (xin[3] + xout[3])/2 * 1e9; |
159 | |
160 | // tranform the the global x coordinate into the local coordinate of the top_volume FTOF |
161 | // defined in the geometry file src/programs/Simulation/hdds/ForwardTOF_HDDS.xml |
162 | // the function transform Coord is defined in src/programs/Simulation/HDGeant/hitutil/hitutil.F |
163 | transformCoord(x,"global",xlocal,"FTOF")transformcoord_(x,"global",xlocal,"FTOF",strlen("global"),strlen ("FTOF")); |
164 | transformCoord(zeroHat,"local",xftof,"FTOF")transformcoord_(zeroHat,"local",xftof,"FTOF",strlen("local"), strlen("FTOF")); |
165 | |
166 | // track vector of this step |
167 | dx[0] = xin[0] - xout[0]; |
168 | dx[1] = xin[1] - xout[1]; |
169 | dx[2] = xin[2] - xout[2]; |
170 | // length of the track of this step |
171 | dr = sqrt(dx[0]*dx[0] + dx[1]*dx[1] + dx[2]*dx[2]); |
172 | // calculate dEdx only if track length is >0.001 cm |
173 | if (dr > 1e-3) { |
174 | dEdx = dEsum/dr; |
175 | } |
176 | else { |
177 | dEdx = 0; |
178 | } |
179 | |
180 | /* post the hit to the truth tree */ |
181 | // in other words: store the GENERATED track information |
182 | |
183 | if ((history == 0) && (plane == 0)) { |
184 | |
185 | // save all tracks from particles that hit the first plane of FTOF |
186 | // save the generated "true" values |
187 | |
188 | int mark = (1<<30) + pointCount; |
189 | |
190 | // getTwig is defined in src/programs/Simulation/HDGeant/bintree.c |
191 | // the two arguments are a pointer to a pointer and an integer |
192 | |
193 | void** twig = getTwig(&forwardTOFTree, mark); |
194 | if (*twig == 0) { |
195 | // make_s_ForwardTOF is defined in src/programs/Analysis/hddm/hddm_s.h |
196 | // and coded in src/programs/Analysis/hddm/hddm_s.c |
197 | // the same holds for make_s_FtofTruthPoints |
198 | |
199 | // make_s_ForwardTOF returns pointer to structure s_ForwardTOF generated memory |
200 | // tof->ftofCoutners and tof-> ftofTruthPoints are initialized already |
201 | |
202 | s_ForwardTOF_t* tof = *twig = make_s_ForwardTOF(); |
203 | s_FtofTruthPoints_t* points = make_s_FtofTruthPoints(1); |
204 | tof->ftofTruthPoints = points; |
205 | points->in[0].primary = (stack == 0); |
206 | points->in[0].track = track; |
207 | points->in[0].x = x[0]; |
208 | points->in[0].y = x[1]; |
209 | points->in[0].z = x[2]; |
210 | points->in[0].t = t; |
211 | points->in[0].px = pin[0]*pin[4]; |
212 | points->in[0].py = pin[1]*pin[4]; |
213 | points->in[0].pz = pin[2]*pin[4]; |
214 | points->in[0].E = pin[3]; |
215 | points->in[0].ptype = ipart; |
216 | points->mult = 1; |
217 | pointCount++; |
218 | } |
219 | } |
220 | |
221 | /* post the hit to the hits tree, mark slab as hit */ |
222 | // in other words now store the simulated detector response |
223 | if (dEsum > 0) { |
224 | int nhit; |
225 | s_FtofNorthTruthHits_t* northHits; |
226 | s_FtofSouthTruthHits_t* southHits; |
227 | |
228 | s_FtofMCHits_t* noMCHits; |
229 | s_FtofMCHits_t* soMCHits; |
230 | |
231 | // getrow and getcolumn are both coded in hddsGeant3.F |
232 | // see above for function getplane() |
233 | |
234 | int row = getrow_(); |
235 | int column = getcolumn_(); |
236 | |
237 | // distance of hit from PMT north w.r.t. center and similar for PMT south |
238 | // this means positive x points north. to get a right handed system y must |
239 | // point vertically up as z is the beam axis. |
240 | // plane==0 horizontal plane, plane==1 vertical plane |
241 | // float dist = xlocal[0]; |
242 | |
243 | float dist = x[1]; // do not use local coordinate for x and y |
244 | if (plane==1) dist = x[0]; |
245 | float dxnorth = BAR_LENGTH/2.-dist; |
246 | float dxsouth = BAR_LENGTH/2.+dist; |
247 | |
248 | // calculate time at the PMT "normalized" to the center, so a hit in the |
249 | // center will have time "t" at both PMTs |
250 | // the speed of signal travel is C_EFFECTIVE |
251 | // propagte time to the end of the bar |
252 | // column = 0 is a full paddle column ==1,2 is a half paddle |
253 | |
254 | float tnorth = (column == 2) ? 0 : t + dxnorth/C_EFFECTIVE; |
255 | float tsouth = (column == 1) ? 0 : t + dxsouth/C_EFFECTIVE; |
256 | |
257 | // calculate energy seen by PM for this track step using attenuation factor |
258 | float dEnorth = (column == 2) ? 0 : dEsum * exp(-dxnorth/ATTEN_LENGTH); |
259 | float dEsouth = (column == 1) ? 0 : dEsum * exp(-dxsouth/ATTEN_LENGTH); |
260 | |
261 | int mark = (plane<<20) + (row<<10) + column; |
262 | void** twig = getTwig(&forwardTOFTree, mark); |
263 | |
264 | if (*twig == 0) { // this paddle has not been hit yet by any particle track |
265 | // get space and store it |
266 | |
267 | s_ForwardTOF_t* tof = *twig = make_s_ForwardTOF(); |
268 | s_FtofCounters_t* counters = make_s_FtofCounters(1); |
269 | counters->mult = 1; |
270 | counters->in[0].plane = plane; |
271 | counters->in[0].bar = row; |
272 | northHits = HDDM_NULL(void*)&hddm_s_nullTarget; |
273 | southHits = HDDM_NULL(void*)&hddm_s_nullTarget; |
274 | noMCHits = HDDM_NULL(void*)&hddm_s_nullTarget; |
Value stored to 'noMCHits' is never read | |
275 | soMCHits = HDDM_NULL(void*)&hddm_s_nullTarget; |
276 | |
277 | // get space for the left/top or right/down PMT data for a total |
278 | // of MAX_HITS possible hits in a single paddle |
279 | // and space for up to MAX_PAD_HITS hits in a paddle to store MC track information |
280 | // Note: column=0 means paddle read out on both ends coumn=1or2 means single ended readout |
281 | |
282 | if (column == 0 || column == 1) { |
283 | counters->in[0].ftofNorthTruthHits = northHits |
284 | = make_s_FtofNorthTruthHits(MAX_HITS); |
285 | } |
286 | if (column == 0 || column == 2) { |
287 | counters->in[0].ftofSouthTruthHits = southHits |
288 | = make_s_FtofSouthTruthHits(MAX_HITS); |
289 | } |
290 | |
291 | tof->ftofCounters = counters; |
292 | counterCount++; |
293 | |
294 | } else { |
295 | |
296 | // this paddle is already registered (was hit before) |
297 | // get the hit list back |
298 | s_ForwardTOF_t* tof = *twig; |
299 | northHits = tof->ftofCounters->in[0].ftofNorthTruthHits; |
300 | southHits = tof->ftofCounters->in[0].ftofSouthTruthHits; |
301 | |
302 | } |
303 | |
304 | if (northHits != HDDM_NULL(void*)&hddm_s_nullTarget) { |
305 | |
306 | // loop over hits in this PM to find correct time slot |
307 | |
308 | for (nhit = 0; nhit < northHits->mult; nhit++) { |
309 | |
310 | if (fabs(northHits->in[nhit].t - t) < TWO_HIT_RESOL) { |
311 | break; |
312 | } |
313 | } |
314 | |
315 | // this hit is within the time frame of a previous hit |
316 | // combine the times of this weighted by the energy of the hit |
317 | |
318 | if (nhit < northHits->mult) { /* merge with former hit */ |
319 | northHits->in[nhit].t = |
320 | (northHits->in[nhit].t * northHits->in[nhit].dE |
321 | + tnorth * dEnorth) |
322 | / (northHits->in[nhit].dE += dEnorth); |
323 | |
324 | // now add MC tracking information |
325 | // first get MC pointer of this paddle |
326 | |
327 | noMCHits = northHits->in[nhit].ftofMCHits; |
328 | unsigned int nMChit = noMCHits->mult; |
329 | if (nMChit<MAX_PAD_HITS) { |
330 | noMCHits->in[nMChit].x = x[0]; |
331 | noMCHits->in[nMChit].y = x[1]; |
332 | noMCHits->in[nMChit].z = x[2]; |
333 | noMCHits->in[nMChit].E = pin[3]; |
334 | noMCHits->in[nMChit].px = pin[0]*pin[4]; |
335 | noMCHits->in[nMChit].py = pin[1]*pin[4]; |
336 | noMCHits->in[nMChit].pz = pin[2]*pin[4]; |
337 | noMCHits->in[nMChit].ptype = ipart; |
338 | noMCHits->in[nMChit].itrack = track; |
339 | noMCHits->in[nMChit].dist = dist; |
340 | noMCHits->mult++; |
341 | } |
342 | |
343 | } else if (nhit < MAX_HITS){ // hit in new time window |
344 | northHits->in[nhit].t = tnorth; |
345 | northHits->in[nhit].dE = dEnorth; |
346 | northHits->mult++; |
347 | |
348 | // create memory for MC track hit information |
349 | northHits->in[nhit].ftofMCHits = noMCHits = make_s_FtofMCHits(MAX_PAD_HITS); |
350 | |
351 | noMCHits->in[0].x = x[0]; |
352 | noMCHits->in[0].y = x[1]; |
353 | noMCHits->in[0].z = x[2]; |
354 | noMCHits->in[0].E = pin[3]; |
355 | noMCHits->in[0].px = pin[0]*pin[4]; |
356 | noMCHits->in[0].py = pin[1]*pin[4]; |
357 | noMCHits->in[0].pz = pin[2]*pin[4]; |
358 | noMCHits->in[0].ptype = ipart; |
359 | noMCHits->in[0].itrack = track; |
360 | noMCHits->in[0].dist = dist; |
361 | noMCHits->mult = 1; |
362 | |
363 | } else { |
364 | fprintf(stderrstderr,"HDGeant error in hitForwardTOF (file hitFTOF.c): "); |
365 | fprintf(stderrstderr,"max hit count %d exceeded, truncating!\n",MAX_HITS); |
366 | } |
367 | } |
368 | |
369 | if (southHits != HDDM_NULL(void*)&hddm_s_nullTarget) { |
370 | |
371 | // loop over hits in this PM to find correct time slot |
372 | |
373 | for (nhit = 0; nhit < southHits->mult; nhit++) { |
374 | |
375 | if (fabs(southHits->in[nhit].t - t) < TWO_HIT_RESOL) { |
376 | break; |
377 | } |
378 | } |
379 | |
380 | // this hit is within the time frame of a previous hit |
381 | // combine the times of this weighted by the energy of the hit |
382 | |
383 | if (nhit < southHits->mult) { /* merge with former hit */ |
384 | southHits->in[nhit].t = |
385 | (southHits->in[nhit].t * southHits->in[nhit].dE |
386 | + tsouth * dEsouth) |
387 | / (southHits->in[nhit].dE += dEsouth); |
388 | |
389 | soMCHits = southHits->in[nhit].ftofMCHits; |
390 | |
391 | // now add MC tracking information |
392 | unsigned int nMChit = soMCHits->mult; |
393 | if (nMChit<MAX_PAD_HITS) { |
394 | |
395 | soMCHits->in[nMChit].x = x[0]; |
396 | soMCHits->in[nMChit].y = x[1]; |
397 | soMCHits->in[nMChit].z = x[2]; |
398 | soMCHits->in[nMChit].E = pin[3]; |
399 | soMCHits->in[nMChit].px = pin[0]*pin[4]; |
400 | soMCHits->in[nMChit].py = pin[1]*pin[4]; |
401 | soMCHits->in[nMChit].pz = pin[2]*pin[4]; |
402 | soMCHits->in[nMChit].ptype = ipart; |
403 | soMCHits->in[nMChit].itrack = track; |
404 | soMCHits->in[nMChit].dist = dist; |
405 | soMCHits->mult++; |
406 | } |
407 | |
408 | } else if (nhit < MAX_HITS){ // hit in new time window |
409 | southHits->in[nhit].t = tsouth; |
410 | southHits->in[nhit].dE = dEsouth; |
411 | southHits->mult++; |
412 | |
413 | // create memory space for MC track hit information |
414 | southHits->in[nhit].ftofMCHits = soMCHits = make_s_FtofMCHits(MAX_PAD_HITS); |
415 | |
416 | soMCHits->in[0].x = x[0]; |
417 | soMCHits->in[0].y = x[1]; |
418 | soMCHits->in[0].z = x[2]; |
419 | soMCHits->in[0].E = pin[3]; |
420 | soMCHits->in[0].px = pin[0]*pin[4]; |
421 | soMCHits->in[0].py = pin[1]*pin[4]; |
422 | soMCHits->in[0].pz = pin[2]*pin[4]; |
423 | soMCHits->in[0].ptype = ipart; |
424 | soMCHits->in[0].itrack = track; |
425 | soMCHits->in[0].dist = dist; |
426 | soMCHits->mult = 1; |
427 | |
428 | } else { |
429 | fprintf(stderrstderr,"HDGeant error in hitForwardTOF (file hitFTOF.c): "); |
430 | fprintf(stderrstderr,"max hit count %d exceeded, truncating!\n",MAX_HITS); |
431 | } |
432 | } |
433 | } |
434 | } |
435 | |
436 | /* entry point from fortran */ |
437 | |
438 | void hitforwardtof_ (float* xin, float* xout, |
439 | float* pin, float* pout, float* dEsum, |
440 | int* track, int* stack, int* history, int* ipart) |
441 | { |
442 | hitForwardTOF(xin,xout,pin,pout,*dEsum,*track,*stack,*history, *ipart); |
443 | } |
444 | |
445 | |
446 | /* pick and package the hits for shipping */ |
447 | // this function is called by loadoutput() (coded in hddmOutput.c) |
448 | // which in turn is called by GUOUT at the end of each event |
449 | |
450 | s_ForwardTOF_t* pickForwardTOF () |
451 | { |
452 | s_ForwardTOF_t* box; |
453 | s_ForwardTOF_t* item; |
454 | |
455 | if ((counterCount == 0) && (pointCount == 0)) |
456 | { |
457 | return HDDM_NULL(void*)&hddm_s_nullTarget; |
458 | } |
459 | |
460 | box = make_s_ForwardTOF(); |
461 | box->ftofCounters = make_s_FtofCounters(counterCount); |
462 | box->ftofTruthPoints = make_s_FtofTruthPoints(pointCount); |
463 | |
464 | while (item = (s_ForwardTOF_t*) pickTwig(&forwardTOFTree)){ |
465 | s_FtofCounters_t* counters = item->ftofCounters; |
466 | int counter; |
467 | s_FtofTruthPoints_t* points = item->ftofTruthPoints; |
468 | int point; |
469 | |
470 | for (counter=0; counter < counters->mult; ++counter) { |
471 | s_FtofNorthTruthHits_t* northHits = counters->in[counter].ftofNorthTruthHits; |
472 | s_FtofSouthTruthHits_t* southHits = counters->in[counter].ftofSouthTruthHits; |
473 | |
474 | /* compress out the hits below threshold */ |
475 | // cut off parameter is THRESH_MEV |
476 | int iok,i; |
477 | int mok=0; |
478 | // loop over all hits in a counter for the left/up PMT |
479 | for (iok=i=0; i < northHits->mult; i++) { |
480 | |
481 | // check threshold |
482 | if (northHits->in[i].dE >= THRESH_MEV/1e3) { |
483 | |
484 | if (iok < i) { |
485 | northHits->in[iok] = northHits->in[i]; |
486 | } |
487 | ++mok; |
488 | ++iok; |
489 | } |
490 | } |
491 | |
492 | if (iok) { |
493 | northHits->mult = iok; |
494 | } else if (northHits != HDDM_NULL(void*)&hddm_s_nullTarget){ // no hits left over for this PMT |
495 | counters->in[counter].ftofNorthHits = HDDM_NULL(void*)&hddm_s_nullTarget; |
496 | FREE(northHits)free(northHits); |
497 | } |
498 | |
499 | // now same loop for the right/bottom PMT of a paddle |
500 | for (iok=i=0; i < southHits->mult; i++){ |
501 | |
502 | if (southHits->in[i].dE >= THRESH_MEV/1e3){ |
503 | if (iok < i){ |
504 | southHits->in[iok] = southHits->in[i]; |
505 | } |
506 | ++mok; |
507 | ++iok; |
508 | } |
509 | } |
510 | |
511 | if (iok){ |
512 | southHits->mult = iok; |
513 | } |
514 | else if (southHits != HDDM_NULL(void*)&hddm_s_nullTarget) { |
515 | counters->in[counter].ftofSouthHits = HDDM_NULL(void*)&hddm_s_nullTarget; |
516 | FREE(southHits)free(southHits); |
517 | } |
518 | if (mok){ // total number of time independent FTOF hits in this counter |
519 | int m = box->ftofCounters->mult++; |
520 | // add the hit list of this counter to the list |
521 | box->ftofCounters->in[m] = counters->in[counter]; |
522 | } |
523 | } // end of loop over all counters |
524 | |
525 | if (counters != HDDM_NULL(void*)&hddm_s_nullTarget) { |
526 | FREE(counters)free(counters); |
527 | } |
528 | |
529 | // keep also the MC generated primary track particles |
530 | for (point=0; point < points->mult; ++point) { |
531 | int m = box->ftofTruthPoints->mult++; |
532 | box->ftofTruthPoints->in[m] = points->in[point]; |
533 | } |
534 | if (points != HDDM_NULL(void*)&hddm_s_nullTarget) { |
535 | FREE(points)free(points); |
536 | } |
537 | FREE(item)free(item); |
538 | } |
539 | |
540 | // reset the counters |
541 | counterCount = pointCount = 0; |
542 | |
543 | // free the hit list memory used by this event |
544 | if ((box->ftofCounters != HDDM_NULL(void*)&hddm_s_nullTarget) && |
545 | (box->ftofCounters->mult == 0)) { |
546 | FREE(box->ftofCounters)free(box->ftofCounters); |
547 | box->ftofCounters = HDDM_NULL(void*)&hddm_s_nullTarget; |
548 | } |
549 | if ((box->ftofTruthPoints != HDDM_NULL(void*)&hddm_s_nullTarget) && |
550 | (box->ftofTruthPoints->mult == 0)) { |
551 | FREE(box->ftofTruthPoints)free(box->ftofTruthPoints); |
552 | box->ftofTruthPoints = HDDM_NULL(void*)&hddm_s_nullTarget; |
553 | } |
554 | if ((box->ftofCounters->mult == 0) && |
555 | (box->ftofTruthPoints->mult == 0)) { |
556 | FREE(box)free(box); |
557 | box = HDDM_NULL(void*)&hddm_s_nullTarget; |
558 | } |
559 | return box; |
560 | } |
561 |