Bug Summary

File:alld2/home/sdobbs/Software/jana/jana_0.8.2/Linux_CentOS7.7-x86_64-gcc4.8.5/include/JANA/JEventLoop.h
Warning:line 398, column 7
Null pointer passed to 2nd parameter expecting 'nonnull'

Annotated Source Code

Press '?' to see keyboard shortcuts

clang -cc1 -cc1 -triple x86_64-unknown-linux-gnu -analyze -disable-free -main-file-name DBCALShower_factory_KLOE.cc -analyzer-store=region -analyzer-opt-analyze-nested-blocks -analyzer-checker=core -analyzer-checker=apiModeling -analyzer-checker=unix -analyzer-checker=deadcode -analyzer-checker=cplusplus -analyzer-checker=security.insecureAPI.UncheckedReturn -analyzer-checker=security.insecureAPI.getpw -analyzer-checker=security.insecureAPI.gets -analyzer-checker=security.insecureAPI.mktemp -analyzer-checker=security.insecureAPI.mkstemp -analyzer-checker=security.insecureAPI.vfork -analyzer-checker=nullability.NullPassedToNonnull -analyzer-checker=nullability.NullReturnedFromNonnull -analyzer-output plist -w -setup-static-analyzer -mrelocation-model pic -pic-level 2 -fhalf-no-semantic-interposition -mframe-pointer=none -fmath-errno -fno-rounding-math -mconstructor-aliases -munwind-tables -target-cpu x86-64 -tune-cpu generic -fno-split-dwarf-inlining -debugger-tuning=gdb -resource-dir /w/halld-scifs17exp/home/sdobbs/clang/llvm-project/install/lib/clang/12.0.0 -D HAVE_CCDB -D HAVE_RCDB -D HAVE_EVIO -D HAVE_TMVA=1 -D RCDB_MYSQL=1 -D RCDB_SQLITE=1 -D SQLITE_USE_LEGACY_STRUCT=ON -I .Linux_CentOS7.7-x86_64-gcc4.8.5/libraries/BCAL -I libraries/BCAL -I . -I libraries -I libraries/include -I /w/halld-scifs17exp/home/sdobbs/clang/halld_recon/Linux_CentOS7.7-x86_64-gcc4.8.5/include -I external/xstream/include -I /usr/include/tirpc -I /group/halld/Software/builds/Linux_CentOS7.7-x86_64-gcc4.8.5/root/root-6.08.06/include -I /w/halld-scifs17exp/halld2/home/sdobbs/Software/jana/jana_0.8.2/Linux_CentOS7.7-x86_64-gcc4.8.5/include -I /group/halld/Software/builds/Linux_CentOS7.7-x86_64-gcc4.8.5/ccdb/ccdb_1.06.06/include -I /group/halld/Software/builds/Linux_CentOS7.7-x86_64-gcc4.8.5/rcdb/rcdb_0.06.00/cpp/include -I /usr/include/mysql -I /group/halld/Software/builds/Linux_CentOS7.7-x86_64-gcc4.8.5/sqlitecpp/SQLiteCpp-2.2.0^bs130/include -I /group/halld/Software/builds/Linux_CentOS7.7-x86_64-gcc4.8.5/sqlite/sqlite-3.13.0^bs130/include -I /group/halld/Software/builds/Linux_CentOS7.7-x86_64-gcc4.8.5/hdds/hdds-4.9.0/Linux_CentOS7.7-x86_64-gcc4.8.5/src -I /group/halld/Software/builds/Linux_CentOS7.7-x86_64-gcc4.8.5/xerces-c/xerces-c-3.1.4/include -I /group/halld/Software/builds/Linux_CentOS7.7-x86_64-gcc4.8.5/evio/evio-4.4.6/Linux-x86_64/include -internal-isystem /usr/lib/gcc/x86_64-redhat-linux/4.8.5/../../../../include/c++/4.8.5 -internal-isystem /usr/lib/gcc/x86_64-redhat-linux/4.8.5/../../../../include/c++/4.8.5/x86_64-redhat-linux -internal-isystem /usr/lib/gcc/x86_64-redhat-linux/4.8.5/../../../../include/c++/4.8.5/backward -internal-isystem /usr/local/include -internal-isystem /w/halld-scifs17exp/home/sdobbs/clang/llvm-project/install/lib/clang/12.0.0/include -internal-externc-isystem /include -internal-externc-isystem /usr/include -O2 -std=c++11 -fdeprecated-macro -fdebug-compilation-dir /home/sdobbs/work/clang/halld_recon/src -ferror-limit 19 -fgnuc-version=4.2.1 -fcxx-exceptions -fexceptions -vectorize-loops -vectorize-slp -analyzer-output=html -faddrsig -o /tmp/scan-build-2021-01-21-110224-160369-1 -x c++ libraries/BCAL/DBCALShower_factory_KLOE.cc

libraries/BCAL/DBCALShower_factory_KLOE.cc

→
1// $Id$
2//
3// File: DBCALShower_factory_KLOE.cc
4// Created: Tue Jul 3 18:25:12 EDT 2007
5// Creator: Matthew Shepherd
6//
7
8#include <cassert>
9#include <math.h>
10#include <map>
11
12#include "BCAL/DBCALHit.h"
13#include "BCAL/DBCALTDCHit.h"
14#include "BCAL/DBCALPoint.h"
15#include "BCAL/DBCALGeometry.h"
16#include "BCAL/DBCALShower_factory_KLOE.h"
17
18#include "DANA/DApplication.h"
19
20#include "units.h"
21
22using namespace std;
23
24//------------------
25// init
26//------------------
27jerror_t DBCALShower_factory_KLOE::init()
28{
29 // this should be lower than cut in mcsmear
30 ethr_cell=0.0001; // MIN ENERGY THRESD OF cell in GeV
31
32 CLUST_THRESH = 0.03; // MIN ENERGY THRESD OF CLUSTER IN GEV (make this match the value used in the other algorithm)
33
34 elyr = 1;
35 xlyr = 2;
36 ylyr = 3;
37 zlyr = 4;
38 tlyr = 5;
39
40 // these four parameters are used in merging clusters
41 MERGE_THRESH_DIST = 40.0; // CENTROID DISTANCE THRESHOLD
42 MERGE_THRESH_TIME = 2.5; // CENTROID TIME THRESHOLD
43 MERGE_THRESH_ZDIST = 30.0; // FIBER DISTANCE THRESHOLD
44 MERGE_THRESH_XYDIST = 40.0; // CENTROID TRANSVERSE DISTANCE THRESHOLD
45
46 // this parameter is used to break clusters based on rms time
47 BREAK_THRESH_TRMS= 5.0; // T RMS THRESHOLD
48
49 gPARMS->SetDefaultParameter( "BCALRECON:CLUST_THRESH", CLUST_THRESH );
50 gPARMS->SetDefaultParameter( "BCALRECON:MERGE_THRESH_DIST", MERGE_THRESH_DIST );
51 gPARMS->SetDefaultParameter( "BCALRECON:MERGE_THRESH_TIME", MERGE_THRESH_TIME );
52 gPARMS->SetDefaultParameter( "BCALRECON:MERGE_THRESH_ZDIST", MERGE_THRESH_ZDIST );
53 gPARMS->SetDefaultParameter( "BCALRECON:MERGE_THRESH_XYDIST", MERGE_THRESH_XYDIST );
54 gPARMS->SetDefaultParameter( "BCALRECON:BREAK_THRESH_TRMS", BREAK_THRESH_TRMS );
55
56 /* this is unused code
57 if( !DBCALGeometry::summingOn() ){
58
59 // these are energy calibration parameters -- no summing of cells
60
61 m_scaleZ_p0 = 0.9597;
62 m_scaleZ_p1 = 0.000454875;
63 m_scaleZ_p2 = -2.29912e-06;
64 m_scaleZ_p3 = 1.49757e-09;
65
66 m_nonlinZ_p0 = -0.00154122;
67 m_nonlinZ_p1 = 6.73594e-05;
68 m_nonlinZ_p2 = 0;
69 m_nonlinZ_p3 = 0;
70
71 }
72 else{
73 */
74
75 // these are energy calibration parameters -- 1.2.3.4 summing
76
77 //last updated for svn revision 9233
78 m_scaleZ_p0 = 0.99798;
79 m_scaleZ_p1 = 0.000361096;
80 m_scaleZ_p2 = -2.17338e-06;
81 m_scaleZ_p3 = 1.32201e-09;
82
83 m_nonlinZ_p0 = -0.0201272;
84 m_nonlinZ_p1 = 0.000103649;
85 m_nonlinZ_p2 = 0;
86 m_nonlinZ_p3 = 0;
87 //}
88
89 return NOERROR;
90}
91
92//------------------
93// brun
94//------------------
95jerror_t DBCALShower_factory_KLOE::brun(JEventLoop *loop, int32_t runnumber)
96{
97 //get target position
98 DApplication* app = dynamic_cast<DApplication*>(loop->GetJApplication());
99 DGeometry* geom = app->GetDGeometry(runnumber);
100 geom->GetTargetZ(m_z_target_center);
101
102 vector<const DBCALGeometry*> bcalGeomVect;
103 loop->Get( bcalGeomVect );
104 bcalGeom = bcalGeomVect[0];
105
106 //////////////////////////////////////////////////////////////////
107 // Calculate Cell Position
108 //////////////////////////////////////////////////////////////////
109
110 ATTEN_LENGTH = bcalGeom->GetBCAL_attenutation_length();
111 C_EFFECTIVE = bcalGeom->GetBCAL_c_effective();
112
113 fiberLength = bcalGeom->GetBCAL_length(); // fiber length in cm
114 zOffset = bcalGeom->GetBCAL_center();
115
116 //the following uses some sad notation in which modmin=0 and modmax=48, when in fact there are 48 modules labelled either 0-47 or 1-48 depending on one's whim, although if we are using the methods from DBCALGeometry (e.g. cellId()), we must start counting from 1 and if we are accessing arrays we must of course start from 0
117 int modmin = 0;
118 int modmax = bcalGeom->GetBCAL_Nmodules();
119 int rowmin1=0;
120 int rowmax1= bcalGeom->GetBCAL_NInnerLayers();
121 int rowmin2= rowmax1;
122 int rowmax2= bcalGeom->GetBCAL_NOuterLayers()+rowmin2;
123 int colmin1=0;
124 int colmax1=bcalGeom->GetBCAL_NInnerSectors();
125 int colmin2=0;
126 int colmax2=bcalGeom->GetBCAL_NOuterSectors();
127
128 float r_inner= bcalGeom->GetBCAL_inner_rad();
129
130 for (int i = (rowmin1+1); i < (rowmax1+1); i++){
131 //this loop starts from 1, so we can use i in cellId with no adjustment
132 int cellId = bcalGeom->cellId(1,i,1); //this gives us a cellId that we can use to get the radius. the module and sector numbers are irrelevant
133 //rt is radius of center of layer - BCAL inner radius
134 rt[i]=bcalGeom->r(cellId)-r_inner;
135 }
136
137 for (int i = (rowmin2+1); i < (rowmax2+1); i++){
138 int cellId = bcalGeom->cellId(1,i,1);
139 rt[i]=bcalGeom->r(cellId)-r_inner;
140 }
141
142 //these are r and phi positions of readout cells
143 float r[modulemax_bcal48][layermax_bcal10][colmax_bcal4];
144 float phi[modulemax_bcal48][layermax_bcal10][colmax_bcal4];
145
146 // Now start to extract cell position information from Geometry class
147 for (int k = modmin; k < modmax; k++){
148 for (int i = rowmin1; i < rowmax1; i++){
149 for (int j = colmin1; j < colmax1; j++){
150 //in this case the loops start at 0, so we have to add 1 to the indices when calling cellId(). hooray!
151 //use DBCALGeometry to get r/phi position of each cell
152 int cellId = bcalGeom->cellId(k+1,i+1,j+1);
153 r[k][i][j]=bcalGeom->r(cellId);
154 phi[k][i][j]=bcalGeom->phi(cellId);
155 //set x and y positions
156 xx[k][i][j]=r[k][i][j]*cos(phi[k][i][j]);
157 yy[k][i][j]=r[k][i][j]*sin(phi[k][i][j]);
158 }
159 }
160
161 for (int i = rowmin2; i < rowmax2; i++){
162 for (int j = colmin2; j < colmax2; j++){
163 int cellId = bcalGeom->cellId(k+1,i+1,j+1);
164 r[k][i][j]=bcalGeom->r(cellId);
165 phi[k][i][j]=bcalGeom->phi(cellId);
166
167 xx[k][i][j]=r[k][i][j]*cos(phi[k][i][j]);
168 yy[k][i][j]=r[k][i][j]*sin(phi[k][i][j]);
169 }
170 }
171 }
172
173 ////////////////////////////////////////////////////////////////////////////
174 // Now the cell information are already contained in xx and yy arrays.
175 // xx and yy arrays are private members of this class
176 ////////////////////////////////////////////////////////////////////////////
177
178 return NOERROR;
179}
180
181//------------------
182// evnt
183//------------------
184jerror_t DBCALShower_factory_KLOE::evnt(JEventLoop *loop, uint64_t eventnumber)
185{
186 // Call core KLOE reconstruction routines
187 CellRecon(loop);
1
Calling 'DBCALShower_factory_KLOE::CellRecon'
→
188 CeleToArray();
189 PreCluster(loop);
190 ClusNorm();
191 ClusAnalysis();
192 Trakfit();
193
194 //Loop over reconstructed clusters and make DBCALShower objects out of them
195 vector<DBCALShower*> clusters;
196 int id = 0;
197 for (int i = 1; i < (clstot+1); i++){
198
199 int j=clspoi[i];
200
201 if( e_cls[j] < CLUST_THRESH ) continue;
202
203 // Time to cook a final shower
204 DBCALShower *shower = new DBCALShower;
205
206 //The algorithm has so far both clustered together hits and determined
207 //the position (x,y,z,t,E) of these clusters.
208 //For now, we will use only the clustering and recompute the position
209 //here. The reason we do this is because the clustering process is
210 //unaware of the errors on measurements of z (hits with TDC information
211 //will have much better z resolution). By taking into account this
212 //information we can greatly improve the z-resolution of clusters.
213
214 /*
215 //this is the old way of setting shower properties
216 shower->id = id++;
217 shower->E_raw = e_cls[j];
218 shower->x = x_cls[j];
219 shower->y = y_cls[j];
220 shower->z = z_cls[j] + zOffset;
221 shower->t = t_cls[j];
222 shower->N_cell = ncltot[j];
223
224 shower->xErr = eapx[1][j];
225 shower->yErr = eapx[2][j];
226 shower->zErr = eapx[3][j];
227
228 shower->tErr = 0.5 * sqrt( trms_a[j] * trms_a[j] +
229 trms_b[j] * trms_b[j] );
230 */
231
232 // Trace back to the DBCALPoint objects used in this shower and
233 // add them as associated objects.
234 vector<const DBCALPoint*> pointsInShower;
235 FindPointsInShower(j, loop, pointsInShower);
236
237 //Determine cluster (x,y,z,t) by averaging (x,y,z,t) of constituent
238 //DBCALPoints.
239
240 //For now just do the most naive averaging (weighting by E) to get
241 //the cluster properties (x,y,t). For z, average with weight of
242 //1/sig_z^2.
243 //Should consider a different weighting scheme (weighting by E^2) or average different quantities (cylindrical or spherical coordinates instead of rectangular)
244 double E=0,x=0,y=0,z=0,t=0;
245 int N_cell=0;
246 double sig_x=0,sig_y=0,sig_z=0,sig_t=0;
247 double sum_z_wt=0;
248 for(unsigned int j=0; j<pointsInShower.size(); j++){
249 double cell_E = pointsInShower[j]->E();
250 double cell_r = pointsInShower[j]->r();
251 double cell_phi = pointsInShower[j]->phi();
252 E += cell_E;
253 x += cell_E*cell_r*cos(cell_phi);
254 sig_x += cell_E*cell_r*cos(cell_phi)*cell_r*cos(cell_phi);
255 y += cell_E*cell_r*sin(cell_phi);
256 sig_y += cell_E*cell_r*sin(cell_phi)*cell_r*sin(cell_phi);
257
258
259 double z_wt = 1/(pointsInShower[j]->sigZ()*pointsInShower[j]->sigZ());
260 double cell_z = pointsInShower[j]->z();
261 double cell_t = pointsInShower[j]->t();
262
263 sum_z_wt += z_wt;
264 z += z_wt*cell_z;
265 sig_z += z_wt*cell_z*cell_z;
266
267 t += cell_E*cell_t;
268 sig_t += cell_E*cell_t*cell_t;
269 N_cell++;
270
271 shower->AddAssociatedObject(pointsInShower[j]);
272 }
273
274 x /= E;
275 sig_x /= E;
276 sig_x = sqrt(sig_x - x*x)/sqrt(N_cell); //really this should be n_effective rather than n, change this later
277 y /= E;
278 sig_y /= E;
279 sig_y = sqrt(sig_y - y*y)/sqrt(N_cell);
280 z /= sum_z_wt;
281 sig_z /= sum_z_wt;
282 sig_z = sqrt(sig_z - z*z)/sqrt(N_cell);
283 t /= E;
284 sig_t /= E;
285 sig_t = sqrt(sig_t - t*t)/sqrt(N_cell);
286
287 shower->id = id++;
288 shower->E_raw = E;
289 shower->x = x;
290 shower->y = y;
291 shower->z = z + m_z_target_center;
292 shower->t = t;
293 shower->N_cell = N_cell;
294
295 // shower->xErr = sig_x;
296 // shower->yErr = sig_y;
297 // shower->zErr = sig_z;
298
299 // shower->tErr = sig_t;
300
301 // calibrate energy:
302 // Energy calibration has a z dependence -- the
303 // calibration comes from fitting E_rec / E_gen to scale * E_gen^nonlin
304 // for slices of z. These fit parameters (scale and nonlin) are then plotted
305 // as a function of z and fit.
306
307 float r = sqrt( shower->x * shower->x + shower->y * shower->y );
308
309 float zEntry = ( shower->z - m_z_target_center ) * ( bcalGeom->GetBCAL_inner_rad() / r );
310
311 float scale = m_scaleZ_p0 + m_scaleZ_p1*zEntry +
312 m_scaleZ_p2*(zEntry*zEntry) + m_scaleZ_p3*(zEntry*zEntry*zEntry);
313 float nonlin = m_nonlinZ_p0 + m_nonlinZ_p1*zEntry +
314 m_nonlinZ_p2*(zEntry*zEntry) + m_nonlinZ_p3*(zEntry*zEntry*zEntry);
315
316 shower->E = pow( (shower->E_raw ) / scale, 1 / ( 1 + nonlin ) );
317
318 //copy xyz errors into covariance matrix
319 // shower->xyzCovariance.ResizeTo(3,3);
320 // shower->xyzCovariance[0][0] = shower->xErr*shower->xErr;
321 // shower->xyzCovariance[1][1] = shower->yErr*shower->yErr;
322 // shower->xyzCovariance[2][2] = shower->zErr*shower->zErr;
323
324 _data.push_back(shower);
325 }
326
327 return NOERROR;
328}
329
330
331//------------------
332// FindPointsInShower()
333//------------------
334void DBCALShower_factory_KLOE::FindPointsInShower(int indx, JEventLoop *loop, vector<const DBCALPoint*> &pointsInShower)
335{
336 /// This is called after the clusters have been completely formed. Our
337 /// job is simply to find the DBCALPoint objects used to form a given cluster.
338 /// This is so the DBCALPoint objects can be added to the DBCALShower object
339 /// as AssociatedObjects.
340
341 // The variable indx indexes the next[] array as the starting cell for the
342 // cluster. It also indexes the narr[][] array which holds the module, layer,
343 // sector(column) values for the hits.
344 //
345 // Here, we need to loop over next[] elements starting at next[indx] until
346 // we find the one pointing to element "indx" (i.e. the start of the list of
347 // cells in the cluster.) For each of these, we must find the LAST
348 // member in bcalhits to have the same module, layer, number indicating
349 // that is a DBCALHit to be added.
350
351
352 vector<const DBCALPoint*> points;
353 loop->Get(points);
354
355 int start_indx = indx;
356 do{
357 int module = narr[1][indx];
358 int layer = narr[2][indx];
359 int sector = narr[3][indx];
360
361 // Loop over BCAL hits, trying to find this one
362 for(unsigned int i=0; i<points.size(); i++){
363 if(points[i]->module() !=module)continue;
364 if(points[i]->layer() !=layer)continue;
365 if(points[i]->sector() !=sector)continue;
366 pointsInShower.push_back(points[i]);
367 }
368
369
370 indx = next[indx];
371 }while(indx != start_indx);
372
373}
374
375
376//------------------
377// CellRecon()
378//------------------
379void DBCALShower_factory_KLOE::CellRecon(JEventLoop *loop)
380{
381 //**********************************************************************
382 // The main purpose of this function is extracting information
383 // from DBCALPoint
384 // objects to form several arrays, which will be used later in the function
385 // CeleToArray()
386 // The four arrays ecel_a,tcel_a,ecel_b,tcel_b hold the energies and times
387 // of individual hits (this it the "attenuated" energy as measured at the
388 // end of the module, in GeV-equivalent units) (however times
389 // here should be *after* timewalk correction) (a=upstream, b=downstream).
390 // The five arrays xcel,ycel,zcel,tcel,ecel hold information about the
391 // position, time, and energy of an event in the calorimeter corresponding
392 // to a DBCALPoint (one hit at each end of the detector).
393 // The final two arrays tcell_anor and tcell_bnor contain the same
394 // information as tcel_a and tcell_b.
395 // These arrays are 3-D arrays and are rather logically are indexed by
396 // [module][layer][column]
397 //**********************************************************************
398
399 //First reset the arrays ecel_a,tcel_a,ecel_b,tcel_b to clear out garbage
400 //information from the previous events
401 memset( ecel_a, 0, modulemax_bcal48 * layermax_bcal10 *
402 colmax_bcal4 * sizeof( float ) );
403 memset( tcel_a, 0, modulemax_bcal48 * layermax_bcal10 *
404 colmax_bcal4 * sizeof( float ) );
405 memset( ecel_b, 0, modulemax_bcal48 * layermax_bcal10 *
406 colmax_bcal4 * sizeof( float ) );
407 memset( tcel_b, 0, modulemax_bcal48 * layermax_bcal10 *
408 colmax_bcal4 * sizeof( float ) );
409 //the other seven arrays will also be filled with garbage values from
410 //previous events, HOWEVER
411 //we don't need to zero out these arrays, as long as the ecel_a,tcel_a,ecel_b,tcel_b arrays are zeroed out
412 //this is because in CeleToArray(), all cells with ecel_a=0 are skipped
413 //and if ecel_a has been to set to a nonzero value for a particular cell
414 //then the other arrays will also have been set properly and not full of garbage
415
416 vector<const DBCALPoint*> points;
417 loop->Get(points);
2
←
Calling 'JEventLoop::Get'
→
418 if(points.size() <=0) return;
419
420 for (vector<const DBCALPoint*>::const_iterator point_iter = points.begin();
421 point_iter != points.end();
422 ++point_iter) {
423 const DBCALPoint &point = **point_iter;
424 int module = point.module();
425 int layer = point.layer();
426 int sector = point.sector();
427 double r = point.r();
428 double phi = point.phi();
429 double x = r*cos(phi);
430 double y = r*sin(phi);
431
432 xcel[module-1][layer-1][sector-1] = x;
433 ycel[module-1][layer-1][sector-1] = y;
434 //This factory expects z values relative to the center of the BCAL (z=212 cm)
435 //but DBCALPoint_factory gives us z values relative to the center of the target
436 zcel[module-1][layer-1][sector-1] = point.z()+m_z_target_center-zOffset;
437 tcel[module-1][layer-1][sector-1] = point.t();
438 ecel[module-1][layer-1][sector-1] = point.E();
439
440 //we require knowledge of the times and energies of the individual upstream and downstream hits
441 //to get these we need the associated objects
442 double EUp=0,EDown=0,tUp=0,tDown=0;
443 vector<const DBCALUnifiedHit*> assoc_hits;
444 point.Get(assoc_hits);
445 for (unsigned int i=0; i<assoc_hits.size(); i++) {
446 if (assoc_hits[i]->end == DBCALGeometry::kUpstream) {
447 EUp = assoc_hits[i]->E;
448 tUp = assoc_hits[i]->t;
449 }
450 if (assoc_hits[i]->end == DBCALGeometry::kDownstream) {
451 EDown = assoc_hits[i]->E;
452 tDown = assoc_hits[i]->t;
453 }
454
455 }
456
457 ecel_a[module-1][layer-1][sector-1] = EUp;
458 //for some reason we need to record the hit time in two different arrays
459 tcel_a[module-1][layer-1][sector-1] = tUp;
460 tcell_anor[module-1][layer-1][sector-1] = tUp;
461
462 ecel_b[module-1][layer-1][sector-1] = EDown;
463 tcel_b[module-1][layer-1][sector-1] = tDown;
464 tcell_bnor[module-1][layer-1][sector-1] = tDown;
465 }
466}
467
468
469//------------------
470// CeleToArray()
471//------------------
472void DBCALShower_factory_KLOE::CeleToArray(void)
473{
474 // THis code is adpapted from kloe code by Chuncheng Xu on June29,2005
475 // The following part is taken from kaloe's clurec_lib.f's subroutine
476 // cele_to_arr.
477 //
478 // It's original purpose is to extract the data information from array RW
479 // and IW into array NARR(j,i), CELDATA(j,i), nclus(i),next(i),
480 // and e_cel(i), x_cel(i),y_cel(i), z_cel(i), t_cel(i),ta_cel(i)
481 // and tb_cel(i)
482 // In the above explanationn, i is the index for cell, and different j
483 // is for different component for such cell.
484 // for example, NARR(1,i), NARR(2,i),NARR(3,i) are for module number,
485 // layer number and sector number in halld bcal simulation
486
487
488 // Now we assume that the information is stored in arrays ECEL_A(k,i,j)
489 // , ECEL_B(k,i,j), TCEL_A(k,i,j), TCEL_B(k,i,j), XCEL(k,I,J),YCEL(k,I,J),
490 // ZCEL(K,I,J), TCEL(K,I,J), ECEL(K,I,J),TCELL_ANOR(K,I,J),TCELL_BNOR(K,I,J),
491 // which are passed into here by event.inc through common block.
492
493 // these values e_cel(i), x_cel(i),y_cel(i), z_cel(i), t_cel(i),ta_cel(i)
494 // and tb_cel(i)) are extremly useful in later's clusterization
495 // and they are passed by common block to precluster subroutine too. (through
496 // clurec_cal.inc)
497 //
498 // we hope this is a good "bridge" from raw data to precluster.
499 // This way we can keep good use of their strategy of clusterization
500 // as much as possible, although we know that Fortran has it's bad fame
501 // of not so easy to be reused.
502
503
504 // Essentially this function takes the cell information (e.g. E,x,y,z,t)
505 // from the 3D arrays filled in CellRecon() and puts it into 1D arrays.
506 // These 1D arrays are indexed by an identifier. The module/layer/sector
507 // associated with this identifier can be found using the narr[][] array,
508 // as described above.
509
510 celtot=0;
511
512 for (int k = 0; k < modulemax_bcal48; k++){
513 for (int i = 0; i < layermax_bcal10; i++){
514 for (int j = 0; j < colmax_bcal4; j++){
515
516 float ea = ecel_a[k][i][j];
517 float eb = ecel_b[k][i][j];
518 float ta = tcel_a[k][i][j];
519 float tb = tcel_b[k][i][j];
520
521 if( (min(ea,eb)>ethr_cell) & (fabs(ta-tb)<35.) & (ta!=0.) & (tb!=0.)) {
522 celtot=celtot+1;
523 } else {
524 continue;
525 }
526
527
528 if(celtot>cellmax_bcal48*10*4) {
529 break;
530 }
531
532 narr[1][celtot]=k+1; // these numbers will
533 narr[2][celtot]=i+1; // be used by preclusters
534 narr[3][celtot]=j+1; // which will start from index of 1
535 // rather than from 0.
536
537 // why 0.145? -- these variables are used as weights
538 celdata[1][celtot]=ea/0.145;
539 celdata[2][celtot]=eb/0.145;
540
541 nclus[celtot] = celtot;
542 next[celtot] = celtot;
543
544 e_cel[celtot] = ecel[k][i][j];
545 x_cel[celtot] = xcel[k][i][j];
546 y_cel[celtot] = ycel[k][i][j];
547 z_cel[celtot] = zcel[k][i][j];
548 t_cel[celtot] = tcel[k][i][j];
549
550 ta_cel[celtot]=tcell_anor[k][i][j];
551 tb_cel[celtot]=tcell_bnor[k][i][j];
552 }
553 }
554 }
555}
556
557
558
559//------------------
560// PreCluster()
561//------------------
562void DBCALShower_factory_KLOE::PreCluster(JEventLoop *loop)
563{
564 //what this function does: for each cell with a hit it finds the maximum
565 //energy neighbor and Connect()'s the two. Two cells are neighbors if they
566 //are within one column of each other and within one row. The situation is
567 //slightly more complicated for two cells on opposite sides of the boundary
568 //between inner cells and outer and is described in more detail below, but
569 //essentially works out the same way. The purpose of Connect() is described
570 //in that function itself.
571
572 int k=1; // NUMBER OF NEARBY ROWS &/OR TO LOOK FOR MAX E CELL
573
574 // extract the BCAL Geometry
575 //vector<const DBCALGeometry*> bcalGeomVect;
576 //loop->Get( bcalGeomVect );
577 //const DBCALGeometry& bcalGeom = *(bcalGeomVect[0]);
578
579 // calculate cell position
580
581 int modmin = 0;
582 int modmax = bcalGeom->GetBCAL_Nmodules();
583
584
585 //these values make sense as actually minima/maxima if it is implied that rowmin1=1,colmin1=1,colmin2=1
586 int rowmax1= bcalGeom->GetBCAL_NInnerLayers();
587 int rowmin2= rowmax1+1;
588 //int rowmax2= bcalGeom->NBCALLAYSOUT+rowmin2-1;
589 int colmax1=bcalGeom->GetBCAL_NInnerSectors();
590 int colmax2=bcalGeom->GetBCAL_NOuterSectors();
591
592 float r_middle= bcalGeom->GetBCAL_middle_rad();
593
594 //radial size of the outermost inner layer
595 float thick_inner=bcalGeom->rSize(bcalGeom->cellId(1,bcalGeom->GetBCAL_NInnerLayers(),1));
596 //radial size of the innermost outer layer
597 float thick_outer=bcalGeom->rSize(bcalGeom->cellId(1,bcalGeom->GetBCAL_NInnerLayers()+1,1));
598
599 // this is the radial distance between the center of the innermost outer layer and the outermost inner layer
600 float dis_in_out=bcalGeom->r(bcalGeom->cellId(1,bcalGeom->GetBCAL_NInnerLayers()+1,1))-bcalGeom->r(bcalGeom->cellId(1,bcalGeom->GetBCAL_NInnerLayers(),1));
601
602 float degree_permodule=360.0/(modmax-modmin);
603 float half_degree_permodule=degree_permodule/2.0;
604
605 //roughly the width of a single cell in the outermost inner layer
606 float width_1=2.0*(r_middle-thick_inner/2.0)*
607 sin(half_degree_permodule*3.141593/180)/colmax1;
608 //roughly the width of a single cell in the innermost outer layer
609 float width_2=2.0*(r_middle+thick_outer/2.0)*
610 sin(half_degree_permodule*3.141593/180)/colmax2;
611
612 //disthres is roughly the azimuthal distance between the center of an outer cell and the center of the most distant inner cell bordering an adjacent outer cell (a picture would be nice wouldn't it)
613 //this value is used for determining if two cells that straddle the boundary between inner and outer layers should be considered as neighboring
614 float disthres=width_2*1.5-width_1*0.5+0.0001;
615
616 for (int i = 1; i < (celtot+1); i++){
617
618 int maxnn=0; //cell index of the maximum energy neighbor (if one is found)
619 float emin=0.; //energy of maximum energy neighbor
620
621
622 for (int j = 1; j < (celtot+1); j++){
623 if ( (j!=i) & (nclus[j]!=nclus[i]) & (e_cel[j]>emin)) {
624
625 int k1= narr[1][i];
626 int k2= narr[1][j];
627 int i1= narr[2][i];
628 int i2= narr[2][j];
629
630
631 int modiff = k1-k2;
632 int amodif = abs(modiff);
633
634 // the following if is to check module and row distance.
635 if ( (abs(i1-i2)<=k) & ((amodif<=1) || (amodif==47)) ) {
636 // further check col distance
637 int j1= narr[3][i];
638 int j2= narr[3][j];
639
640 if(amodif==0) { // same module
641 // further check col distance if both are inner layers
642 if ( (i1<=rowmax1) & (i2<=rowmax1) & (abs(j2-j1)<=k) ) {
643 emin=e_cel[j];
644 maxnn=j;
645 }
646
647 // further check col distance if both are outer layers
648
649 if ( (i1>=rowmin2) & (i2>=rowmin2) & (abs(j2-j1)<=k) ) {
650 emin=e_cel[j];
651 maxnn=j;
652 }
653 }
654
655 if(amodif>0) { // different module
656 if( (modiff==1) || (modiff==-47) ) {
657 if ( (i1<=rowmax1) & (i2<=rowmax1) ){
658 if(abs((j1+colmax1)-j2)<=k){
659 emin=e_cel[j];
660 maxnn=j;
661 }
662 }
663
664 if ( (i1>=rowmin2) & (i2>=rowmin2) ) {
665 if(abs((j1+colmax2)-j2)<=k){
666 emin=e_cel[j];
667 maxnn=j;
668 }
669 }
670 }
671
672 if ( (modiff==-1) || (modiff==47) ) {
673
674 if ( (i1<=rowmax1) & (i2<=rowmax1) ){
675 if(abs((j2+colmax1)-j1)<=k){
676 emin=e_cel[j];
677 maxnn=j;
678 }
679 }
680
681 if ( (i1>=rowmin2) & (i2>=rowmin2) ){
682 if(abs((j2+colmax2)-j1)<=k){
683 emin=e_cel[j];
684 maxnn=j;
685 }
686 }
687 }
688 }
689
690 // further check col distance if one is inner layer, another is outer
691 // so that the two may be between the boundary of two different size
692 // of cells.
693 if( ( (i1 == rowmax1) & (i2 == rowmin2) ) ||
694 ( (i1 == rowmin2) & (i2 == rowmax1) ) ) {
695
696 float delta_xx=xx[k1-1][i1-1][j1-1]-xx[k2-1][i2-1][j2-1];
697 float delta_yy=yy[k1-1][i1-1][j1-1]-yy[k2-1][i2-1][j2-1];
698
699 //distance between centers of two cells
700 float dis = sqrt( delta_xx * delta_xx + delta_yy * delta_yy );
701 //dis_in_out is the distance in radial direction, so we now isolate distance in direction perpendicular to radius
702 dis = sqrt( dis*dis - dis_in_out * dis_in_out );
703 //disthres is described above
704 if( dis < disthres ){
705 emin = e_cel[j];
706 maxnn = j;
707 }
708 }
709 }
710 }
711 } // finish second loop
712
713 if(maxnn>0){
714
715 Connect(maxnn,i);
716 }
717 } // finish first loop
718}
719
720
721
722//------------------
723// Connect(n,m);
724//------------------
725void DBCALShower_factory_KLOE::Connect(int n,int m)
726{
727 //----------------------------------------------------------------------
728 // Purpose and Methods :CONNECTS CELL M TO THE NEAREST MAX NEIGHBOR N
729 // Created 31-JUL-1992 WON KIM
730 //----------------------------------------------------------------------
731
732 // This little piece of code wasn't so easy to decipher, but I *think* I
733 // understand what it's doing. The idea is the following:
734 //
735 // (Prior to entering this routine)
736 // The sparsified list of hits is copied into some 1-D arrays with
737 // dimension cellmax_bcal+1. This includes the nclus[] and next[]
738 // arrays. The nclus[] array keeps the cluster number which is initalized
739 // to the sparsified hit number. Thus, every (double-ended) hit cell
740 // is it's own cluster. A loop over pairs of hits is done above to find
741 // nearest neighbors that should be merged into the same cluster.
742 //
743 // The next[] array contains an index to the "next" cell in the cluster
744 // for the given hit cell. This is a circular list such that if one follows
745 // the "next[next[next[...]]]" values, the last element points back to
746 // the first element of the cluster. As such, these are also initialized
747 // to index themselves as all single hits are considered a cluster of
748 // one element prior to calling this Connect() routine.
749 //
750 // When we are called, the value of "m" will be less than than value
751 // of "n". The cluster number (kept in nclus[]) is therefore updated
752 // for all members of nclus[m] to be the same as nclus[n]. Furthermore,
753 // the hit cell "m" is appended to the cluster nclus[n]. This also means
754 // that the cluster numbers become non-sequential since all elements of
755 // nclus whose value is "m" will be changed such that their values are
756 // "n".
757 //
758 // 5/10/2011 DL
759
760 if(nclus[n]!=nclus[m]){
761 int j=m;
762 nclus[j]=nclus[n];
763 while(next[j]!=m){
764 j=next[j];
765 nclus[j]=nclus[n];
766 }
767 next[j]=next[n];
768 next[n]=m;
769 }
770}
771
772
773//------------------
774// ClusNorm()
775//------------------
776void DBCALShower_factory_KLOE::ClusNorm(void)
777{
778 // fast initialization of arrays:
779 memset( e_cls, 0, ( clsmax_bcal48*10*4 + 1 ) * sizeof( float ) );
780 memset( x_cls, 0, ( clsmax_bcal48*10*4 + 1 ) * sizeof( float ) );
781 memset( y_cls, 0, ( clsmax_bcal48*10*4 + 1 ) * sizeof( float ) );
782 memset( z_cls, 0, ( clsmax_bcal48*10*4 + 1 ) * sizeof( float ) );
783 memset( t_cls, 0, ( clsmax_bcal48*10*4 + 1 ) * sizeof( float ) );
784 memset( ea_cls, 0, ( clsmax_bcal48*10*4 + 1 ) * sizeof( float ) );
785 memset( eb_cls, 0, ( clsmax_bcal48*10*4 + 1 ) * sizeof( float ) );
786 memset( ta_cls, 0, ( clsmax_bcal48*10*4 + 1 ) * sizeof( float ) );
787 memset( tb_cls, 0, ( clsmax_bcal48*10*4 + 1 ) * sizeof( float ) );
788 memset( tsqr_a, 0, ( clsmax_bcal48*10*4 + 1 ) * sizeof( float ) );
789 memset( tsqr_b, 0, ( clsmax_bcal48*10*4 + 1 ) * sizeof( float ) );
790 memset( trms_a, 0, ( clsmax_bcal48*10*4 + 1 ) * sizeof( float ) );
791 memset( trms_b, 0, ( clsmax_bcal48*10*4 + 1 ) * sizeof( float ) );
792 memset( e2_a, 0, ( clsmax_bcal48*10*4 + 1 ) * sizeof( float ) );
793 memset( e2_b, 0, ( clsmax_bcal48*10*4 + 1 ) * sizeof( float ) );
794 memset( clspoi, 0, ( clsmax_bcal48*10*4 + 1 ) * sizeof( float ) );
795 memset( ncltot, 0, ( clsmax_bcal48*10*4 + 1 ) * sizeof( float ) );
796 memset( ntopol, 0, ( clsmax_bcal48*10*4 + 1 ) * sizeof( float ) );
797
798 // Part of what is being done here is to further sparsify the
799 // data into a list of clusters. This starts to fill arrays
800 // with dimension clsmax_bcal+1. One important thing is that
801 // the clspoi[] array is being filled with the cluster number
802 // of what should be the index of the first cell hit in the
803 // cluster.
804 //
805 // 5/10/2011 DL
806
807 clstot=0;
808
809 for (int ix = 1; ix < (celtot+1); ix++){
810
811 //----------------------------------------------------------------------
812 // KEEP TALLY OF CLUSTERS
813 //----------------------------------------------------------------------
814
815 int n=nclus[ix];
816 int j=0;
817
818 for (int i = 1; i < (clstot+1); i++){
819 if(n==clspoi[i]) j=i;
820 }
821
822 if(j==0) {
823 clstot=clstot+1;
824 clspoi[clstot]=n;
825 }
826
827 //----------------------------------------------------------------------
828 // NORMALIZED QUANTITIES OF THE TOTAL CLUSTER
829 //----------------------------------------------------------------------
830 if(e_cel[ix]<0.000000001)continue; // just for protection
831
832 x_cls[n]=(e_cls[n]*x_cls[n]+e_cel[ix]*x_cel[ix])
833 /(e_cls[n]+e_cel[ix]);
834
835 y_cls[n]=(e_cls[n]*y_cls[n]+e_cel[ix]*y_cel[ix])
836 /(e_cls[n]+e_cel[ix]);
837
838 z_cls[n]=(e_cls[n]*z_cls[n]+e_cel[ix]*z_cel[ix])
839 /(e_cls[n]+e_cel[ix]);
840
841 t_cls[n]=(e_cls[n]*t_cls[n]+e_cel[ix]*t_cel[ix])
842 /(e_cls[n]+e_cel[ix]);
843
844 e_cls[n]=e_cls[n]+e_cel[ix];
845
846 // write(*,*)' x-y-z-T-e done'
847 //----------------------------------------------------------------------
848 // NORMALIZED QUANTITIES FOR EACH SIDE
849 //----------------------------------------------------------------------
850
851 ta_cls[n]=(ea_cls[n]*ta_cls[n]+celdata[1][ix]*ta_cel[ix])
852 /(ea_cls[n]+celdata[1][ix]);
853 tsqr_a[n]=(ea_cls[n]*tsqr_a[n]+celdata[1][ix]*ta_cel[ix]*
854 ta_cel[ix])/(ea_cls[n]+celdata[1][ix]);
855 ea_cls[n]=ea_cls[n]+celdata[1][ix];
856 e2_a[n]=e2_a[n]+celdata[1][ix]*celdata[1][ix];
857 tb_cls[n]=(eb_cls[n]*tb_cls[n]+celdata[2][ix]*tb_cel[ix])/
858 (eb_cls[n]+celdata[2][ix]);
859 tsqr_b[n]=(eb_cls[n]*tsqr_b[n]+celdata[2][ix]*tb_cel[ix]*
860 tb_cel[ix])/(eb_cls[n]+celdata[2][ix]);
861 eb_cls[n]=eb_cls[n]+celdata[2][ix];
862 e2_b[n]=e2_b[n]+celdata[2][ix]*celdata[2][ix];
863
864 //----------------------------------------------------------------------
865 // TOTAL CELLS AND CELLS IN OTHER MODULES
866 //----------------------------------------------------------------------
867
868 ncltot[n]++;
869
870 if( narr[1][n] != narr[1][ix] || narr[2][n] != narr[2][ix] )
871 ntopol[n]++;
872
873 //----------------------------------------------------------------------
874
875 }
876
877 for (int n = 1; n < (clstot+1); n++){
878
879 int ix=clspoi[n];
880 if( ncltot[ix] > 1) {
881
882 float effnum = ea_cls[ix] * ea_cls[ix] / e2_a[ix];
883 trms_a[ix] = effnum / ( effnum - 1 ) *
884 ( tsqr_a[ix] - ta_cls[ix] * ta_cls[ix] );
885
886 effnum = eb_cls[ix] * eb_cls[ix] / e2_b[ix];
887 trms_b[ix] = effnum / ( effnum - 1 ) *
888 ( tsqr_b[ix] - tb_cls[ix] * tb_cls[ix] );
889
890 if( trms_a[ix] <= 0.0 ) trms_a[ix] = 0.;
891 if( trms_b[ix] <= 0.0 ) trms_b[ix] = 0.;
892 trms_a[ix] = sqrt( trms_a[ix] );
893 trms_b[ix] = sqrt( trms_b[ix] );
894 }
895 else {
896 trms_a[ix] = 0.;
897 trms_b[ix] = 0.;
898 }
899 }
900}
901
902//------------------
903// ClusAnalysis()
904//------------------
905void DBCALShower_factory_KLOE::ClusAnalysis()
906{
907 // track when clusters change to cut down
908 // on excess calles to expensive ClusNorm
909 bool newClust = false;
910
911 //----------------------------------------------------------------------
912 // Check for overlapping clusters
913 //----------------------------------------------------------------------
914
915 for (int i = 0; i < 2; i++){
916 for (int j = 1; j < (clstot+1); j++){
917 int ix=clspoi[j];
918 if(e_cls[ix]>0.0){
919 float dist=sqrt(trms_a[ix]*trms_a[ix]+trms_b[ix]*trms_b[ix]);
920 if(dist>BREAK_THRESH_TRMS) {
921 Clus_Break(ix);
922 newClust = true;
923 }
924 }
925 }
926
927 if( newClust ){
928
929 ClusNorm();
930 newClust = false;
931 }
932 }
933
934 //----------------------------------------------------------------------
935 // merge clusters likely to be from the same shower
936 //----------------------------------------------------------------------
937
938 int icls[3];
939 for (int i = 1; i < clstot; i++){
940 icls[1]=0;
941 icls[2]=0;
942 for (int j = (i+1); j < (clstot+1); j++){
943
944 int ix=clspoi[i];
945 int iy=clspoi[j];
946
947
948 if ( (e_cls[ix]>0.0) & (e_cls[iy]>0.0) ) {
949
950 float delta_x=x_cls[ix]-x_cls[iy];
951 float delta_y=y_cls[ix]-y_cls[iy];
952 float delta_z=z_cls[ix]-z_cls[iy];
953 float dist=sqrt(delta_x*delta_x+delta_y*delta_y+delta_z*delta_z);
954
955 float tdif=fabs(t_cls[ix]-t_cls[iy]);
956
957 // float distz=abs(z_cls[ix]-z_cls[iy]);
958 // cout<<"dist="<<dist<<" distz="<<distz<<" tdif="<<tdif<<"\n";
959
960 if ( (dist<MERGE_THRESH_DIST) & (tdif<MERGE_THRESH_TIME) ){
961 float zdif=fabs(z_cls[ix]-z_cls[iy]);
962 float distran=sqrt(delta_x*delta_x+delta_y*delta_y);
963
964 if ( (zdif<MERGE_THRESH_ZDIST) & (distran<MERGE_THRESH_XYDIST) ){
965 if(e_cls[ix]>=e_cls[iy]) {
966 icls[1]=ix;
967 icls[2]=iy;
968 }
969 else {
970 icls[1]=iy;
971 icls[2]=ix;
972 }
973 }
974 }
975
976 }
977
978 if(min(icls[1],icls[2])>0){
979
980 Connect(icls[1],icls[2]);
981 newClust = true;
982 }
983 }
984 }
985
986 if( newClust ){
987
988 ClusNorm();
989 }
990}
991
992//------------------
993// Clus_Break(ix);
994//------------------
995void DBCALShower_factory_KLOE::Clus_Break(int nclust)
996{
997 int nseed[5],selnum,selcel[cellmax_bcal48*10*4+1];
998 float tdif,tdif_a,tdif_b,tseed[5];
999
1000 //----------------------------------------------------------------------
1001 for (int i =0; i < 5; i++){
1002 nseed[i]=0;
1003 tseed[i]=0;
1004 }
1005 //----------------------------------------------------------------------
1006 // Divide cluster cells into four quadrant groups
1007 //----------------------------------------------------------------------
1008 int n=nclust;
1009 tdif_a=ta_cel[n]-ta_cls[nclust];
1010 tdif_b=tb_cel[n]-tb_cls[nclust];
1011// selnum=0;
1012
1013 //----------------------------------------------------------------------
1014 if(tdif_a>0.0) {
1015 if(tdif_b>0){
1016 selnum=1;
1017 }
1018 else {
1019 selnum=2;
1020 }
1021 }
1022
1023 else {
1024 if(tdif_b>0.0) {
1025 selnum=3;
1026 }
1027 else {
1028 selnum=4;
1029 }
1030 }
1031
1032 //------------------------------------------------------------------------
1033
1034 if(selnum>0) {
1035 float tdif=sqrt(tdif_a*tdif_a+tdif_b*tdif_b);
1036 if(tdif>tseed[selnum]){
1037 nseed[selnum]=n;
1038 tseed[selnum]=tdif;
1039 }
1040 selcel[n]=selnum;
1041 }
1042
1043 //-------------------------------------------------------------------------
1044
1045 while(next[n]!=nclust) {
1046 n=next[n];
1047 tdif_a=ta_cel[n]-ta_cls[nclust];
1048 tdif_b=tb_cel[n]-tb_cls[nclust];
1049// selnum=0;
1050
1051 //**************************************************************************
1052
1053 if(tdif_a>0.0) {
1054
1055 if(tdif_b>0.0) {
1056 selnum=1;
1057 }
1058 else {
1059 selnum=2;
1060 }
1061 }
1062
1063 else {
1064 if(tdif_b>0.0) {
1065 selnum=3;
1066 }
1067 else {
1068 selnum=4;
1069 }
1070 }
1071
1072 //***************************************************************************
1073
1074 //...........................................................................
1075
1076
1077 if(selnum>0){
1078 tdif=sqrt(tdif_a*tdif_a+tdif_b*tdif_b);
1079
1080 if(tdif>tseed[selnum]){
1081 nseed[selnum]=n;
1082 tseed[selnum]=tdif;
1083 }
1084
1085 selcel[n]=selnum;
1086 }
1087
1088 //...........................................................................
1089
1090 } //this bracket is related to "while" above.
1091
1092
1093 //----------------------------------------------------------------------
1094 // If successful, divide cluster chain into the new cluster chains
1095 //----------------------------------------------------------------------
1096
1097 for (int i =1; i < 5; i++){
1098
1099 if(nseed[i]>0) {
1100
1101 nclus[nseed[i]]=nseed[i];
1102 next[nseed[i]]=nseed[i];
1103
1104 for (int j =1; j < (celtot+1); j++){
1105 if ( (nclus[j]==nclust) & (j!=nseed[i]) ){
1106 if(selcel[j]==i) {
1107 nclus[j]=j;
1108 next[j]=j;
1109 Connect(nseed[i],j);
1110 }
1111 }
1112 }
1113 }
1114 }
1115}
1116
1117
1118//------------------
1119// Trakfit()
1120//------------------
1121
1122// routine modified by MRS to do expensive filling of layer information
1123// contains code previously in ClusNorm which is called multiple times
1124// per event, but there is no need to call Trakfit multiple times per event
1125void DBCALShower_factory_KLOE::Trakfit( void )
1126{
1127
1128 float emin=0.0001;
1129
1130 memset( clslyr, 0, ( clsmax_bcal48*10*4 + 1 ) *
1131 ( layermax_bcal10 + 1 ) * 6 * sizeof( float ) );
1132
1133 memset( apx, 0, ( clsmax_bcal48*10*4 + 1 ) * 4 * sizeof( float ) );
1134 memset( eapx, 0, ( clsmax_bcal48*10*4 + 1 ) * 4 * sizeof( float ) );
1135 memset( ctrk, 0, ( clsmax_bcal48*10*4 + 1 ) * 4 * sizeof( float ) );
1136 memset( ectrk, 0, ( clsmax_bcal48*10*4 + 1 ) * 4 * sizeof( float ) );
1137
1138 for (int ix = 1; ix < (celtot+1); ix++){
1139
1140 int n = nclus[ix];
1141
1142 //----------------------------------------------------------------------
1143 // NORMALIZED QUANTITIES OF THE TOTAL CLUSTER PER LAYER
1144 //----------------------------------------------------------------------
1145
1146 int lyr = narr[2][ix];
1147
1148 // write(*,*)'X, E',E_CEL(ix),CLSLYR(XLYR,LYR,N)
1149 clslyr[xlyr][lyr][n]=(clslyr[elyr][lyr][n]*clslyr[xlyr][lyr][n]
1150 +e_cel[ix]*x_cel[ix])/(clslyr[elyr][lyr][n]+e_cel[ix]);
1151
1152 clslyr[ylyr][lyr][n]=(clslyr[elyr][lyr][n]*clslyr[ylyr][lyr][n]
1153 +e_cel[ix]*y_cel[ix])/(clslyr[elyr][lyr][n]+e_cel[ix]);
1154 // write(*,*)' Y'
1155
1156 clslyr[zlyr][lyr][n]=(clslyr[elyr][lyr][n]*clslyr[zlyr][lyr][n]
1157 +e_cel[ix]*z_cel[ix])/(clslyr[elyr][lyr][n]+e_cel[ix]);
1158
1159 // write(*,*)' z'
1160 clslyr[tlyr][lyr][n]=(clslyr[elyr][lyr][n]*clslyr[tlyr][lyr][n]
1161 +e_cel[ix]*t_cel[ix])/(clslyr[elyr][lyr][n]+e_cel[ix]);
1162
1163 // write(*,*)'E'
1164 clslyr[elyr][lyr][n]=clslyr[elyr][lyr][n]+e_cel[ix];
1165
1166 // write(*,*)' slopes done'
1167 }
1168
1169 memset( nlrtot, 0, ( clsmax_bcal48*10*4 + 1 ) * sizeof( float ) );
1170
1171 for (int n = 1; n < ( clstot + 1 ); n++){
1172
1173 int ix=clspoi[n];
1174
1175 for (int i = 1; i < (layermax_bcal10+1); i++){
1176
1177 if( clslyr[elyr][i][ix] > 0.0 ) nlrtot[ix]++;
1178 }
1179
1180 for (int i = 0; i < ( layermax_bcal10 + 1 ); i++){
1181
1182 x[i]=0.0;
1183 y[i]=0.0;
1184 z[i]=0.0;
1185 e[i]=0.0;
1186 sigx[i]=0.0; // cm
1187 sigy[i]=0.0; // cm
1188 sigz[i]=0.0;
1189 }
1190
1191 int nltot=0;
1192
1193 for (int il = 1; il < (layermax_bcal10+1); il++){
1194
1195 if(clslyr[elyr][il][ix]>emin) {
1196
1197 nltot=nltot+1;
1198 x[nltot]= clslyr[xlyr][il][ix];
1199 y[nltot]= clslyr[ylyr][il][ix];
1200 z[nltot]= clslyr[zlyr][il][ix];
1201 e[nltot]= clslyr[elyr][il][ix];
1202
1203 sigy[nltot] = 1.0/e[nltot];
1204 sigx[nltot] = 1.0/e[nltot];
1205 sigz[nltot] = 1.0/sqrt(e[nltot]);
1206 }
1207 }
1208
1209 // The following error bar is the estimation of error bar
1210 // based on the experience of my fortran code running
1211 // If the structure of BCAL changes drastically, you have to make changes
1212 // accordingly. Xu Chuncheng , Jan 9, 2006
1213
1214 // **** this seems strange since it appears to overwrite what is above
1215
1216 sigx[1]=0.5;
1217 sigx[2]=0.5;
1218 sigx[3]=0.5;
1219 sigx[4]=0.5;
1220 sigx[5]=0.5;
1221 sigx[6]=0.8;
1222 sigx[7]=0.9;
1223 sigx[8]=1.2;
1224 sigx[9]=1.3;
1225
1226 sigy[1]=0.5;
1227 sigy[2]=0.5;
1228 sigy[3]=0.5;
1229 sigy[4]=0.5;
1230 sigy[5]=0.5;
1231 sigy[6]=0.8;
1232 sigy[7]=0.9;
1233 sigy[8]=1.2;
1234 sigy[9]=1.3;
1235
1236
1237 sigz[1]=0.5;
1238 sigz[2]=0.5;
1239 sigz[3]=0.5;
1240 sigz[4]=0.5;
1241 sigz[5]=0.5;
1242 sigz[6]=0.8;
1243 sigz[7]=0.9;
1244 sigz[8]=1.2;
1245 sigz[9]=1.3;
1246
1247 if( nltot > 1 ){
1248
1249 Fit_ls();
1250
1251 for (int i = 1; i < 4; i++){
1252
1253 ctrk[i][ix]=ctrk_ix[i];
1254 ectrk[i][ix]=ectrk_ix[i];
1255 apx[i][ix]=apx_ix[i];
1256 eapx[i][ix]=eapx_ix[i];
1257 }
1258 }
1259 else{
1260
1261 // we have 1 or less layers hit -- no fit
1262
1263 apx[1][ix] = x[1];
1264 apx[2][ix] = y[1];
1265 apx[3][ix] = z[1];
1266 eapx[1][ix] = sigx[1];
1267 eapx[2][ix] = sigy[1];
1268 eapx[3][ix] = sigz[1];
1269 ectrk[1][ix] = 0.0;
1270 ectrk[2][ix] = 0.0;
1271 ectrk[3][ix] = 0.0;
1272 ctrk[1][ix] = 999.0;
1273 ctrk[2][ix] = 999.0;
1274 ctrk[3][ix] = 999.0;
1275 }
1276 }
1277}
1278
1279//------------------
1280// Fit_ls()
1281//------------------
1282void DBCALShower_factory_KLOE::Fit_ls()
1283{
1284 float a,b,c;
1285 float d,e,f,chi2,q,norm;
1286 float siga,sigb,sigc,sigd,sige,sigf;
1287 float sigb2,sigd2,sigf2;
1288
1289 // fitting for X=a+bt
1290 Linefit(1,1,a,b,siga,sigb,chi2,q);
1291 // fitting for Y=c+dt
1292 Linefit(2,1,c,d,sigc,sigd,chi2,q);
1293 // fitting for Z=e+ft
1294 Linefit(3,1,e,f,sige,sigf,chi2,q);
1295 sigb2=sigb*sigb;
1296 sigd2=sigd*sigd;
1297 sigf2=sigf*sigf;
1298
1299 apx_ix[1]=a;
1300 apx_ix[2]=c;
1301 apx_ix[3]=e;
1302 eapx_ix[1]=siga;
1303 eapx_ix[2]=sigc;
1304 eapx_ix[3]=sige;
1305
1306 // write(*,*) "b,d,f = ",b,d,f
1307 // cout<<"b= "<<b<<" d="<<d<<" f="<<f<<"\n";
1308
1309 // write(*,*) "sigb,sigd,sigf = ",sigb,sigd,sigf
1310
1311 norm=sqrt(b*b+d*d+f*f);
1312
1313 ctrk_ix[1]=b/norm;
1314 ctrk_ix[2]=d/norm;
1315 ctrk_ix[3]=f/norm;
1316
1317 float norm3=norm*norm*norm;
1318
1319 ectrk_ix[1]=sqrt((d*d+f*f)*(d*d+f*f)*sigb2+b*b*d*d*sigd2+b*b*f*f*sigf2)/norm3;
1320 ectrk_ix[2]=sqrt((b*b+f*f)*(b*b+f*f)*sigd2+d*d*b*b*sigb2+d*d*f*f*sigf2)/norm3;
1321 ectrk_ix[3]=sqrt((b*b+d*d)*(b*b+d*d)*sigf2+f*f*b*b*sigb2+f*f*d*d*sigd2)/norm3;
1322
1323 return;
1324}
1325
1326
1327//------------------
1328// Linefit(x,y,ndata,sig,mwt,a,b,siga,sigb,chi2,q)
1329//------------------
1330void DBCALShower_factory_KLOE::Linefit(int ixyz,int mwt,float &a,
1331 float &b,float &siga,float &sigb,float &chi2,float &q)
1332{
1333
1334 // This programme is taken from Garth's book
1335 // "Numerical Recipes in Fortran"
1336 // and we tested it. It works well. Chuncheng Xu,2005
1337
1338
1339 float sig[layermax_bcal10+1],etemp;
1340 float xtemp[layermax_bcal10+1],ytemp[layermax_bcal10+1];
1341 // uses gammq
1342
1343 // Given a set of data points X(1:ndata),Y(1:ndata) with individual standard
1344 // deviations sig(1:ndata), fit them to a strait line y=a+bx by minimum
1345 // chisq. Returned are a,b and their respective probable uncertainties
1346 // siga and sigb, the chisq, and the goodness-of-fit probability q(that the
1347 // fit would have chisq this large or larger). if mwt=0 on input, then the
1348 // the standard deviations are assumed to be unavalable:q is returned as
1349 // 1.0 and the normalization of chi2 is to the unit standard deviation on all
1350 // points.
1351
1352
1353 float sigdat,ss,st2,sx,sxoss,sy,t,wt;
1354 sx=0.0; // Initialize sums to zero
1355 sy=0.0;
1356 st2=0.0;
1357 b=0.0;
1358
1359 int ndata=0;
1360
1361
1362
1363 if(ixyz==1) {
1364 for (int i = 1; i < (layermax_bcal10+1); i++){
1365 xtemp[i]=rt[i];
1366 ytemp[i]=x[i];
1367 sig[i]=sigx[i];
1368 etemp=e[i];
1369 if(etemp>0.0001)ndata=ndata+1;
1370
1371 }
1372 }
1373 else if(ixyz==2) {
1374 for (int i = 1; i < (layermax_bcal10+1); i++){
1375 xtemp[i]=rt[i];
1376 ytemp[i]=y[i];
1377 sig[i]=sigy[i];
1378 etemp=e[i];
1379 if(etemp>0.000001)ndata=ndata+1;
1380 }
1381 }
1382 else if(ixyz==3) {
1383 for (unsigned int i = 1; i < (layermax_bcal10+1); i++){
1384 xtemp[i]=rt[i];
1385 ytemp[i]=z[i];
1386 sig[i]=sigz[i];
1387 etemp=e[i];
1388 if(etemp>0.000001)ndata=ndata+1;
1389 }
1390 }
1391
1392 if(mwt!=0) { // Accumulate sums
1393 ss=0.0;
1394 for (int i = 1; i < (ndata+1); i++){
1395 wt=1.0/(sig[i]*sig[i]);
1396 ss=ss+wt;
1397 sx=sx+xtemp[i]*wt;
1398 sy=sy+ytemp[i]*wt;
1399 }
1400 }
1401
1402 else {
1403 for (int i = 1; i < (ndata+1); i++){
1404 sx=sx+xtemp[i];
1405 sy=sy+ytemp[i];
1406 }
1407 ss=float(ndata);
1408 }
1409
1410 sxoss=sx/ss;
1411
1412 if(mwt!=0) {
1413 for (int i = 1; i < (ndata+1); i++){
1414 t=(xtemp[i]-sxoss)/sig[i];
1415 st2=st2+t*t;
1416 b=b+t*ytemp[i]/sig[i];
1417 }
1418
1419 }
1420
1421 else {
1422
1423 for (int i = 1; i < (ndata+1); i++){
1424 t=xtemp[i]-sxoss;
1425 st2=st2+t*t;
1426 b=b+t*ytemp[i];
1427 }
1428 }
1429
1430 b=b/st2; // Solve for a,b,siga and sigb
1431 a=(sy-sx*b)/ss;
1432 siga=sqrt((1.0+sx*sx/(ss*st2))/ss);
1433 sigb=sqrt(1.0/st2);
1434 chi2=0.0; // Calculate chisq
1435
1436 if(mwt==0) {
1437 for (int i = 1; i < (ndata+1); i++){
1438 chi2=chi2+(ytemp[i]-a-b*xtemp[i])*(ytemp[i]-a-b*xtemp[i]);
1439 }
1440 q=1.0;
1441 sigdat=sqrt(chi2/(ndata-2));
1442 siga=siga*sigdat;
1443 sigb=sigb*sigdat;
1444 }
1445 else {
1446 for (int i = 1; i < (ndata+1); i++){
1447 chi2=chi2+((ytemp[i]-a-b*xtemp[i])/
1448 sig[i])*((ytemp[i]-a-b*xtemp[i])/sig[i]);
1449 }
1450 q=Gammq(0.5*(ndata-2),0.5*chi2);
1451 }
1452
1453}
1454
1455
1456
1457//------------------
1458// Gammq(a_gammq,x_gammq)
1459//------------------
1460float DBCALShower_factory_KLOE::Gammq(float a_gammq,float x_gammq)
1461{
1462
1463 //============================================================================
1464
1465
1466 float gammq;
1467
1468 // uses gcf,gser
1469 //Returns the incomplete gamma function Q(a_gammq,x_gammq)=1-P(a_gammq,x_gammq)
1470
1471
1472 float gammcf=0,gamser;
1473
1474 if(a_gammq==0.0) {
1475 gammq=999.0;
1476 return gammq;
1477 }
1478
1479
1480
1481 if(x_gammq<0. || a_gammq<= 0.0) {
1482 // cout<<"bad arguments in gammq"<<"\n";
1483 return 999.0;// pause 'bad arguments in gammq'
1484 }
1485
1486 if(x_gammq<(a_gammq+1.)) { // Use the series representation
1487 Gser(gamser,a_gammq,x_gammq);
1488 gammq=1.0-gamser; // and take its complement
1489 }
1490 else { // Use continued fraction representation
1491 Gcf(gammcf,a_gammq,x_gammq);
1492 gammq=gammcf;
1493 }
1494 return gammq;
1495}
1496
1497
1498//------------------
1499// Gser(gamser,a_gser,x_gser,gln)
1500//------------------
1501void DBCALShower_factory_KLOE::Gser(float &gamser,float a_gser,float x_gser)
1502{
1503 //============================================================================
1504 int itmax=100;
1505 float eps=3.0e-7;
1506 float gln;
1507
1508 // uses gammln
1509 // Returns the incomplete gamma function P(a,x) evaluated by its
1510 // series representation as gamser. Also returns ln(Gamma(a)) as gln
1511
1512 float ap, del,sum;
1513
1514 gln=Gammln(a_gser);
1515
1516 if(x_gser<=0.0) {
1517 if(x_gser<0.0) cout<<"x_gser<0 in gser"<<"\n";
1518 gamser=0.0;
1519 return;
1520 }
1521
1522 ap=a_gser;
1523 sum=1.0/a_gser;
1524 del=sum;
1525
1526
1527 for (int n = 1; n < (itmax+1); n++){
1528 ap=ap+1.0;
1529 del=del*x_gser/ap;
1530 sum=sum+del;
1531
1532 if(fabs(del)<fabs(sum)*eps) {
1533 gamser=sum*exp(-x_gser+a_gser*log(x_gser)-gln);
1534 return;
1535 }
1536
1537 }
1538
1539 // cout<< "a too large, ITMAX too small in gser"<<"\n";
1540 return;
1541}
1542
1543
1544//------------------
1545// Gcf(gammcf,a,x,gln)
1546//------------------
1547void DBCALShower_factory_KLOE::Gcf(float &gammcf,float a_gcf,float x_gcf)
1548{
1549 //========================================================================
1550
1551 int itmax=100;
1552 float eps=3.0e-7;
1553 float fpmin=1.0e-30;
1554
1555 float gln;
1556
1557 // uses gammln
1558 // Returns the incomplete gamma function Q(a,x) evaluated by
1559 // its continued fraction representation as gammcf. Also returns
1560 // ln(Gamma(a)) as gln
1561
1562 // Parameters: ITMAX is the maximum allowed number of iterations;
1563 // EPS is the relative accuracy; FPMIN is a number near the smallest
1564 // representable floating-point number.
1565
1566 float an,b,c,d,del,h;
1567
1568 gln=Gammln(a_gcf);
1569 b=x_gcf+1.0-a_gcf;
1570 c=1.0/fpmin;
1571 d=1.0/b;
1572 h=d;
1573
1574
1575 for (int i = 1; i < (itmax+1); i++){
1576 an=-i*(i-a_gcf);
1577 b=b+2.0;
1578 d=an*d+b;
1579 if(fabs(d)<fpmin)d=fpmin;
1580 c=b+an/c;
1581 if(fabs(c)<fpmin)c=fpmin;
1582 d=1.0/d;
1583 del=d*c;
1584 h=h*del;
1585 if(fabs(del-1.0)<eps) {
1586 gammcf=exp(-x_gcf+a_gcf*log(x_gcf)-gln)*h; // Put factors in front
1587 return;
1588 }
1589 // cout<< "a too large,ITMAX too small in gcf"<<"\n";
1590 return;
1591 }
1592} //???
1593
1594//------------------
1595// Gammln(xx)
1596//------------------
1597float DBCALShower_factory_KLOE::Gammln(float xx_gln)
1598{
1599 // Returns the value ln[Gamma(xx_gln)] for xx_gln>0.0
1600
1601 float ser,stp,tmp,x_gln,y_gln;
1602 float cof[7];
1603 float gammln;
1604
1605 // Internal arithmetic will be done in double precision, a nicety that
1606 // that you can omit if five-figure accuracy is good enoug
1607 // save cof,stp
1608 // data cof,stp/76.18009172947146d0,-86.50532032941677d0,
1609 // * 24.01409824083091d0,-1.231739572450155d0,.1208650973866179d-2,
1610 // * -.5395239384953d-5,2.5066282746310005d0/
1611
1612 stp=2.5066282746310005;
1613 cof[1]=76.18009172947146;
1614 cof[2]=-86.50532032941677;
1615 cof[3]=24.01409824083091;
1616 cof[4]=-1.231739572450155;
1617 cof[5]=.1208650973866179e-2;
1618 cof[6]=-.5395239384953e-5;
1619
1620 x_gln=xx_gln;
1621 y_gln=x_gln;
1622 tmp=x_gln+5.5;
1623 tmp=(x_gln+0.5)*log(tmp)-tmp;
1624
1625 ser=1.000000000190015;
1626
1627 for (int j = 1; j < 7; j++){
1628 y_gln=y_gln+1.0;
1629 ser=ser+cof[j]/y_gln;
1630 }
1631
1632
1633 gammln=tmp+log(stp*ser/x_gln);
1634 return gammln;
1635}
1636
←

/w/halld-scifs17exp/halld2/home/sdobbs/Software/jana/jana_0.8.2/Linux_CentOS7.7-x86_64-gcc4.8.5/include/JANA/JEventLoop.h

→
1// $Id: JEventLoop.h 1763 2006-05-10 14:29:25Z davidl $
2//
3// File: JEventLoop.h
4// Created: Wed Jun 8 12:30:51 EDT 2005
5// Creator: davidl (on Darwin wire129.jlab.org 7.8.0 powerpc)
6//
7
8#ifndef _JEventLoop_
9#define _JEventLoop_
10
11#include <sys/time.h>
12
13#include <vector>
14#include <list>
15#include <string>
16#include <utility>
17#include <typeinfo>
18#include <string.h>
19#include <map>
20#include <utility>
21using std::vector;
22using std::list;
23using std::string;
24using std::type_info;
25
26#include <JANA/jerror.h>
27#include <JANA/JObject.h>
28#include <JANA/JException.h>
29#include <JANA/JEvent.h>
30#include <JANA/JThread.h>
31#include <JANA/JFactory_base.h>
32#include <JANA/JCalibration.h>
33#include <JANA/JGeometry.h>
34#include <JANA/JResourceManager.h>
35#include <JANA/JStreamLog.h>
36
37// The following is here just so we can use ROOT's THtml class to generate documentation.
38#include "cint.h"
39
40
41// Place everything in JANA namespace
42namespace jana{
43
44
45template<class T> class JFactory;
46class JApplication;
47class JEventProcessor;
48
49
50class JEventLoop{
51 public:
52
53 friend class JApplication;
54
55 enum data_source_t{
56 DATA_NOT_AVAILABLE = 1,
57 DATA_FROM_CACHE,
58 DATA_FROM_SOURCE,
59 DATA_FROM_FACTORY
60 };
61
62 typedef struct{
63 string caller_name;
64 string caller_tag;
65 string callee_name;
66 string callee_tag;
67 double start_time;
68 double end_time;
69 data_source_t data_source;
70 }call_stack_t;
71
72 typedef struct{
73 const char* factory_name;
74 string tag;
75 const char* filename;
76 int line;
77 }error_call_stack_t;
78
79 JEventLoop(JApplication *app); ///< Constructor
80 virtual ~JEventLoop(); ////< Destructor
81 virtual const char* className(void){return static_className();}
82 static const char* static_className(void){return "JEventLoop";}
83
84 JApplication* GetJApplication(void) const {return app;} ///< Get pointer to the JApplication object
85 void RefreshProcessorListFromJApplication(void); ///< Re-copy the list of JEventProcessors from JApplication
86 virtual jerror_t AddFactory(JFactory_base* factory); ///< Add a factory
87 jerror_t RemoveFactory(JFactory_base* factory); ///< Remove a factory
88 JFactory_base* GetFactory(const string data_name, const char *tag="", bool allow_deftag=true); ///< Get a specific factory pointer
89 vector<JFactory_base*> GetFactories(void) const {return factories;} ///< Get all factory pointers
90 void GetFactoryNames(vector<string> &factorynames); ///< Get names of all factories in name:tag format
91 void GetFactoryNames(map<string,string> &factorynames); ///< Get names of all factories in map with key=name, value=tag
92 map<string,string> GetDefaultTags(void) const {return default_tags;}
93 jerror_t ClearFactories(void); ///< Reset all factories in preparation for next event.
94 jerror_t PrintFactories(int sparsify=0); ///< Print a list of all factories.
95 jerror_t Print(const string data_name, const char *tag=""); ///< Print the data of the given type
96
97 JCalibration* GetJCalibration();
98 template<class T> bool GetCalib(string namepath, map<string,T> &vals);
99 template<class T> bool GetCalib(string namepath, vector<T> &vals);
100 template<class T> bool GetCalib(string namepath, T &val);
101
102 JGeometry* GetJGeometry();
103 template<class T> bool GetGeom(string namepath, map<string,T> &vals);
104 template<class T> bool GetGeom(string namepath, T &val);
105
106 JResourceManager* GetJResourceManager(void);
107 string GetResource(string namepath);
108 template<class T> bool GetResource(string namepath, T vals, int event_number=0);
109
110 void Initialize(void); ///< Do initializations just before event processing starts
111 jerror_t Loop(void); ///< Loop over events
112 jerror_t OneEvent(uint64_t event_number); ///< Process a specific single event (if source supports it)
113 jerror_t OneEvent(void); ///< Process a single event
114 inline void Pause(void){pause = 1;} ///< Pause event processing
115 inline void Resume(void){pause = 0;} ///< Resume event processing
116 inline void Quit(void){quit = 1;} ///< Clean up and exit the event loop
117 inline bool GetQuit(void) const {return quit;}
118 void QuitProgram(void);
119
120 // Support for random access of events
121 bool HasRandomAccess(void);
122 void AddToEventQueue(uint64_t event_number){ next_events_to_process.push_back(event_number); }
123 void AddToEventQueue(list<uint64_t> &event_numbers) { next_events_to_process.insert(next_events_to_process.end(), event_numbers.begin(), event_numbers.end()); }
124 list<uint64_t> GetEventQueue(void){ return next_events_to_process; }
125 void ClearEventQueue(void){ next_events_to_process.clear(); }
126
127 template<class T> JFactory<T>* GetSingle(const T* &t, const char *tag="", bool exception_if_not_one=true); ///< Get pointer to first data object from (source or factory).
128 template<class T> JFactory<T>* Get(vector<const T*> &t, const char *tag="", bool allow_deftag=true); ///< Get data object pointers from (source or factory)
129 template<class T> JFactory<T>* GetFromFactory(vector<const T*> &t, const char *tag="", data_source_t &data_source=null_data_source, bool allow_deftag=true); ///< Get data object pointers from factory
130 template<class T> jerror_t GetFromSource(vector<const T*> &t, JFactory_base *factory=NULL__null); ///< Get data object pointers from source.
131 inline JEvent& GetJEvent(void){return event;} ///< Get pointer to the current JEvent object.
132 inline void SetJEvent(JEvent *event){this->event = *event;} ///< Set the JEvent pointer.
133 inline void SetAutoFree(int auto_free){this->auto_free = auto_free;} ///< Set the Auto-Free flag on/off
134 inline pthread_t GetPThreadID(void) const {return pthread_id;} ///< Get the pthread of the thread to which this JEventLoop belongs
135 double GetInstantaneousRate(void) const {return rate_instantaneous;} ///< Get the current event processing rate
136 double GetIntegratedRate(void) const {return rate_integrated;} ///< Get the current event processing rate
137 double GetLastEventProcessingTime(void) const {return delta_time_single;}
138 unsigned int GetNevents(void) const {return Nevents;}
139
140 inline bool CheckEventBoundary(uint64_t event_numberA, uint64_t event_numberB);
141
142 inline bool GetCallStackRecordingStatus(void){ return record_call_stack; }
143 inline void DisableCallStackRecording(void){ record_call_stack = false; }
144 inline void EnableCallStackRecording(void){ record_call_stack = true; }
145 inline void CallStackStart(JEventLoop::call_stack_t &cs, const string &caller_name, const string &caller_tag, const string callee_name, const string callee_tag);
146 inline void CallStackEnd(JEventLoop::call_stack_t &cs);
147 inline vector<call_stack_t> GetCallStack(void){return call_stack;} ///< Get the current factory call stack
148 inline void AddToCallStack(call_stack_t &cs){if(record_call_stack) call_stack.push_back(cs);} ///< Add specified item to call stack record but only if record_call_stack is true
149 inline void AddToErrorCallStack(error_call_stack_t &cs){error_call_stack.push_back(cs);} ///< Add layer to the factory call stack
150 inline vector<error_call_stack_t> GetErrorCallStack(void){return error_call_stack;} ///< Get the current factory error call stack
151 void PrintErrorCallStack(void); ///< Print the current factory call stack
152
153 const JObject* FindByID(JObject::oid_t id); ///< Find a data object by its identifier.
154 template<class T> const T* FindByID(JObject::oid_t id); ///< Find a data object by its type and identifier
155 JFactory_base* FindOwner(const JObject *t); ///< Find the factory that owns a data object by pointer
156 JFactory_base* FindOwner(JObject::oid_t id); ///< Find a factory that owns a data object by identifier
157
158 // User defined references
159 template<class T> void SetRef(T *t); ///< Add a user reference to this JEventLoop (must be a pointer)
160 template<class T> T* GetRef(void); ///< Get a user-defined reference of a specific type
161 template<class T> vector<T*> GetRefsT(void); ///< Get all user-defined refrences of a specicif type
162 vector<pair<const char*, void*> > GetRefs(void){ return user_refs; } ///< Get copy of full list of user-defined references
163 template<class T> void RemoveRef(T *t); ///< Remove user reference from list
164
165 // Convenience methods wrapping JEvent methods of same name
166 uint64_t GetStatus(void){return event.GetStatus();}
167 bool GetStatusBit(uint32_t bit){return event.GetStatusBit(bit);}
168 bool SetStatusBit(uint32_t bit, bool val=true){return event.SetStatusBit(bit, val);}
169 bool ClearStatusBit(uint32_t bit){return event.ClearStatusBit(bit);}
170 void ClearStatus(void){event.ClearStatus();}
171 void SetStatusBitDescription(uint32_t bit, string description){event.SetStatusBitDescription(bit, description);}
172 string GetStatusBitDescription(uint32_t bit){return event.GetStatusBitDescription(bit);}
173 void GetStatusBitDescriptions(map<uint32_t, string> &status_bit_descriptions){return event.GetStatusBitDescriptions(status_bit_descriptions);}
174 string StatusWordToString(void);
175
176 private:
177 JEvent event;
178 vector<JFactory_base*> factories;
179 vector<JEventProcessor*> processors;
180 vector<error_call_stack_t> error_call_stack;
181 vector<call_stack_t> call_stack;
182 JApplication *app;
183 JThread *jthread;
184 bool initialized;
185 bool print_parameters_called;
186 int pause;
187 int quit;
188 int auto_free;
189 pthread_t pthread_id;
190 map<string, string> default_tags;
191 vector<pair<string,string> > auto_activated_factories;
192 bool record_call_stack;
193 string caller_name;
194 string caller_tag;
195 vector<uint64_t> event_boundaries;
196 int32_t event_boundaries_run; ///< Run number boundaries were retrieved from (possbily 0)
197 list<uint64_t> next_events_to_process;
198
199 uint64_t Nevents; ///< Total events processed (this thread)
200 uint64_t Nevents_rate; ///< Num. events accumulated for "instantaneous" rate
201 double delta_time_single; ///< Time spent processing last event
202 double delta_time_rate; ///< Integrated time accumulated "instantaneous" rate (partial number of events)
203 double delta_time; ///< Total time spent processing events (this thread)
204 double rate_instantaneous; ///< Latest instantaneous rate
205 double rate_integrated; ///< Rate integrated over all events
206
207 static data_source_t null_data_source;
208
209 vector<pair<const char*, void*> > user_refs;
210};
211
212
213// The following is here just so we can use ROOT's THtml class to generate documentation.
214#ifdef G__DICTIONARY
215typedef JEventLoop::call_stack_t call_stack_t;
216typedef JEventLoop::error_call_stack_t error_call_stack_t;
217#endif
218
219#if !defined(__CINT__) && !defined(__CLING__)
220
221//-------------
222// GetSingle
223//-------------
224template<class T>
225JFactory<T>* JEventLoop::GetSingle(const T* &t, const char *tag, bool exception_if_not_one)
226{
227 /// This is a convenience method that can be used to get a pointer to the single
228 /// object of type T from the specified factory. It simply calls the Get(vector<...>) method
229 /// and copies the first pointer into "t" (or NULL if something other than 1 object is returned).
230 ///
231 /// This is intended to address the common situation in which there is an interest
232 /// in the event if and only if there is exactly 1 object of type T. If the event
233 /// has no objects of that type or more than 1 object of that type (for the specified
234 /// factory) then an exception of type "unsigned long" is thrown with the value
235 /// being the number of objects of type T. You can supress the exception by setting
236 /// exception_if_not_one to false. In that case, you will have to check if t==NULL to
237 /// know if the call succeeded.
238 vector<const T*> v;
239 JFactory<T> *fac = Get(v, tag);
240
241 if(v.size()!=1){
242 t = NULL__null;
243 if(exception_if_not_one) throw v.size();
244 }
245
246 t = v[0];
247
248 return fac;
249}
250
251//-------------
252// Get
253//-------------
254template<class T>
255JFactory<T>* JEventLoop::Get(vector<const T*> &t, const char *tag, bool allow_deftag)
256{
257 /// Retrieve or generate the array of objects of
258 /// type T for the curent event being processed
259 /// by this thread.
260 ///
261 /// By default, preference is given to reading the
262 /// objects from the data source(e.g. file) before generating
263 /// them in the factory. A flag exists in the factory
264 /// however to change this so that the factory is
265 /// given preference.
266 ///
267 /// Note that regardless of the setting of this flag,
268 /// the data are only either read in or generated once.
269 /// Ownership of the objects will always be with the
270 /// factory so subsequent calls will always return pointers to
271 /// the same data.
272 ///
273 /// If the factory is called on to generate the data,
274 /// it is done by calling the factory's Get() method
275 /// which, in turn, calls the evnt() method.
276 ///
277 /// First, we just call the GetFromFactory() method.
278 /// It will make the initial decision as to whether
279 /// it should look in the source first or not. If
280 /// it returns NULL, then the factory couldn't be
281 /// found so we automatically try the file.
282 ///
283 /// Note that if no factory exists to hold the objects
284 /// from the file, one can be created automatically
285 /// providing the <i>JANA:AUTOFACTORYCREATE</i>
286 /// configuration parameter is set.
287
288 // Check if a tag was specified for this data type to use for the
289 // default.
290 const char *mytag = tag
2.1
'tag' is not equal to NULL
2.1
'tag' is not equal to NULL
2.1
'tag' is not equal to NULL
==NULL__null ? "":tag; // protection against NULL tags
3
←
'?' condition is false
→
291 if(strlen(mytag)==0 && allow_deftag
3.1
'allow_deftag' is true
3.1
'allow_deftag' is true
3.1
'allow_deftag' is true
){
4
←
Taking true branch
→
292 map<string, string>::const_iterator iter = default_tags.find(T::static_className());
293 if(iter!=default_tags.end())tag = iter->second.c_str();
5
←
Assuming the condition is true
→
6
←
Taking true branch
→
7
←
Value assigned to 'tag'
→
294 }
295
296
297 // If we are trying to keep track of the call stack then we
298 // need to add a new call_stack_t object to the the list
299 // and initialize it with the start time and caller/callee
300 // info.
301 call_stack_t cs;
302
303 // Optionally record starting info of call stack entry
304 if(record_call_stack) CallStackStart(cs, caller_name, caller_tag, T::static_className(), tag);
8
←
Assuming field 'record_call_stack' is false
→
9
←
Taking false branch
→
305
306 // Get the data (or at least try to)
307 JFactory<T>* factory=NULL__null;
308 try{
309 factory = GetFromFactory(t, tag, cs.data_source, allow_deftag);
10
←
Passing value via 2nd parameter 'tag'
→
11
←
Calling 'JEventLoop::GetFromFactory'
→
310 if(!factory){
311 // No factory exists for this type and tag. It's possible
312 // that the source may be able to provide the objects
313 // but it will need a place to put them. We can create a
314 // dumb JFactory just to hold the data in case the source
315 // can provide the objects. Before we do though, make sure
316 // the user condones this via the presence of the
317 // "JANA:AUTOFACTORYCREATE" config parameter.
318 string p;
319 try{
320 gPARMS->GetParameter("JANA:AUTOFACTORYCREATE", p);
321 }catch(...){}
322 if(p.size()==0){
323 jout<<std::endl;
324 _DBG__std::cerr<<"/w/halld-scifs17exp/halld2/home/sdobbs/Software/jana/jana_0.8.2/Linux_CentOS7.7-x86_64-gcc4.8.5/include/JANA/JEventLoop.h"
<<":"<<324<<std::endl;
325 jout<<"No factory of type \""<<T::static_className()<<"\" with tag \""<<tag<<"\" exists."<<std::endl;
326 jout<<"If you are reading objects from a file, I can auto-create a factory"<<std::endl;
327 jout<<"of the appropriate type to hold the objects, but this feature is turned"<<std::endl;
328 jout<<"off by default. To turn it on, set the \"JANA:AUTOFACTORYCREATE\""<<std::endl;
329 jout<<"configuration parameter. This can usually be done by passing the"<<std::endl;
330 jout<<"following argument to the program from the command line:"<<std::endl;
331 jout<<std::endl;
332 jout<<" -PJANA:AUTOFACTORYCREATE=1"<<std::endl;
333 jout<<std::endl;
334 jout<<"Note that since the most commonly expected occurance of this situation."<<std::endl;
335 jout<<"is an error, the program will now throw an exception so that the factory."<<std::endl;
336 jout<<"call stack can be printed."<<std::endl;
337 jout<<std::endl;
338 throw exception();
339 }else{
340 AddFactory(new JFactory<T>(tag));
341 jout<<__FILE__"/w/halld-scifs17exp/halld2/home/sdobbs/Software/jana/jana_0.8.2/Linux_CentOS7.7-x86_64-gcc4.8.5/include/JANA/JEventLoop.h"<<":"<<__LINE__341<<" Auto-created "<<T::static_className()<<":"<<tag<<" factory"<<std::endl;
342
343 // Now try once more. The GetFromFactory method will call
344 // GetFromSource since it's empty.
345 factory = GetFromFactory(t, tag, cs.data_source, allow_deftag);
346 }
347 }
348 }catch(exception &e){
349 // Uh-oh, an exception was thrown. Add us to the call stack
350 // and re-throw the exception
351 error_call_stack_t ecs;
352 ecs.factory_name = T::static_className();
353 ecs.tag = tag;
354 ecs.filename = NULL__null;
355 error_call_stack.push_back(ecs);
356 throw e;
357 }
358
359 // If recording the call stack, update the end_time field and add to stack
360 if(record_call_stack) CallStackEnd(cs);
361
362 return factory;
363}
364
365//-------------
366// GetFromFactory
367//-------------
368template<class T>
369JFactory<T>* JEventLoop::GetFromFactory(vector<const T*> &t, const char *tag, data_source_t &data_source, bool allow_deftag)
370{
371 // We need to find the factory providing data type T with
372 // tag given by "tag".
373 vector<JFactory_base*>::iterator iter=factories.begin();
374 JFactory<T> *factory = NULL__null;
375 string className(T::static_className());
376
377 // Check if a tag was specified for this data type to use for the
378 // default.
379 const char *mytag = tag==NULL__null ? "":tag; // protection against NULL tags
12
←
Assuming 'tag' is equal to NULL
→
13
←
Assuming pointer value is null
→
14
←
'?' condition is true
→
380 if(strlen(mytag)==0 && allow_deftag
14.1
'allow_deftag' is true
14.1
'allow_deftag' is true
14.1
'allow_deftag' is true
){
15
←
Taking true branch
→
381 map<string, string>::const_iterator iter = default_tags.find(className);
382 if(iter!=default_tags.end())tag = iter->second.c_str();
16
←
Assuming the condition is false
→
17
←
Taking false branch
→
383 }
384
385 for(; iter!=factories.end(); iter++){
18
←
Calling 'operator!=<jana::JFactory_base **, std::vector<jana::JFactory_base *>>'
→
21
←
Returning from 'operator!=<jana::JFactory_base **, std::vector<jana::JFactory_base *>>'
→
22
←
Loop condition is true. Entering loop body
→
386 // It turns out a long standing bug in g++ makes dynamic_cast return
387 // zero improperly when used on objects created on one side of
388 // a dynamically shared object (DSO) and the cast occurs on the
389 // other side. I saw bug reports ranging from 2001 to 2004. I saw
390 // saw it first-hand on LinuxEL4 using g++ 3.4.5. This is too bad
391 // since it is much more elegant (and safe) to use dynamic_cast.
392 // To avoid this problem which can occur with plugins, we check
393 // the name of the data classes are the same. (sigh)
394 //factory = dynamic_cast<JFactory<T> *>(*iter);
395 if(className == (*iter)->GetDataClassName())factory = (JFactory<T>*)(*iter);
23
←
Taking true branch
→
396 if(factory == NULL__null)continue;
24
←
Assuming 'factory' is not equal to NULL
→
25
←
Taking false branch
→
397 const char *factag = factory->Tag()==NULL__null ? "":factory->Tag();
26
←
Assuming the condition is true
→
27
←
'?' condition is true
→
398 if(!strcmp(factag, tag)){
28
←
Null pointer passed to 2nd parameter expecting 'nonnull'
399 break;
400 }else{
401 factory=NULL__null;
402 }
403 }
404
405 // If factory not found, just return now
406 if(!factory){
407 data_source = DATA_NOT_AVAILABLE;
408 return NULL__null;
409 }
410
411 // OK, we found the factory. If the evnt() routine has already
412 // been called, then just call the factory's Get() routine
413 // to return a copy of the existing data
414 if(factory->evnt_was_called()){
415 factory->CopyFrom(t);
416 data_source = DATA_FROM_CACHE;
417 return factory;
418 }
419
420 // Next option is to get the objects from the data source
421 if(factory->GetCheckSourceFirst()){
422 // If the object type/tag is found in the source, it
423 // will return NOERROR, even if there are zero instances
424 // of it. If it is not available in the source then it
425 // will return OBJECT_NOT_AVAILABLE.
426
427 jerror_t err = GetFromSource(t, factory);
428 if(err == NOERROR){
429 // A return value of NOERROR means the source had the objects
430 // even if there were zero of them.(If the source had no
431 // information about the objects OBJECT_NOT_AVAILABLE would
432 // have been returned.)
433 // The GetFromSource() call will eventually lead to a call to
434 // the GetObjects() method of the concrete class derived
435 // from JEventSource. That routine should copy the object
436 // pointers into the factory using the factory's CopyTo()
437 // method which also sets the evnt_called flag for the factory.
438 // Note also that the "t" vector is then filled with a call
439 // to the factory's CopyFrom() method in JEvent::GetObjects().
440 // All we need to do now is just set the factory pointers in
441 // the newly generated JObjects and return the factory pointer.
442
443 factory->SetFactoryPointers();
444 data_source = DATA_FROM_SOURCE;
445
446 return factory;
447 }
448 }
449
450 // OK. It looks like we have to have the factory make this.
451 // Get pointers to data from the factory.
452 factory->Get(t);
453 factory->SetFactoryPointers();
454 data_source = DATA_FROM_FACTORY;
455
456 return factory;
457}
458
459//-------------
460// GetFromSource
461//-------------
462template<class T>
463jerror_t JEventLoop::GetFromSource(vector<const T*> &t, JFactory_base *factory)
464{
465 /// This tries to get objects from the event source.
466 /// "factory" must be a valid pointer to a JFactory
467 /// object since that will take ownership of the objects
468 /// created by the source.
469 /// This should usually be called from JEventLoop::GetFromFactory
470 /// which is called from JEventLoop::Get. The latter will
471 /// create a dummy JFactory of the proper flavor and tag if
472 /// one does not already exist so if objects exist in the
473 /// file without a corresponding factory to create them, they
474 /// can still be used.
475 if(!factory)throw OBJECT_NOT_AVAILABLE;
476
477 return event.GetObjects(t, factory);
478}
479
480//-------------
481// CallStackStart
482//-------------
483inline void JEventLoop::CallStackStart(JEventLoop::call_stack_t &cs, const string &caller_name, const string &caller_tag, const string callee_name, const string callee_tag)
484{
485 /// This is used to fill initial info into a call_stack_t stucture
486 /// for recording the call stack. It should be matched with a call
487 /// to CallStackEnd. It is normally called from the Get() method
488 /// above, but may also be used by external actors to manipulate
489 /// the call stack (presumably for good and not evil).
490
491 struct itimerval tmr;
492 getitimer(ITIMER_PROFITIMER_PROF, &tmr);
493
494 cs.caller_name = this->caller_name;
495 cs.caller_tag = this->caller_tag;
496 this->caller_name = cs.callee_name = callee_name;
497 this->caller_tag = cs.callee_tag = callee_tag;
498 cs.start_time = tmr.it_value.tv_sec + tmr.it_value.tv_usec/1.0E6;
499}
500
501//-------------
502// CallStackEnd
503//-------------
504inline void JEventLoop::CallStackEnd(JEventLoop::call_stack_t &cs)
505{
506 /// Complete a call stack entry. This should be matched
507 /// with a previous call to CallStackStart which was
508 /// used to fill the cs structure.
509
510 struct itimerval tmr;
511 getitimer(ITIMER_PROFITIMER_PROF, &tmr);
512 cs.end_time = tmr.it_value.tv_sec + tmr.it_value.tv_usec/1.0E6;
513 caller_name = cs.caller_name;
514 caller_tag = cs.caller_tag;
515 call_stack.push_back(cs);
516}
517
518//-------------
519// CheckEventBoundary
520//-------------
521inline bool JEventLoop::CheckEventBoundary(uint64_t event_numberA, uint64_t event_numberB)
522{
523 /// Check whether the two event numbers span one or more boundaries
524 /// in the calibration/conditions database for the current run number.
525 /// Return true if they do and false if they don't. The first parameter
526 /// "event_numberA" is also checked if it lands on a boundary in which
527 /// case true is also returned. If event_numberB lands on a boundary,
528 /// but event_numberA does not, then false is returned.
529 ///
530 /// This method is not expected to be called by a user. It is, however called,
531 /// everytime a JFactory's Get() method is called.
532
533 // Make sure our copy of the boundaries is up to date
534 if(event.GetRunNumber()!=event_boundaries_run){
535 event_boundaries.clear(); // in case we can't get the JCalibration pointer
536 JCalibration *jcalib = GetJCalibration();
537 if(jcalib)jcalib->GetEventBoundaries(event_boundaries);
538 event_boundaries_run = event.GetRunNumber();
539 }
540
541 // Loop over boundaries
542 for(unsigned int i=0; i<event_boundaries.size(); i++){
543 uint64_t eb = event_boundaries[i];
544 if((eb - event_numberA)*(eb - event_numberB) < 0.0 || eb==event_numberA){ // think about it ....
545 // events span a boundary or is on a boundary. Return true
546 return true;
547 }
548 }
549
550 return false;
551}
552
553//-------------
554// FindByID
555//-------------
556template<class T>
557const T* JEventLoop::FindByID(JObject::oid_t id)
558{
559 /// This is a templated method that can be used in place
560 /// of the non-templated FindByID(oid_t) method if one knows
561 /// the class of the object with the specified id.
562 /// This method is faster than calling the non-templated
563 /// FindByID and dynamic_cast-ing the JObject since
564 /// this will only search the objects of factories that
565 /// produce the desired data type.
566 /// This method will cast the JObject pointer to one
567 /// of the specified type. To use this method,
568 /// a type is specified in the call as follows:
569 ///
570 /// const DMyType *t = loop->FindByID<DMyType>(id);
571
572 // Loop over factories looking for ones that provide
573 // specified data type.
574 for(unsigned int i=0; i<factories.size(); i++){
575 if(factories[i]->GetDataClassName() != T::static_className())continue;
576
577 // This factory provides data of type T. Search it for
578 // the object with the specified id.
579 const JObject *my_obj = factories[i]->GetByID(id);
580 if(my_obj)return dynamic_cast<const T*>(my_obj);
581 }
582
583 return NULL__null;
584}
585
586//-------------
587// GetCalib (map)
588//-------------
589template<class T>
590bool JEventLoop::GetCalib(string namepath, map<string,T> &vals)
591{
592 /// Get the JCalibration object from JApplication for the run number of
593 /// the current event and call its Get() method to get the constants.
594
595 // Note that we could do this by making "vals" a generic type T thus, combining
596 // this with the vector version below. However, doing this explicitly will make
597 // it easier for the user to understand how to call us.
598
599 vals.clear();
600
601 JCalibration *calib = GetJCalibration();
602 if(!calib){
603 _DBG_std::cerr<<"/w/halld-scifs17exp/halld2/home/sdobbs/Software/jana/jana_0.8.2/Linux_CentOS7.7-x86_64-gcc4.8.5/include/JANA/JEventLoop.h"
<<":"<<603<<" "<<"Unable to get JCalibration object for run "<<event.GetRunNumber()<<std::endl;
604 return true;
605 }
606
607 return calib->Get(namepath, vals, event.GetEventNumber());
608}
609
610//-------------
611// GetCalib (vector)
612//-------------
613template<class T> bool JEventLoop::GetCalib(string namepath, vector<T> &vals)
614{
615 /// Get the JCalibration object from JApplication for the run number of
616 /// the current event and call its Get() method to get the constants.
617
618 vals.clear();
619
620 JCalibration *calib = GetJCalibration();
621 if(!calib){
622 _DBG_std::cerr<<"/w/halld-scifs17exp/halld2/home/sdobbs/Software/jana/jana_0.8.2/Linux_CentOS7.7-x86_64-gcc4.8.5/include/JANA/JEventLoop.h"
<<":"<<622<<" "<<"Unable to get JCalibration object for run "<<event.GetRunNumber()<<std::endl;
623 return true;
624 }
625
626 return calib->Get(namepath, vals, event.GetEventNumber());
627}
628
629//-------------
630// GetCalib (single)
631//-------------
632template<class T> bool JEventLoop::GetCalib(string namepath, T &val)
633{
634 /// This is a convenience method for getting a single entry. It
635 /// simply calls the vector version and returns the first entry.
636 /// It returns true if the vector version returns true AND there
637 /// is at least one entry in the vector. No check is made for there
638 /// there being more than one entry in the vector.
639
640 vector<T> vals;
641 bool ret = GetCalib(namepath, vals);
642 if(vals.empty()) return true;
643 val = vals[0];
644
645 return ret;
646}
647
648//-------------
649// GetGeom (map)
650//-------------
651template<class T>
652bool JEventLoop::GetGeom(string namepath, map<string,T> &vals)
653{
654 /// Get the JGeometry object from JApplication for the run number of
655 /// the current event and call its Get() method to get the constants.
656
657 // Note that we could do this by making "vals" a generic type T thus, combining
658 // this with the vector version below. However, doing this explicitly will make
659 // it easier for the user to understand how to call us.
660
661 vals.clear();
662
663 JGeometry *geom = GetJGeometry();
664 if(!geom){
665 _DBG_std::cerr<<"/w/halld-scifs17exp/halld2/home/sdobbs/Software/jana/jana_0.8.2/Linux_CentOS7.7-x86_64-gcc4.8.5/include/JANA/JEventLoop.h"
<<":"<<665<<" "<<"Unable to get JGeometry object for run "<<event.GetRunNumber()<<std::endl;
666 return true;
667 }
668
669 return geom->Get(namepath, vals);
670}
671
672//-------------
673// GetGeom (atomic)
674//-------------
675template<class T> bool JEventLoop::GetGeom(string namepath, T &val)
676{
677 /// Get the JCalibration object from JApplication for the run number of
678 /// the current event and call its Get() method to get the constants.
679
680 JGeometry *geom = GetJGeometry();
681 if(!geom){
682 _DBG_std::cerr<<"/w/halld-scifs17exp/halld2/home/sdobbs/Software/jana/jana_0.8.2/Linux_CentOS7.7-x86_64-gcc4.8.5/include/JANA/JEventLoop.h"
<<":"<<682<<" "<<"Unable to get JGeometry object for run "<<event.GetRunNumber()<<std::endl;
683 return true;
684 }
685
686 return geom->Get(namepath, val);
687}
688
689//-------------
690// SetRef
691//-------------
692template<class T>
693void JEventLoop::SetRef(T *t)
694{
695 pair<const char*, void*> p(typeid(T).name(), (void*)t);
696 user_refs.push_back(p);
697}
698
699//-------------
700// GetResource
701//-------------
702template<class T> bool JEventLoop::GetResource(string namepath, T vals, int event_number)
703{
704 JResourceManager *resource_manager = GetJResourceManager();
705 if(!resource_manager){
706 string mess = string("Unable to get the JResourceManager object (namepath=\"")+namepath+"\")";
707 throw JException(mess);
708 }
709
710 return resource_manager->Get(namepath, vals, event_number);
711}
712
713//-------------
714// GetRef
715//-------------
716template<class T>
717T* JEventLoop::GetRef(void)
718{
719 /// Get a user-defined reference (a pointer)
720 for(unsigned int i=0; i<user_refs.size(); i++){
721 if(user_refs[i].first == typeid(T).name()) return (T*)user_refs[i].second;
722 }
723
724 return NULL__null;
725}
726
727//-------------
728// GetRefsT
729//-------------
730template<class T>
731vector<T*> JEventLoop::GetRefsT(void)
732{
733 vector<T*> refs;
734 for(unsigned int i=0; i<user_refs.size(); i++){
735 if(user_refs[i].first == typeid(T).name()){
736 refs.push_back((T*)user_refs[i].second);
737 }
738 }
739
740 return refs;
741}
742
743//-------------
744// RemoveRef
745//-------------
746template<class T>
747void JEventLoop::RemoveRef(T *t)
748{
749 vector<pair<const char*, void*> >::iterator iter;
750 for(iter=user_refs.begin(); iter!= user_refs.end(); iter++){
751 if(iter->second == (void*)t){
752 user_refs.erase(iter);
753 return;
754 }
755 }
756 _DBG_std::cerr<<"/w/halld-scifs17exp/halld2/home/sdobbs/Software/jana/jana_0.8.2/Linux_CentOS7.7-x86_64-gcc4.8.5/include/JANA/JEventLoop.h"
<<":"<<756<<" "<<" Attempt to remove user reference not in event loop!" << std::endl;
757}
758
759
760#endif //__CINT__ __CLING__
761
762} // Close JANA namespace
763
764
765
766#endif // _JEventLoop_
767
←

/usr/lib/gcc/x86_64-redhat-linux/4.8.5/../../../../include/c++/4.8.5/bits/stl_iterator.h

1// Iterators -*- C++ -*-
2
3// Copyright (C) 2001-2013 Free Software Foundation, Inc.
4//
5// This file is part of the GNU ISO C++ Library. This library is free
6// software; you can redistribute it and/or modify it under the
7// terms of the GNU General Public License as published by the
8// Free Software Foundation; either version 3, or (at your option)
9// any later version.
10
11// This library is distributed in the hope that it will be useful,
12// but WITHOUT ANY WARRANTY; without even the implied warranty of
13// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14// GNU General Public License for more details.
15
16// Under Section 7 of GPL version 3, you are granted additional
17// permissions described in the GCC Runtime Library Exception, version
18// 3.1, as published by the Free Software Foundation.
19
20// You should have received a copy of the GNU General Public License and
21// a copy of the GCC Runtime Library Exception along with this program;
22// see the files COPYING3 and COPYING.RUNTIME respectively. If not, see
23// <http://www.gnu.org/licenses/>.
24
25/*
26 *
27 * Copyright (c) 1994
28 * Hewlett-Packard Company
29 *
30 * Permission to use, copy, modify, distribute and sell this software
31 * and its documentation for any purpose is hereby granted without fee,
32 * provided that the above copyright notice appear in all copies and
33 * that both that copyright notice and this permission notice appear
34 * in supporting documentation. Hewlett-Packard Company makes no
35 * representations about the suitability of this software for any
36 * purpose. It is provided "as is" without express or implied warranty.
37 *
38 *
39 * Copyright (c) 1996-1998
40 * Silicon Graphics Computer Systems, Inc.
41 *
42 * Permission to use, copy, modify, distribute and sell this software
43 * and its documentation for any purpose is hereby granted without fee,
44 * provided that the above copyright notice appear in all copies and
45 * that both that copyright notice and this permission notice appear
46 * in supporting documentation. Silicon Graphics makes no
47 * representations about the suitability of this software for any
48 * purpose. It is provided "as is" without express or implied warranty.
49 */
50
51/** @file bits/stl_iterator.h
52 * This is an internal header file, included by other library headers.
53 * Do not attempt to use it directly. @headername{iterator}
54 *
55 * This file implements reverse_iterator, back_insert_iterator,
56 * front_insert_iterator, insert_iterator, __normal_iterator, and their
57 * supporting functions and overloaded operators.
58 */
59
60#ifndef _STL_ITERATOR_H1
61#define _STL_ITERATOR_H1 1
62
63#include <bits/cpp_type_traits.h>
64#include <ext/type_traits.h>
65#include <bits/move.h>
66
67namespace std _GLIBCXX_VISIBILITY(default)__attribute__ ((__visibility__ ("default")))
68{
69_GLIBCXX_BEGIN_NAMESPACE_VERSION
70
71 /**
72 * @addtogroup iterators
73 * @{
74 */
75
76 // 24.4.1 Reverse iterators
77 /**
78 * Bidirectional and random access iterators have corresponding reverse
79 * %iterator adaptors that iterate through the data structure in the
80 * opposite direction. They have the same signatures as the corresponding
81 * iterators. The fundamental relation between a reverse %iterator and its
82 * corresponding %iterator @c i is established by the identity:
83 * @code
84 * &*(reverse_iterator(i)) == &*(i - 1)
85 * @endcode
86 *
87 * <em>This mapping is dictated by the fact that while there is always a
88 * pointer past the end of an array, there might not be a valid pointer
89 * before the beginning of an array.</em> [24.4.1]/1,2
90 *
91 * Reverse iterators can be tricky and surprising at first. Their
92 * semantics make sense, however, and the trickiness is a side effect of
93 * the requirement that the iterators must be safe.
94 */
95 template<typename _Iterator>
96 class reverse_iterator
97 : public iterator<typename iterator_traits<_Iterator>::iterator_category,
98 typename iterator_traits<_Iterator>::value_type,
99 typename iterator_traits<_Iterator>::difference_type,
100 typename iterator_traits<_Iterator>::pointer,
101 typename iterator_traits<_Iterator>::reference>
102 {
103 protected:
104 _Iterator current;
105
106 typedef iterator_traits<_Iterator> __traits_type;
107
108 public:
109 typedef _Iterator iterator_type;
110 typedef typename __traits_type::difference_type difference_type;
111 typedef typename __traits_type::pointer pointer;
112 typedef typename __traits_type::reference reference;
113
114 /**
115 * The default constructor value-initializes member @p current.
116 * If it is a pointer, that means it is zero-initialized.
117 */
118 // _GLIBCXX_RESOLVE_LIB_DEFECTS
119 // 235 No specification of default ctor for reverse_iterator
120 reverse_iterator() : current() { }
121
122 /**
123 * This %iterator will move in the opposite direction that @p x does.
124 */
125 explicit
126 reverse_iterator(iterator_type __x) : current(__x) { }
127
128 /**
129 * The copy constructor is normal.
130 */
131 reverse_iterator(const reverse_iterator& __x)
132 : current(__x.current) { }
133
134 /**
135 * A %reverse_iterator across other types can be copied if the
136 * underlying %iterator can be converted to the type of @c current.
137 */
138 template<typename _Iter>
139 reverse_iterator(const reverse_iterator<_Iter>& __x)
140 : current(__x.base()) { }
141
142 /**
143 * @return @c current, the %iterator used for underlying work.
144 */
145 iterator_type
146 base() const
147 { return current; }
148
149 /**
150 * @return A reference to the value at @c --current
151 *
152 * This requires that @c --current is dereferenceable.
153 *
154 * @warning This implementation requires that for an iterator of the
155 * underlying iterator type, @c x, a reference obtained by
156 * @c *x remains valid after @c x has been modified or
157 * destroyed. This is a bug: http://gcc.gnu.org/PR51823
158 */
159 reference
160 operator*() const
161 {
162 _Iterator __tmp = current;
163 return *--__tmp;
164 }
165
166 /**
167 * @return A pointer to the value at @c --current
168 *
169 * This requires that @c --current is dereferenceable.
170 */
171 pointer
172 operator->() const
173 { return &(operator*()); }
174
175 /**
176 * @return @c *this
177 *
178 * Decrements the underlying iterator.
179 */
180 reverse_iterator&
181 operator++()
182 {
183 --current;
184 return *this;
185 }
186
187 /**
188 * @return The original value of @c *this
189 *
190 * Decrements the underlying iterator.
191 */
192 reverse_iterator
193 operator++(int)
194 {
195 reverse_iterator __tmp = *this;
196 --current;
197 return __tmp;
198 }
199
200 /**
201 * @return @c *this
202 *
203 * Increments the underlying iterator.
204 */
205 reverse_iterator&
206 operator--()
207 {
208 ++current;
209 return *this;
210 }
211
212 /**
213 * @return A reverse_iterator with the previous value of @c *this
214 *
215 * Increments the underlying iterator.
216 */
217 reverse_iterator
218 operator--(int)
219 {
220 reverse_iterator __tmp = *this;
221 ++current;
222 return __tmp;
223 }
224
225 /**
226 * @return A reverse_iterator that refers to @c current - @a __n
227 *
228 * The underlying iterator must be a Random Access Iterator.
229 */
230 reverse_iterator
231 operator+(difference_type __n) const
232 { return reverse_iterator(current - __n); }
233
234 /**
235 * @return *this
236 *
237 * Moves the underlying iterator backwards @a __n steps.
238 * The underlying iterator must be a Random Access Iterator.
239 */
240 reverse_iterator&
241 operator+=(difference_type __n)
242 {
243 current -= __n;
244 return *this;
245 }
246
247 /**
248 * @return A reverse_iterator that refers to @c current - @a __n
249 *
250 * The underlying iterator must be a Random Access Iterator.
251 */
252 reverse_iterator
253 operator-(difference_type __n) const
254 { return reverse_iterator(current + __n); }
255
256 /**
257 * @return *this
258 *
259 * Moves the underlying iterator forwards @a __n steps.
260 * The underlying iterator must be a Random Access Iterator.
261 */
262 reverse_iterator&
263 operator-=(difference_type __n)
264 {
265 current += __n;
266 return *this;
267 }
268
269 /**
270 * @return The value at @c current - @a __n - 1
271 *
272 * The underlying iterator must be a Random Access Iterator.
273 */
274 reference
275 operator[](difference_type __n) const
276 { return *(*this + __n); }
277 };
278
279 //@{
280 /**
281 * @param __x A %reverse_iterator.
282 * @param __y A %reverse_iterator.
283 * @return A simple bool.
284 *
285 * Reverse iterators forward many operations to their underlying base()
286 * iterators. Others are implemented in terms of one another.
287 *
288 */
289 template<typename _Iterator>
290 inline bool
291 operator==(const reverse_iterator<_Iterator>& __x,
292 const reverse_iterator<_Iterator>& __y)
293 { return __x.base() == __y.base(); }
294
295 template<typename _Iterator>
296 inline bool
297 operator<(const reverse_iterator<_Iterator>& __x,
298 const reverse_iterator<_Iterator>& __y)
299 { return __y.base() < __x.base(); }
300
301 template<typename _Iterator>
302 inline bool
303 operator!=(const reverse_iterator<_Iterator>& __x,
304 const reverse_iterator<_Iterator>& __y)
305 { return !(__x == __y); }
306
307 template<typename _Iterator>
308 inline bool
309 operator>(const reverse_iterator<_Iterator>& __x,
310 const reverse_iterator<_Iterator>& __y)
311 { return __y < __x; }
312
313 template<typename _Iterator>
314 inline bool
315 operator<=(const reverse_iterator<_Iterator>& __x,
316 const reverse_iterator<_Iterator>& __y)
317 { return !(__y < __x); }
318
319 template<typename _Iterator>
320 inline bool
321 operator>=(const reverse_iterator<_Iterator>& __x,
322 const reverse_iterator<_Iterator>& __y)
323 { return !(__x < __y); }
324
325 template<typename _Iterator>
326 inline typename reverse_iterator<_Iterator>::difference_type
327 operator-(const reverse_iterator<_Iterator>& __x,
328 const reverse_iterator<_Iterator>& __y)
329 { return __y.base() - __x.base(); }
330
331 template<typename _Iterator>
332 inline reverse_iterator<_Iterator>
333 operator+(typename reverse_iterator<_Iterator>::difference_type __n,
334 const reverse_iterator<_Iterator>& __x)
335 { return reverse_iterator<_Iterator>(__x.base() - __n); }
336
337 // _GLIBCXX_RESOLVE_LIB_DEFECTS
338 // DR 280. Comparison of reverse_iterator to const reverse_iterator.
339 template<typename _IteratorL, typename _IteratorR>
340 inline bool
341 operator==(const reverse_iterator<_IteratorL>& __x,
342 const reverse_iterator<_IteratorR>& __y)
343 { return __x.base() == __y.base(); }
344
345 template<typename _IteratorL, typename _IteratorR>
346 inline bool
347 operator<(const reverse_iterator<_IteratorL>& __x,
348 const reverse_iterator<_IteratorR>& __y)
349 { return __y.base() < __x.base(); }
350
351 template<typename _IteratorL, typename _IteratorR>
352 inline bool
353 operator!=(const reverse_iterator<_IteratorL>& __x,
354 const reverse_iterator<_IteratorR>& __y)
355 { return !(__x == __y); }
356
357 template<typename _IteratorL, typename _IteratorR>
358 inline bool
359 operator>(const reverse_iterator<_IteratorL>& __x,
360 const reverse_iterator<_IteratorR>& __y)
361 { return __y < __x; }
362
363 template<typename _IteratorL, typename _IteratorR>
364 inline bool
365 operator<=(const reverse_iterator<_IteratorL>& __x,
366 const reverse_iterator<_IteratorR>& __y)
367 { return !(__y < __x); }
368
369 template<typename _IteratorL, typename _IteratorR>
370 inline bool
371 operator>=(const reverse_iterator<_IteratorL>& __x,
372 const reverse_iterator<_IteratorR>& __y)
373 { return !(__x < __y); }
374
375 template<typename _IteratorL, typename _IteratorR>
376#if __cplusplus201103L >= 201103L
377 // DR 685.
378 inline auto
379 operator-(const reverse_iterator<_IteratorL>& __x,
380 const reverse_iterator<_IteratorR>& __y)
381 -> decltype(__y.base() - __x.base())
382#else
383 inline typename reverse_iterator<_IteratorL>::difference_type
384 operator-(const reverse_iterator<_IteratorL>& __x,
385 const reverse_iterator<_IteratorR>& __y)
386#endif
387 { return __y.base() - __x.base(); }
388 //@}
389
390 // 24.4.2.2.1 back_insert_iterator
391 /**
392 * @brief Turns assignment into insertion.
393 *
394 * These are output iterators, constructed from a container-of-T.
395 * Assigning a T to the iterator appends it to the container using
396 * push_back.
397 *
398 * Tip: Using the back_inserter function to create these iterators can
399 * save typing.
400 */
401 template<typename _Container>
402 class back_insert_iterator
403 : public iterator<output_iterator_tag, void, void, void, void>
404 {
405 protected:
406 _Container* container;
407
408 public:
409 /// A nested typedef for the type of whatever container you used.
410 typedef _Container container_type;
411
412 /// The only way to create this %iterator is with a container.
413 explicit
414 back_insert_iterator(_Container& __x) : container(&__x) { }
415
416 /**
417 * @param __value An instance of whatever type
418 * container_type::const_reference is; presumably a
419 * reference-to-const T for container<T>.
420 * @return This %iterator, for chained operations.
421 *
422 * This kind of %iterator doesn't really have a @a position in the
423 * container (you can think of the position as being permanently at
424 * the end, if you like). Assigning a value to the %iterator will
425 * always append the value to the end of the container.
426 */
427#if __cplusplus201103L < 201103L
428 back_insert_iterator&
429 operator=(typename _Container::const_reference __value)
430 {
431 container->push_back(__value);
432 return *this;
433 }
434#else
435 back_insert_iterator&
436 operator=(const typename _Container::value_type& __value)
437 {
438 container->push_back(__value);
439 return *this;
440 }
441
442 back_insert_iterator&
443 operator=(typename _Container::value_type&& __value)
444 {
445 container->push_back(std::move(__value));
446 return *this;
447 }
448#endif
449
450 /// Simply returns *this.
451 back_insert_iterator&
452 operator*()
453 { return *this; }
454
455 /// Simply returns *this. (This %iterator does not @a move.)
456 back_insert_iterator&
457 operator++()
458 { return *this; }
459
460 /// Simply returns *this. (This %iterator does not @a move.)
461 back_insert_iterator
462 operator++(int)
463 { return *this; }
464 };
465
466 /**
467 * @param __x A container of arbitrary type.
468 * @return An instance of back_insert_iterator working on @p __x.
469 *
470 * This wrapper function helps in creating back_insert_iterator instances.
471 * Typing the name of the %iterator requires knowing the precise full
472 * type of the container, which can be tedious and impedes generic
473 * programming. Using this function lets you take advantage of automatic
474 * template parameter deduction, making the compiler match the correct
475 * types for you.
476 */
477 template<typename _Container>
478 inline back_insert_iterator<_Container>
479 back_inserter(_Container& __x)
480 { return back_insert_iterator<_Container>(__x); }
481
482 /**
483 * @brief Turns assignment into insertion.
484 *
485 * These are output iterators, constructed from a container-of-T.
486 * Assigning a T to the iterator prepends it to the container using
487 * push_front.
488 *
489 * Tip: Using the front_inserter function to create these iterators can
490 * save typing.
491 */
492 template<typename _Container>
493 class front_insert_iterator
494 : public iterator<output_iterator_tag, void, void, void, void>
495 {
496 protected:
497 _Container* container;
498
499 public:
500 /// A nested typedef for the type of whatever container you used.
501 typedef _Container container_type;
502
503 /// The only way to create this %iterator is with a container.
504 explicit front_insert_iterator(_Container& __x) : container(&__x) { }
505
506 /**
507 * @param __value An instance of whatever type
508 * container_type::const_reference is; presumably a
509 * reference-to-const T for container<T>.
510 * @return This %iterator, for chained operations.
511 *
512 * This kind of %iterator doesn't really have a @a position in the
513 * container (you can think of the position as being permanently at
514 * the front, if you like). Assigning a value to the %iterator will
515 * always prepend the value to the front of the container.
516 */
517#if __cplusplus201103L < 201103L
518 front_insert_iterator&
519 operator=(typename _Container::const_reference __value)
520 {
521 container->push_front(__value);
522 return *this;
523 }
524#else
525 front_insert_iterator&
526 operator=(const typename _Container::value_type& __value)
527 {
528 container->push_front(__value);
529 return *this;
530 }
531
532 front_insert_iterator&
533 operator=(typename _Container::value_type&& __value)
534 {
535 container->push_front(std::move(__value));
536 return *this;
537 }
538#endif
539
540 /// Simply returns *this.
541 front_insert_iterator&
542 operator*()
543 { return *this; }
544
545 /// Simply returns *this. (This %iterator does not @a move.)
546 front_insert_iterator&
547 operator++()
548 { return *this; }
549
550 /// Simply returns *this. (This %iterator does not @a move.)
551 front_insert_iterator
552 operator++(int)
553 { return *this; }
554 };
555
556 /**
557 * @param __x A container of arbitrary type.
558 * @return An instance of front_insert_iterator working on @p x.
559 *
560 * This wrapper function helps in creating front_insert_iterator instances.
561 * Typing the name of the %iterator requires knowing the precise full
562 * type of the container, which can be tedious and impedes generic
563 * programming. Using this function lets you take advantage of automatic
564 * template parameter deduction, making the compiler match the correct
565 * types for you.
566 */
567 template<typename _Container>
568 inline front_insert_iterator<_Container>
569 front_inserter(_Container& __x)
570 { return front_insert_iterator<_Container>(__x); }
571
572 /**
573 * @brief Turns assignment into insertion.
574 *
575 * These are output iterators, constructed from a container-of-T.
576 * Assigning a T to the iterator inserts it in the container at the
577 * %iterator's position, rather than overwriting the value at that
578 * position.
579 *
580 * (Sequences will actually insert a @e copy of the value before the
581 * %iterator's position.)
582 *
583 * Tip: Using the inserter function to create these iterators can
584 * save typing.
585 */
586 template<typename _Container>
587 class insert_iterator
588 : public iterator<output_iterator_tag, void, void, void, void>
589 {
590 protected:
591 _Container* container;
592 typename _Container::iterator iter;
593
594 public:
595 /// A nested typedef for the type of whatever container you used.
596 typedef _Container container_type;
597
598 /**
599 * The only way to create this %iterator is with a container and an
600 * initial position (a normal %iterator into the container).
601 */
602 insert_iterator(_Container& __x, typename _Container::iterator __i)
603 : container(&__x), iter(__i) {}
604
605 /**
606 * @param __value An instance of whatever type
607 * container_type::const_reference is; presumably a
608 * reference-to-const T for container<T>.
609 * @return This %iterator, for chained operations.
610 *
611 * This kind of %iterator maintains its own position in the
612 * container. Assigning a value to the %iterator will insert the
613 * value into the container at the place before the %iterator.
614 *
615 * The position is maintained such that subsequent assignments will
616 * insert values immediately after one another. For example,
617 * @code
618 * // vector v contains A and Z
619 *
620 * insert_iterator i (v, ++v.begin());
621 * i = 1;
622 * i = 2;
623 * i = 3;
624 *
625 * // vector v contains A, 1, 2, 3, and Z
626 * @endcode
627 */
628#if __cplusplus201103L < 201103L
629 insert_iterator&
630 operator=(typename _Container::const_reference __value)
631 {
632 iter = container->insert(iter, __value);
633 ++iter;
634 return *this;
635 }
636#else
637 insert_iterator&
638 operator=(const typename _Container::value_type& __value)
639 {
640 iter = container->insert(iter, __value);
641 ++iter;
642 return *this;
643 }
644
645 insert_iterator&
646 operator=(typename _Container::value_type&& __value)
647 {
648 iter = container->insert(iter, std::move(__value));
649 ++iter;
650 return *this;
651 }
652#endif
653
654 /// Simply returns *this.
655 insert_iterator&
656 operator*()
657 { return *this; }
658
659 /// Simply returns *this. (This %iterator does not @a move.)
660 insert_iterator&
661 operator++()
662 { return *this; }
663
664 /// Simply returns *this. (This %iterator does not @a move.)
665 insert_iterator&
666 operator++(int)
667 { return *this; }
668 };
669
670 /**
671 * @param __x A container of arbitrary type.
672 * @return An instance of insert_iterator working on @p __x.
673 *
674 * This wrapper function helps in creating insert_iterator instances.
675 * Typing the name of the %iterator requires knowing the precise full
676 * type of the container, which can be tedious and impedes generic
677 * programming. Using this function lets you take advantage of automatic
678 * template parameter deduction, making the compiler match the correct
679 * types for you.
680 */
681 template<typename _Container, typename _Iterator>
682 inline insert_iterator<_Container>
683 inserter(_Container& __x, _Iterator __i)
684 {
685 return insert_iterator<_Container>(__x,
686 typename _Container::iterator(__i));
687 }
688
689 // @} group iterators
690
691_GLIBCXX_END_NAMESPACE_VERSION
692} // namespace
693
694namespace __gnu_cxx _GLIBCXX_VISIBILITY(default)__attribute__ ((__visibility__ ("default")))
695{
696_GLIBCXX_BEGIN_NAMESPACE_VERSION
697
698 // This iterator adapter is @a normal in the sense that it does not
699 // change the semantics of any of the operators of its iterator
700 // parameter. Its primary purpose is to convert an iterator that is
701 // not a class, e.g. a pointer, into an iterator that is a class.
702 // The _Container parameter exists solely so that different containers
703 // using this template can instantiate different types, even if the
704 // _Iterator parameter is the same.
705 using std::iterator_traits;
706 using std::iterator;
707 template<typename _Iterator, typename _Container>
708 class __normal_iterator
709 {
710 protected:
711 _Iterator _M_current;
712
713 typedef iterator_traits<_Iterator> __traits_type;
714
715 public:
716 typedef _Iterator iterator_type;
717 typedef typename __traits_type::iterator_category iterator_category;
718 typedef typename __traits_type::value_type value_type;
719 typedef typename __traits_type::difference_type difference_type;
720 typedef typename __traits_type::reference reference;
721 typedef typename __traits_type::pointer pointer;
722
723 _GLIBCXX_CONSTEXPRconstexpr __normal_iterator() : _M_current(_Iterator()) { }
724
725 explicit
726 __normal_iterator(const _Iterator& __i) : _M_current(__i) { }
727
728 // Allow iterator to const_iterator conversion
729 template<typename _Iter>
730 __normal_iterator(const __normal_iterator<_Iter,
731 typename __enable_if<
732 (std::__are_same<_Iter, typename _Container::pointer>::__value),
733 _Container>::__type>& __i)
734 : _M_current(__i.base()) { }
735
736 // Forward iterator requirements
737 reference
738 operator*() const
739 { return *_M_current; }
740
741 pointer
742 operator->() const
743 { return _M_current; }
744
745 __normal_iterator&
746 operator++()
747 {
748 ++_M_current;
749 return *this;
750 }
751
752 __normal_iterator
753 operator++(int)
754 { return __normal_iterator(_M_current++); }
755
756 // Bidirectional iterator requirements
757 __normal_iterator&
758 operator--()
759 {
760 --_M_current;
761 return *this;
762 }
763
764 __normal_iterator
765 operator--(int)
766 { return __normal_iterator(_M_current--); }
767
768 // Random access iterator requirements
769 reference
770 operator[](const difference_type& __n) const
771 { return _M_current[__n]; }
772
773 __normal_iterator&
774 operator+=(const difference_type& __n)
775 { _M_current += __n; return *this; }
776
777 __normal_iterator
778 operator+(const difference_type& __n) const
779 { return __normal_iterator(_M_current + __n); }
780
781 __normal_iterator&
782 operator-=(const difference_type& __n)
783 { _M_current -= __n; return *this; }
784
785 __normal_iterator
786 operator-(const difference_type& __n) const
787 { return __normal_iterator(_M_current - __n); }
788
789 const _Iterator&
790 base() const
791 { return _M_current; }
792 };
793
794 // Note: In what follows, the left- and right-hand-side iterators are
795 // allowed to vary in types (conceptually in cv-qualification) so that
796 // comparison between cv-qualified and non-cv-qualified iterators be
797 // valid. However, the greedy and unfriendly operators in std::rel_ops
798 // will make overload resolution ambiguous (when in scope) if we don't
799 // provide overloads whose operands are of the same type. Can someone
800 // remind me what generic programming is about? -- Gaby
801
802 // Forward iterator requirements
803 template<typename _IteratorL, typename _IteratorR, typename _Container>
804 inline bool
805 operator==(const __normal_iterator<_IteratorL, _Container>& __lhs,
806 const __normal_iterator<_IteratorR, _Container>& __rhs)
807 { return __lhs.base() == __rhs.base(); }
808
809 template<typename _Iterator, typename _Container>
810 inline bool
811 operator==(const __normal_iterator<_Iterator, _Container>& __lhs,
812 const __normal_iterator<_Iterator, _Container>& __rhs)
813 { return __lhs.base() == __rhs.base(); }
814
815 template<typename _IteratorL, typename _IteratorR, typename _Container>
816 inline bool
817 operator!=(const __normal_iterator<_IteratorL, _Container>& __lhs,
818 const __normal_iterator<_IteratorR, _Container>& __rhs)
819 { return __lhs.base() != __rhs.base(); }
820
821 template<typename _Iterator, typename _Container>
822 inline bool
823 operator!=(const __normal_iterator<_Iterator, _Container>& __lhs,
824 const __normal_iterator<_Iterator, _Container>& __rhs)
825 { return __lhs.base() != __rhs.base(); }
19
←
Assuming the condition is true
→
20
←
Returning the value 1, which participates in a condition later
→
826
827 // Random access iterator requirements
828 template<typename _IteratorL, typename _IteratorR, typename _Container>
829 inline bool
830 operator<(const __normal_iterator<_IteratorL, _Container>& __lhs,
831 const __normal_iterator<_IteratorR, _Container>& __rhs)
832 { return __lhs.base() < __rhs.base(); }
833
834 template<typename _Iterator, typename _Container>
835 inline bool
836 operator<(const __normal_iterator<_Iterator, _Container>& __lhs,
837 const __normal_iterator<_Iterator, _Container>& __rhs)
838 { return __lhs.base() < __rhs.base(); }
839
840 template<typename _IteratorL, typename _IteratorR, typename _Container>
841 inline bool
842 operator>(const __normal_iterator<_IteratorL, _Container>& __lhs,
843 const __normal_iterator<_IteratorR, _Container>& __rhs)
844 { return __lhs.base() > __rhs.base(); }
845
846 template<typename _Iterator, typename _Container>
847 inline bool
848 operator>(const __normal_iterator<_Iterator, _Container>& __lhs,
849 const __normal_iterator<_Iterator, _Container>& __rhs)
850 { return __lhs.base() > __rhs.base(); }
851
852 template<typename _IteratorL, typename _IteratorR, typename _Container>
853 inline bool
854 operator<=(const __normal_iterator<_IteratorL, _Container>& __lhs,
855 const __normal_iterator<_IteratorR, _Container>& __rhs)
856 { return __lhs.base() <= __rhs.base(); }
857
858 template<typename _Iterator, typename _Container>
859 inline bool
860 operator<=(const __normal_iterator<_Iterator, _Container>& __lhs,
861 const __normal_iterator<_Iterator, _Container>& __rhs)
862 { return __lhs.base() <= __rhs.base(); }
863
864 template<typename _IteratorL, typename _IteratorR, typename _Container>
865 inline bool
866 operator>=(const __normal_iterator<_IteratorL, _Container>& __lhs,
867 const __normal_iterator<_IteratorR, _Container>& __rhs)
868 { return __lhs.base() >= __rhs.base(); }
869
870 template<typename _Iterator, typename _Container>
871 inline bool
872 operator>=(const __normal_iterator<_Iterator, _Container>& __lhs,
873 const __normal_iterator<_Iterator, _Container>& __rhs)
874 { return __lhs.base() >= __rhs.base(); }
875
876 // _GLIBCXX_RESOLVE_LIB_DEFECTS
877 // According to the resolution of DR179 not only the various comparison
878 // operators but also operator- must accept mixed iterator/const_iterator
879 // parameters.
880 template<typename _IteratorL, typename _IteratorR, typename _Container>
881#if __cplusplus201103L >= 201103L
882 // DR 685.
883 inline auto
884 operator-(const __normal_iterator<_IteratorL, _Container>& __lhs,
885 const __normal_iterator<_IteratorR, _Container>& __rhs)
886 -> decltype(__lhs.base() - __rhs.base())
887#else
888 inline typename __normal_iterator<_IteratorL, _Container>::difference_type
889 operator-(const __normal_iterator<_IteratorL, _Container>& __lhs,
890 const __normal_iterator<_IteratorR, _Container>& __rhs)
891#endif
892 { return __lhs.base() - __rhs.base(); }
893
894 template<typename _Iterator, typename _Container>
895 inline typename __normal_iterator<_Iterator, _Container>::difference_type
896 operator-(const __normal_iterator<_Iterator, _Container>& __lhs,
897 const __normal_iterator<_Iterator, _Container>& __rhs)
898 { return __lhs.base() - __rhs.base(); }
899
900 template<typename _Iterator, typename _Container>
901 inline __normal_iterator<_Iterator, _Container>
902 operator+(typename __normal_iterator<_Iterator, _Container>::difference_type
903 __n, const __normal_iterator<_Iterator, _Container>& __i)
904 { return __normal_iterator<_Iterator, _Container>(__i.base() + __n); }
905
906_GLIBCXX_END_NAMESPACE_VERSION
907} // namespace
908
909#if __cplusplus201103L >= 201103L
910
911namespace std _GLIBCXX_VISIBILITY(default)__attribute__ ((__visibility__ ("default")))
912{
913_GLIBCXX_BEGIN_NAMESPACE_VERSION
914
915 /**
916 * @addtogroup iterators
917 * @{
918 */
919
920 // 24.4.3 Move iterators
921 /**
922 * Class template move_iterator is an iterator adapter with the same
923 * behavior as the underlying iterator except that its dereference
924 * operator implicitly converts the value returned by the underlying
925 * iterator's dereference operator to an rvalue reference. Some
926 * generic algorithms can be called with move iterators to replace
927 * copying with moving.
928 */
929 template<typename _Iterator>
930 class move_iterator
931 {
932 protected:
933 _Iterator _M_current;
934
935 typedef iterator_traits<_Iterator> __traits_type;
936
937 public:
938 typedef _Iterator iterator_type;
939 typedef typename __traits_type::iterator_category iterator_category;
940 typedef typename __traits_type::value_type value_type;
941 typedef typename __traits_type::difference_type difference_type;
942 // NB: DR 680.
943 typedef _Iterator pointer;
944 typedef value_type&& reference;
945
946 move_iterator()
947 : _M_current() { }
948
949 explicit
950 move_iterator(iterator_type __i)
951 : _M_current(__i) { }
952
953 template<typename _Iter>
954 move_iterator(const move_iterator<_Iter>& __i)
955 : _M_current(__i.base()) { }
956
957 iterator_type
958 base() const
959 { return _M_current; }
960
961 reference
962 operator*() const
963 { return std::move(*_M_current); }
964
965 pointer
966 operator->() const
967 { return _M_current; }
968
969 move_iterator&
970 operator++()
971 {
972 ++_M_current;
973 return *this;
974 }
975
976 move_iterator
977 operator++(int)
978 {
979 move_iterator __tmp = *this;
980 ++_M_current;
981 return __tmp;
982 }
983
984 move_iterator&
985 operator--()
986 {
987 --_M_current;
988 return *this;
989 }
990
991 move_iterator
992 operator--(int)
993 {
994 move_iterator __tmp = *this;
995 --_M_current;
996 return __tmp;
997 }
998
999 move_iterator
1000 operator+(difference_type __n) const
1001 { return move_iterator(_M_current + __n); }
1002
1003 move_iterator&
1004 operator+=(difference_type __n)
1005 {
1006 _M_current += __n;
1007 return *this;
1008 }
1009
1010 move_iterator
1011 operator-(difference_type __n) const
1012 { return move_iterator(_M_current - __n); }
1013
1014 move_iterator&
1015 operator-=(difference_type __n)
1016 {
1017 _M_current -= __n;
1018 return *this;
1019 }
1020
1021 reference
1022 operator[](difference_type __n) const
1023 { return std::move(_M_current[__n]); }
1024 };
1025
1026 // Note: See __normal_iterator operators note from Gaby to understand
1027 // why there are always 2 versions for most of the move_iterator
1028 // operators.
1029 template<typename _IteratorL, typename _IteratorR>
1030 inline bool
1031 operator==(const move_iterator<_IteratorL>& __x,
1032 const move_iterator<_IteratorR>& __y)
1033 { return __x.base() == __y.base(); }
1034
1035 template<typename _Iterator>
1036 inline bool
1037 operator==(const move_iterator<_Iterator>& __x,
1038 const move_iterator<_Iterator>& __y)
1039 { return __x.base() == __y.base(); }
1040
1041 template<typename _IteratorL, typename _IteratorR>
1042 inline bool
1043 operator!=(const move_iterator<_IteratorL>& __x,
1044 const move_iterator<_IteratorR>& __y)
1045 { return !(__x == __y); }
1046
1047 template<typename _Iterator>
1048 inline bool
1049 operator!=(const move_iterator<_Iterator>& __x,
1050 const move_iterator<_Iterator>& __y)
1051 { return !(__x == __y); }
1052
1053 template<typename _IteratorL, typename _IteratorR>
1054 inline bool
1055 operator<(const move_iterator<_IteratorL>& __x,
1056 const move_iterator<_IteratorR>& __y)
1057 { return __x.base() < __y.base(); }
1058
1059 template<typename _Iterator>
1060 inline bool
1061 operator<(const move_iterator<_Iterator>& __x,
1062 const move_iterator<_Iterator>& __y)
1063 { return __x.base() < __y.base(); }
1064
1065 template<typename _IteratorL, typename _IteratorR>
1066 inline bool
1067 operator<=(const move_iterator<_IteratorL>& __x,
1068 const move_iterator<_IteratorR>& __y)
1069 { return !(__y < __x); }
1070
1071 template<typename _Iterator>
1072 inline bool
1073 operator<=(const move_iterator<_Iterator>& __x,
1074 const move_iterator<_Iterator>& __y)
1075 { return !(__y < __x); }
1076
1077 template<typename _IteratorL, typename _IteratorR>
1078 inline bool
1079 operator>(const move_iterator<_IteratorL>& __x,
1080 const move_iterator<_IteratorR>& __y)
1081 { return __y < __x; }
1082
1083 template<typename _Iterator>
1084 inline bool
1085 operator>(const move_iterator<_Iterator>& __x,
1086 const move_iterator<_Iterator>& __y)
1087 { return __y < __x; }
1088
1089 template<typename _IteratorL, typename _IteratorR>
1090 inline bool
1091 operator>=(const move_iterator<_IteratorL>& __x,
1092 const move_iterator<_IteratorR>& __y)
1093 { return !(__x < __y); }
1094
1095 template<typename _Iterator>
1096 inline bool
1097 operator>=(const move_iterator<_Iterator>& __x,
1098 const move_iterator<_Iterator>& __y)
1099 { return !(__x < __y); }
1100
1101 // DR 685.
1102 template<typename _IteratorL, typename _IteratorR>
1103 inline auto
1104 operator-(const move_iterator<_IteratorL>& __x,
1105 const move_iterator<_IteratorR>& __y)
1106 -> decltype(__x.base() - __y.base())
1107 { return __x.base() - __y.base(); }
1108
1109 template<typename _Iterator>
1110 inline auto
1111 operator-(const move_iterator<_Iterator>& __x,
1112 const move_iterator<_Iterator>& __y)
1113 -> decltype(__x.base() - __y.base())
1114 { return __x.base() - __y.base(); }
1115
1116 template<typename _Iterator>
1117 inline move_iterator<_Iterator>
1118 operator+(typename move_iterator<_Iterator>::difference_type __n,
1119 const move_iterator<_Iterator>& __x)
1120 { return __x + __n; }
1121
1122 template<typename _Iterator>
1123 inline move_iterator<_Iterator>
1124 make_move_iterator(_Iterator __i)
1125 { return move_iterator<_Iterator>(__i); }
1126
1127 template<typename _Iterator, typename _ReturnType
1128 = typename conditional<__move_if_noexcept_cond
1129 <typename iterator_traits<_Iterator>::value_type>::value,
1130 _Iterator, move_iterator<_Iterator>>::type>
1131 inline _ReturnType
1132 __make_move_if_noexcept_iterator(_Iterator __i)
1133 { return _ReturnType(__i); }
1134
1135 // @} group iterators
1136
1137_GLIBCXX_END_NAMESPACE_VERSION
1138} // namespace
1139
1140#define _GLIBCXX_MAKE_MOVE_ITERATOR(_Iter)std::make_move_iterator(_Iter) std::make_move_iterator(_Iter)
1141#define _GLIBCXX_MAKE_MOVE_IF_NOEXCEPT_ITERATOR(_Iter)std::__make_move_if_noexcept_iterator(_Iter) \
1142 std::__make_move_if_noexcept_iterator(_Iter)
1143#else
1144#define _GLIBCXX_MAKE_MOVE_ITERATOR(_Iter)std::make_move_iterator(_Iter) (_Iter)
1145#define _GLIBCXX_MAKE_MOVE_IF_NOEXCEPT_ITERATOR(_Iter)std::__make_move_if_noexcept_iterator(_Iter) (_Iter)
1146#endif // C++11
1147
1148#endif