StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
L2jetAlgo2012.cxx
1 #include <stdio.h>
2 #include <string.h>
3 #include <stdlib.h>
4 #include <time.h>
5 #include <math.h>
6 #include <fakeRtsLog.h>
7 
8 /***********************************************************
9  * $Id: L2jetAlgo2012.cxx,v 1.4 2012/03/21 18:18:03 jml Exp $
10  * \author Jan Balewski, IUCF, 2006
11  ***********************************************************
12  * Descripion:
13  Reco of mono- & di-jets in L2 using BTOW+ETOW
14  depends on L2-DB class
15  ***********************************************************
16  * version updated 2/17/2006 W.W. Jacobs
17  * 2008 to 2009 change by B.S. Page 11/12/2008
18  */
19 
20 
21 
22 #ifdef IS_REAL_L2 //in l2-ana environment
23  #include "../L2algoUtil/L2EmcDb2012.h"
24  #include "../L2algoUtil/L2Histo.h"
25 #else
26  #include "StTriggerUtilities/L2Emulator/L2algoUtil/L2EmcDb2012.h"
27  #include "StTriggerUtilities/L2Emulator/L2algoUtil/L2Histo.h"
28 #endif
29 
30 #include "L2jetAlgo2012.h"
31 #include "L2jetResults2012.h"
32 #include "Map_DeltaPhiJets.h"
33 
34 //=================================================
35 //=================================================
36 L2jetAlgo2012::L2jetAlgo2012(const char* name, const char *uid, L2EmcDb2012* db, char* outDir, int resOff, bool writeHighResult)
37  : L2VirtualAlgo2012( name, uid, db, outDir, true, true, resOff) {
38  /* called one per days
39  all memory allocation must be done here
40  */
41  namechar='L';
42  if (writeHighResult) namechar='H';
43 
44  createHisto();
45  run_number=-1;
46  LOG(DBG,"L2jetAlgo2012 instantiated, logPath='%s'\n",mOutDir1.c_str());
47  eve_Jet[0]=new L2Jet;
48  eve_Jet[1]=new L2Jet;
49 }
50 
51 /*========================================
52  ======================================== */
53 bool
54 L2jetAlgo2012::paramsChanged( int *rc_ints, float *rc_floats) {
55  int i;
56  for(i=0;i<5;i++)
57  if(rc_ints[i]!=raw_ints[i]) goto foundProblem;
58 
59  for(i=0;i<5;i++)
60  if(fabs(rc_floats[i]-raw_floats[i])>0.00001) goto foundProblem;
61 
62  return false;
63 
64  foundProblem:
65  if (mLogFile) fprintf(mLogFile,"L2jet-ghost initRun - inconsistent params, ABORT initialization\n");
66  return true;
67 }
68 
69 /*========================================
70  ======================================== */
71 int
72 L2jetAlgo2012::initRunUser( int runNo, int *rc_ints, float *rc_floats) {
73 
74  // clear content of all histograms
75  int i;
76  for (i=0; i<mxHA;i++) if(hA[i])hA[i]->reset();
77 
78  /* .... clear content, set threshold @ max as default */
79  memset(db_btowL2PhiBin, 0, sizeof(db_btowL2PhiBin));
80  memset(db_btowL2PatchBin,0, sizeof(db_btowL2PatchBin));
81 
82  memset(db_etowL2PhiBin, 0 ,sizeof(db_etowL2PhiBin));
83  memset(db_etowL2PatchBin,0 ,sizeof(db_etowL2PatchBin));
84 
85  run_startUnix=time(0);
86  run_number =runNo; // serves as a flag this run is initialized
87  raw_ints =rc_ints;
88  raw_floats =rc_floats;
89  run_nEventOneJet=run_nEventDiJet= run_nEventRnd=0;
90 
91  mEventsInRun=0;
92 
93  if( mLogFile==0) LOG(ERR," L2jetAlgo2012() UNABLE to open run summary log file, continue anyhow\n");
94 
95  // unpack params from run control GUI
96  par_useBtowEast= (rc_ints[0]&1 )>0;
97  par_useBtowWest= (rc_ints[0]&2)>0;
98  par_useEndcap = rc_ints[1]&1;
99  int par_adcThr = rc_ints[2]; // needed only in initRun()
100  par_dbg = rc_ints[3];
101  par_RndAcceptPrescale = rc_ints[4];
102  par_minPhiBinDiff= 7; //for 1X1 JP was 5 for .6x.6 JP BP
103 
104  par_oneJetThr = rc_floats[0];
105  par_diJetThrHigh= rc_floats[1];
106  par_diJetThrLow = rc_floats[2];
107  par_diJetThr_2 = rc_floats[3]; // wwj 2/11
108  par_diJetThr_3 = rc_floats[4]; // wwj 2/11
109  par_diJetThr_4 = rc_floats[5]; // wwj 2/11
110  par_diJetThr_5 = rc_floats[6]; // wwj 2/11
111  par_diJetEtaHigh = rc_floats[7]; // wwj 2/11
112  par_diJetEtaLow = rc_floats[8]; // wwj 2/11
113 
114  // set other hardcoded or calculated params
115 
116  // to monitor hot towers in E+B Emc
117  float monTwEtThr=2.0; // (GeV) Et threshold,WARN, edit histo title by hand
118  par_hotTwEtThr = monTwEtThr;
119 
120  if(par_dbg) LOG(DBG,"Brian version is running!\n");
121 
122  if (mLogFile) {
123  fprintf(mLogFile,"L2jet algorithm initRun(%d), compiled: %s , %s\n params:\n",run_number,__DATE__,__TIME__);
124  fprintf(mLogFile," - use BTOW: East=%d West=%d, Endcap=%d L2ResOffset=%d\n", par_useBtowEast, par_useBtowWest,par_useEndcap ,mResultOffset);
125  fprintf(mLogFile," - threshold: ADC-ped> %d \n", par_adcThr);
126  fprintf(mLogFile," - min phi opening angle Jet1<->Jet2: %d in L2phiBins\n",par_minPhiBinDiff);
127  fprintf(mLogFile," - diJet Et thrHigh= %.2f (GeV) thrLow= %.2f (GeV)\n", par_diJetThrHigh, par_diJetThrLow);
128  fprintf(mLogFile," - new diJet Et thr_2= %.2f (GeV) ; thr_3= %.2f (Gev) ; thr_4= %.2f (GeV) ; thr_5= %.2f (GeV)\n",par_diJetThr_2, par_diJetThr_3, par_diJetThr_4, par_diJetThr_5); // wwj 2/10
129  fprintf(mLogFile," - oneJet Et thr = %.2f (GeV) ; rndAccPrescale=%d \n",par_oneJetThr,par_RndAcceptPrescale);
130  fprintf(mLogFile," - new diJet Eta_low= %.2f ; diJet Eta_High= %.2f \n", par_diJetEtaLow, par_diJetEtaHigh); // wwj 2/10
131  fprintf(mLogFile," - debug=%d, hot tower threshold: Et> %.1f GeV ( only monitoring)\n",par_dbg, monTwEtThr);
132  }
133 
134  // verify consistency of input params
135  int kBad=0;
136  kBad+=0x0001 * ( !par_useBtowEast & !par_useBtowWest & !par_useEndcap);
137  // par_pedOff no longer exists: kBad+=0x0002 * (par_adcThr<par_pedOff);
138  kBad+=0x0004 * (par_adcThr>16);
139  kBad+=0x0008 * (par_minPhiBinDiff<5);
140  kBad+=0x0010 * (par_minPhiBinDiff>=15);
141  kBad+=0x0020 * (par_oneJetThr<3.);
142  kBad+=0x0040 * (par_oneJetThr>12.);
143  kBad+=0x0080 * (par_diJetThrLow<1.9); // LOOK: wwj 2/11 using lowest threshold
144  kBad+=0x0100 * (par_diJetThrHigh<par_diJetThrLow);
145  kBad+=0x0200 * (par_diJetThrHigh>12.);
146  //0x0400 not used
147  kBad+=0x0800 * ( par_RndAcceptPrescale<0 );
148  kBad+=0x1000 * (par_diJetThr_2<par_diJetThrLow || par_diJetThr_5>par_diJetThrHigh); // wwj 2/16
149  kBad+=0x2000 * (par_diJetThr_2>par_diJetThr_3 || par_diJetThr_3>par_diJetThr_4 || par_diJetThr_4>par_diJetThr_5); // wwj 2/16
150  kBad+=0x4000 * (par_diJetEtaLow>par_diJetEtaHigh); // wwj 2/16
151  if (mLogFile) {
152  fprintf(mLogFile,"L2jet initRun() params checked for consistency, Error flag=0x%04x\n",kBad);
153  if(kBad) fprintf(mLogFile,"L2jet initRun() ABORT\n");
154  }
155 
156  if(kBad) {
157  run_number=-66;
158  if (mLogFile) {
159  fprintf(mLogFile,"L2jet algorithm init: crashB in internal logic\n");
160  fclose(mLogFile);
161  return kBad;
162  }
163  }
164 
165  char tit[100];
166  sprintf(tit,"# BTOW towers>ped+%d (input); x: # of towers/event",par_adcThr);
167  hA[47]->setTitle(tit);
168 
169  sprintf(tit,"# ETOW towers>ped+%d (input); x: # of towers/event",par_adcThr);
170  hA[48]->setTitle(tit);
171 
172 
173  const int mxEtaBinsE=12,mxEtaBinsB=40;
174 
175  // rebuild local lookup tables
176 
177  int etowEtaBin2Patch[mxEtaBinsE]={14,14,13,13,12,12,11,11,11,10,10,10};
178 
179  int nB=0, nE=0; /* counts # of unmasekd towers */
180  int nBg=0, nEg=0; /* counts # of reasonable calibrated towers */
181 
182  for(i=0; i<EmcDbIndexMax; i++) {
183  const L2EmcDb2012::EmcCDbItem *x=mDb->getByIndex(i);
184  if(mDb->isEmpty(x)) continue; /* dropped not mapped channels */
185  if(x->fail) continue; /* dropped masked channels */
186  if(x->gain<=0) continue; /* dropped uncalibrated towers , tmp */
187  /* if(x->sec!=1) continue; tmp, to test patch mapping */
188 
189  /* WARN, calculate index to match RDO order in the ADC data banks */
190  int ietaP, iphiP;
191  if (mDb->isBTOW(x) ) {
192  /*....... B A R R E L .................*/
193  nB++;
194  if(x->eta<0 || x->eta>mxEtaBinsB) goto crashIt_1;
195  if(!par_useBtowEast && x->eta<=20) continue;
196  if(!par_useBtowWest && x->eta>=21) continue;
197  ietaP= (x->eta-1)/4; /* correct */
198  int iphiTw=(x->sec-1)*10 + x->sub-'a';
199  // allign in phi TP @ L2 w/ L0
200  iphiTw--;
201  if(iphiTw<0) iphiTw=119;
202  // now cr0x1e, mod1, subm2, is beginning of the first BTOW TP
203  iphiP= iphiTw/4 ; /* correct */
204  // if(ietaP==0 && iphiP==5)
205  db_btowL2PhiBin[x->rdo]=iphiP;
206  db_btowL2PatchBin[x->rdo]=ietaP+ iphiP*cl2jetMaxEtaBins;
207  nBg++;
208  } else if(mDb->isETOW(x) && par_useEndcap) {
209  /*....... E N D C A P ........................*/
210  nE++;
211  int iphiTw= (x->sec-1)*5 + x->sub-'A';
212  // allign in phi TP @ L2 w/ L0
213  iphiTw--;
214  if(iphiTw<0) iphiTw=59;
215  // now subsector 01TB is beginning of the first(==0) ETOW TP in phi
216  iphiP= iphiTw/2 ; /* correct */
217  if(x->eta<0 || x->eta>mxEtaBinsE) goto crashIt_1;
218  ietaP=etowEtaBin2Patch[x->eta-1];
219  db_etowL2PhiBin[x->rdo]=iphiP;
220  db_etowL2PatchBin[x->rdo]=ietaP+ iphiP*cl2jetMaxEtaBins;
221  nEg++;
222  }
223 
224  }
225 
226  if (mLogFile) {
227  fprintf(mLogFile,"L2jet algorithm: found working/calibrated: %d/%d=ETOW & %d/%d=BTOW, based on ASCII DB\n",nE,nEg,nB,nBg);
228  }
229 
230  return 0; //OK
231 
232  crashIt_1: /* fatal initialization error */
233  run_number=-55;
234  if (mLogFile) {
235  fprintf(mLogFile,"L2jet algorithm init: crashC in internal logic\n");
236  fclose(mLogFile);
237  }
238  return -6;
239 
240 }
241 
242 /*========================================
243  ======================================== */
244 void L2jetAlgo2012::computeUser(int token){
245  //put in to satisfy l2virtualalgo2012
246 }
247 
248 /*========================================
249  ======================================== */
250 bool
251 L2jetAlgo2012::decisionUser(int token, int *myL2Result){
252  /* STRICT TIME BUDGET START ....*/
253 
254  //Ross says: Eventually you should take these out - timing is done in the virtual algo class as well.
255  unsigned long mEveTimeStart,mEveTimeStop, mEveTimeDiff;
256 rdtscl_macro(mEveTimeStart);
257 
258  //Ross says: I will check how the following line is done in L2bemcGamma2012
259  //following if statement was eve_ID!=inpEveId BP
260  if(1) {//UUUU
261  //.................... this event has NOT been processed
262 
263  /*
264  Chris doesn't want us to write anything out
265  during event processing ...
266  */
267 
268  bool bemcIn = true;
269  bool eemcIn = true;
270 
271  //eve_ID=inpEveId; // every events is processed only once
272  mEventsInRun++;
273  clearEvent(); /* price=13 kTicks */
274  int runTimeSec=time(0)- run_startUnix;
275  hA[10]->fill(0);
276  hA[12]->fill(runTimeSec);
277 
278  if(par_dbg>1) LOG(DBG,"\n......... in L2Jet_doEvent(ID=%d)... bIn=%d eIn=%d\n",eve_ID,bemcIn,eemcIn);//bemcIn and eemcIn were passed to the old doEvent
279 
280 
281  //Ross asks: will the jet code run on just barrel or just endcap?
282  // if that's the case, will it save time to move those sections of the
283  // code to computeUser() rather than decisionUser()?
284  if (bemcIn || eemcIn){//VVVV this algo has nothing to do w/o any Ecal data
285 
286  //===== step 1: unpack raw data blocks ==============
287 
288  int nBtowTw=0, nEtowTw=0;
289  //Ross says: I've changed projectAdc to use the L2eventStream2012 class
290  /*......... BBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBB */
291 
292  //we can move this to compute.
293  if(bemcIn==1 && (par_useBtowEast||par_useBtowWest) ) {
294 
295  nBtowTw=projectAdc( mEveStream_btow[token].get_hits(), mEveStream_btow[token].get_hitSize(),
296  db_btowL2PhiBin, db_btowL2PatchBin,
297  hA[20] );
298  }
299 
300  /*........... EEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEE */
301  //and this as well.
302  if(eemcIn==1 && par_useEndcap ) {
303  nEtowTw=projectAdc( mEveStream_etow[token].get_hits(), mEveStream_etow[token].get_hitSize(),
304  db_etowL2PhiBin, db_etowL2PatchBin,
305  hA[30] );
306  }
307 
308  //============Do 2D scan==========
309  float totEneGeV=true2Dscan();
310  int itotEneGeV=(int)totEneGeV;
311  int iK;
312  for(iK=0; iK<mxJ; iK++)
313  {
314  weightedEtaPhi(iK);
315  L2Jet *K=eve_Jet[iK];
316  K->eneGeV=K->iene;
317  K->phiRad=0.21*(6.0-K->fphiBin);
318  while(K->phiRad<0) K->phiRad+=6.2832;
319  while(K->phiRad>6.2832) K->phiRad-=6.2832;
320  }
321 
322  if(eve_Jet[0]->eneGeV <eve_Jet[1]->eneGeV) {// swap jets fo E1>E2
323  L2Jet *Jx=eve_Jet[0];
324  eve_Jet[0]=eve_Jet[1];
325  eve_Jet[1]=Jx;
326  }
327 
328  if(par_dbg>2) LOG(DBG,"doEvent iphiBin1=%d iene1=%f , iphiBin2=%d iene2=%f rms1PhiBin=%f \n",eve_Jet[0]->iphiBin,eve_Jet[0]->iene,eve_Jet[1]->iphiBin,eve_Jet[1]->iene,eve_Jet[0]->rmsPhiBin);
329 
330  //====== step 4: make trigger decisions====
331 
332  // now put in eta vs. Jet Threshold "contour" cuts ... this vesion cuts on eta_1 + eta_2 .. wwj 2/16
333  // looks like fetaBin goes 0.0->15.0 for eta -1. -> +2. ... so "eta" = fetaBin*0.2 -1.0 ... wwj 2/10
334 
335  float rjetEta_0=eve_Jet[0]->fetaBin*0.2 -1.0; // wwj 2/10
336  float rjetEta_1=eve_Jet[1]->fetaBin*0.2 -1.0; // wwj 2/10
337 
338  bool acceptDiJet_EE=((rjetEta_0+rjetEta_1) > par_diJetEtaHigh && eve_Jet[0]->eneGeV > par_diJetThr_4 && eve_Jet[1]->eneGeV > par_diJetThrLow); // wwj 2/10
339 
340  bool acceptDiJet_EB=((rjetEta_0+rjetEta_1) > par_diJetEtaLow && eve_Jet[0]->eneGeV > par_diJetThr_5 && eve_Jet[1]->eneGeV > par_diJetThr_2); // wwj 2/10
341 
342 
343  // bool acceptDiJet=( eve_Jet[0]->eneGeV > par_diJetThrHigh) && ( eve_Jet[1]->eneGeV > par_diJetThrLow);
344  bool acceptDiJet_BB=(eve_Jet[0]->eneGeV > par_diJetThrHigh && eve_Jet[1]->eneGeV > par_diJetThr_3); // wwj 2/10
345 
346  bool acceptDiJet=acceptDiJet_EE || acceptDiJet_EB || acceptDiJet_BB; // wwj 2/10
347 
348  bool acceptOneJet=( eve_Jet[0]->eneGeV > par_oneJetThr) ;
349 
350  bool acceptRnd=mRandomAccept>0; // provided by Virtual algo
351  mAccept=acceptDiJet || acceptOneJet || acceptRnd;
352 
353  //====== step 5: update various monitorig histos
354 
355  // histogramming reco Et1-Et2, no cuts
356  int iet1 =(int)eve_Jet[0]->eneGeV;
357  int iet2 =(int)eve_Jet[1]->eneGeV;
358  int ieta1=(int)eve_Jet[0]->fetaBin;
359  int ieta2=(int)eve_Jet[1]->fetaBin;
360  int iphi1=(int)eve_Jet[0]->fphiBin;
361  int iphi2=(int)eve_Jet[1]->fphiBin;
362 
363  hA[40]->fill(iet1,iet2);
364 
365  hA[41]->fill(ieta1,iphi1);
366  hA[42]->fill(ieta2,iphi2);
367  hA[43]->fill(iphi1,iphi2);
368  hA[44]->fill(iet1);
369  hA[45]->fill(iet2);
370  hA[46]->fill(itotEneGeV);
371  hA[47]->fill(nBtowTw);
372  hA[48]->fill(nEtowTw);
373 
374  // sivers delta zeta, the map is still worng
375  int kphi1=int(eve_Jet[0]->phiRad*10.);
376  int kphi2=int(eve_Jet[1]->phiRad*10.);
377  int idelZeta=map_DelPhiJets[kphi1*MxPhiRad10 + kphi2];
378 
379  if( mAccept) hA[10]->fill(8);
380 
381  if(acceptOneJet ){
382  hA[10]->fill(4);
383  run_nEventOneJet++;
384  hA[13]->fill(runTimeSec);
385  hA[50]->fill(iet1);
386  hA[51]->fill(ieta1,iphi1);
387  hA[52]->fill(ieta1);
388  hA[53]->fill(iphi1);
389  }
390 
391  if(acceptDiJet ){
392  hA[10]->fill(5);
393  run_nEventDiJet++;
394  hA[14]->fill(runTimeSec);
395  hA[60]->fill(iet1,iet2);
396  hA[61]->fill(ieta1,iphi1);
397  hA[62]->fill(ieta2,iphi2);
398  hA[63]->fill(iphi1,iphi2);
399  hA[64]->fill(iet1);
400  hA[65]->fill(iet2);
401  hA[66]->fill(ieta1);
402  hA[67]->fill(ieta2);
403  hA[68]->fill(iphi1);
404  hA[69]->fill(iphi2);
405  hA[70]->fill(idelZeta);
406  hA[71]->fill(ieta1,idelZeta);
407  hA[72]->fill(ieta1,ieta2);
408  hA[73]->fill((iphi1+iphi2)/2,idelZeta);
409  hA[74]->fill(itotEneGeV);
410  }
411 
412  if(acceptDiJet_EE ){
413  run_nEventDiJet++; // all these histos #80-#94 added ... wwj 2/10
414  hA[80]->fill(iet1,iet2);
415  hA[81]->fill(ieta1,iphi1);
416  hA[82]->fill(ieta2,iphi2);
417  hA[83]->fill(iphi1,iphi2);
418  hA[84]->fill(iet1);
419  hA[85]->fill(iet2);
420  hA[86]->fill(ieta1);
421  hA[87]->fill(ieta2);
422  hA[88]->fill(iphi1);
423  hA[89]->fill(iphi2);
424  hA[90]->fill(idelZeta);
425  hA[91]->fill(ieta1,idelZeta);
426  hA[92]->fill(ieta1,ieta2);
427  hA[93]->fill((iphi1+iphi2)/2,idelZeta);
428  hA[94]->fill(itotEneGeV);
429  }
430 
431  if(acceptDiJet_EB ){
432  run_nEventDiJet++; // all these histos #100-#114 added ... wwj 2/10
433  hA[100]->fill(iet1,iet2);
434  hA[101]->fill(ieta1,iphi1);
435  hA[102]->fill(ieta2,iphi2);
436  hA[103]->fill(iphi1,iphi2);
437  hA[104]->fill(iet1);
438  hA[105]->fill(iet2);
439  hA[106]->fill(ieta1);
440  hA[107]->fill(ieta2);
441  hA[108]->fill(iphi1);
442  hA[109]->fill(iphi2);
443  hA[110]->fill(idelZeta);
444  hA[111]->fill(ieta1,idelZeta);
445  hA[112]->fill(ieta1,ieta2);
446  hA[113]->fill((iphi1+iphi2)/2,idelZeta);
447  hA[114]->fill(itotEneGeV);
448  }
449 
450  if(acceptDiJet_BB ){
451  run_nEventDiJet++; // all these histos #120-#134 added ... wwj 2/10
452  hA[120]->fill(iet1,iet2);
453  hA[121]->fill(ieta1,iphi1);
454  hA[122]->fill(ieta2,iphi2);
455  hA[123]->fill(iphi1,iphi2);
456  hA[124]->fill(iet1);
457  hA[125]->fill(iet2);
458  hA[126]->fill(ieta1);
459  hA[127]->fill(ieta2);
460  hA[128]->fill(iphi1);
461  hA[129]->fill(iphi2);
462  hA[130]->fill(idelZeta);
463  hA[131]->fill(ieta1,idelZeta);
464  hA[132]->fill(ieta1,ieta2);
465  hA[133]->fill((iphi1+iphi2)/2,idelZeta);
466  hA[134]->fill(itotEneGeV);
467  }
468 
469  if(mRandomAccept){
470  hA[10]->fill(6);
471  run_nEventRnd++;
472  hA[15]->fill(runTimeSec);
473  }
474 
475 
476  //====== step 6: fill L2Result (except time)
477  L2jetResults2012 out; // all output bits lump together
478  memset(&out,0,sizeof(out)); // clear content
479 
480 
481  // add acceptDiJet_EE, acceptDiJet_EB and acceptDiJet_BB to below ... wwj 2/10
482  out.int0.version=L2JET_RESULTS_VERSION;
483  out.int0.decision=
484  ( par_useBtowEast <<0 ) +
485  ( par_useBtowWest <<1 ) +
486  ( par_useEndcap <<2 ) +
487  ( (bemcIn>0) <<3 ) +
488  ( (eemcIn>0) <<4 ) +
489  ( mRandomAccept <<5 ) +
490  ( acceptOneJet <<6 ) +
491  ( acceptDiJet <<7 ) +
492  ( acceptDiJet_EE <<8 ) +
493  ( acceptDiJet_EB <<9 ) +
494  ( acceptDiJet_BB <<10 ) ; // e.g., added last three lines here ... wwj 2/10
495  out.int0.dumm=namechar;
496 
497  out.int1.iTotEne=(unsigned short)(totEneGeV*100.); // now 1=10 MeV
498  out.int2.nBtowTw=nBtowTw;
499  out.int2.nEtowTw=nEtowTw;
500 
501  out.jet1.jPhi=(int)(eve_Jet[0]->phiRad*28.65); //so phi/deg=2*jPhi
502  out.jet1.jEta=(int)(eve_Jet[0]->fetaBin*10.);
503  out.jet1.iEne=(unsigned short)(eve_Jet[0]->eneGeV*100.); // now 1=10 MeV
504 
505  out.jet2.jPhi=(int)(eve_Jet[1]->phiRad*28.65);
506  out.jet2.jEta=(int)(eve_Jet[1]->fetaBin*10.);
507  out.jet2.iEne=(unsigned short)(eve_Jet[1]->eneGeV*100.); // now 1=10 MeV
508 
509  out.jet1.rmsEta=(unsigned short)(eve_Jet[0]->rmsEtaBin*20.); //Fix: 1=0.01 physical eta
510  out.jet1.rmsPhi=(unsigned short)(eve_Jet[0]->rmsPhiBin*120.); //1=0.1 deg
511 
512  out.jet2.rmsEta=(unsigned short)(eve_Jet[1]->rmsEtaBin*20.); //Fix: 1=0.01 physical eta
513  out.jet2.rmsPhi=(unsigned short)(eve_Jet[1]->rmsPhiBin*120.);// 1=0.1 deg
514 
515 
516  rdtscl_macro(mEveTimeStop);
517  mEveTimeDiff=mEveTimeStop-mEveTimeStart;
518  int kTick=mEveTimeDiff/1000;
519 
520  //mhT->fill(kTick); taken out 12/11/08 BP
521 
522  out.int0.kTick= kTick>255 ? 255 : kTick;
523 
524  //calculate and match the check sum
525  out.int1.checkSum=-L2jetResults2012_doCheckSum(&out);
526 
527 
528 
529  //===== step 7: write L2Result
530  memcpy(myL2Result,&out,sizeof(L2jetResults2012));
531 
532  // dirty tests, clean it up before real use
533 
534  if(par_dbg){//WWWW
535  L2jetResults2012_print(&out);
536  LOG(DBG," phiRad1=%f phiRad2=%f \n",eve_Jet[0]->phiRad,eve_Jet[1]->phiRad);
537  LOG(DBG,"idelZeta=%d delZeta/deg=%.1f \n\n",idelZeta,idelZeta/31.416*180);
538 
539  //tmp printouts of errors:
540  if( out.jet1.iEne+out.jet2.iEne > out.int1.iTotEne) {
541  }
542  if(iphi1==iphi2) {
543  LOG(ERR,"L2jet-fatal error,neveId=%d, phi1,2=%d,%d\n",mEventsInRun,iphi1,iphi2);
544  dumpPatchEneA();
545  }
546 
547  if( L2jetResults2012_doCheckSum(&out)) {
548  LOG(ERR,"L2jet-fatal error, wrong cSum=%d\n", L2jetResults2012_doCheckSum(&out));
549  L2jetResults2012_print(&out);
550  }
551 
552  } // end of WWWW
553 
554  }// end of VVVV (etow or btow has some data)
555  }// end of UUUU event processing
556 
557 
558  return mAccept;
559 }
560 
561 
562 /*========================================
563  ======================================== */
564 void
565 L2jetAlgo2012::finishRunUser() { /* called once at the end of the run */
566  if(run_number<0) return; // already finished
567 
568  if (mLogFile) {
569  fprintf(mLogFile,"L2-jet algorithm finishRun(%d)\n",run_number);
570  fprintf(mLogFile," - %d events seen by L2 di-jet\n",mEventsInRun);
571  fprintf(mLogFile," - accepted: rnd=%d oneJet=%d diJet=%d \n", run_nEventRnd, run_nEventOneJet, run_nEventDiJet);
572 
573  // print few basic histos
574  hA[10]->printCSV(mLogFile); // event accumulated
575 
576  }
577  finishRunHisto(); // still needs current DB
578  if( mHistFile==0) {
579  LOG(ERR," L2jetAlgo2012: finishRun() UNABLE to open run summary log file, continue anyhow\n");
580  if (mLogFile)
581  fprintf(mLogFile,"L2 di-jet histos NOT saved, I/O error\n");
582  } else { // save histos
583  int j;
584  int nh=0;
585  for(j=0;j<mxHA;j++) {
586  if(hA[j]==0) continue;
587  hA[j]->write(mHistFile);
588  nh++;
589  }
590  }
591 
592  run_number=-2; // clear run #
593 
594 }
595 
596 //=======================================
597 //=======================================
598 void
599 L2jetAlgo2012::createHisto() {
600  memset(hA,0,sizeof(hA));
601 
602  hA[10]=new L2Histo(10,"total event counter; x=cases",9);
603  hA[11]=new L2Histo(11,"L2 time used per input event; x: time (CPU kTics), range=100muSec; y: events ",160);
604 
605  int mxRunDration=2500;
606  hA[12]=new L2Histo(12,"rate of input events; x: time in this run (seconds); y: rate (Hz)", mxRunDration);
607 
608  hA[13]=new L2Histo(13,"rate of accepted one-Jet; x: time in this run (seconds); y: rate (Hz)", mxRunDration);
609  hA[14]=new L2Histo(14,"rate of accepted di-Jet ; x: time in this run (seconds); y: rate (Hz)", mxRunDration);
610  hA[15]=new L2Histo(15,"rate of random accepted ; x: time in this run (seconds); y: rate (Hz)", mxRunDration);
611 
612  // BTOW raw spectra
613  hA[20]=new L2Histo(20,"BTOW tower, Et>2.0 GeV (input); x: BTOW RDO index=chan*30+fiber; y: counts", 4800);
614  hA[21]=new L2Histo(21,"BTOW tower, Et>2.0 GeV (input); x: BTOW softID", 4800);
615  hA[22]=new L2Histo(22,"BTOW tower, Et>2.0 GeV (input); x: eta bin, [-1,+1]; y: phi bin ~sector",40,120);
616 
617  // ETOW raw spectra
618  hA[30]=new L2Histo(30,"ETOW tower, Et>2.0 GeV (input); x: ETOW RDO index=chan*6+fiber; y: counts", 720 );
619  hA[31]=new L2Histo(31,"ETOW tower, Et>2.0 GeV (input); x: i=chan+128*crate", 768);
620  hA[32]=new L2Histo(32,"ETOW tower, Et>2.0 GeV (input); x: 12 - Endcap etaBin ,[+1,+2]; y: phi bin ~sector",12,60);
621 
622  // Di-Jet raw yields
623  hA[40]=new L2Histo(40,"Et Jet1-Jet2 (input); x: Jet1 Et/GeV ; Jet2 Et/GeV",30,30);
624  hA[41]=new L2Histo(41,"diJet1 eta-phi (input); x: iEta [-1,+2] ; y: iPhi ~sector ",15,30);
625  hA[42]=new L2Histo(42,"diJet2 eta-phi (input); x: iEta [-1,+2] ; y: iPhi ~sector",15,30);
626 
627  hA[43]=new L2Histo(43,"diJet phi1-phi2 (input); x: iPhi1 ~sector ; y: iPhi2 ~sector ",30,30);
628 
629  hA[44]=new L2Histo(44,"Jet1 Et (input); x: Et (GeV)", 60);
630  hA[45]=new L2Histo(45,"Jet2 Et (input); x: Et (GeV)", 60);
631  hA[46]=new L2Histo(46,"total Et (input); x: Et (GeV)", 100);
632  hA[47]=new L2Histo(47,"# BTOW towers>thrXX (input); x: # of towers/event", 200);
633  hA[48]=new L2Histo(48,"# ETOW towers>thrXX (input); x: # of towers/event", 100);
634 
635  // ........accepted one-jet events
636  hA[50]=new L2Histo(50,"one-Jet Et (accepted); x: jet Et (GeV)", 60);
637  hA[51]=new L2Histo(51,"one-Jet eta-phi (accepted); x: iEta [-1,+2] ; y: iPhi ~sector ",15,30);
638  hA[52]=new L2Histo(52,"one-Jet eta (accepted); x: iEta [-1,+2]", 15);
639  hA[53]=new L2Histo(53,"one-Jet phi (accepted); x: iPhi ~sector", 30);
640 
641  // Di-Jet accepted
642  hA[60]=new L2Histo(60,"Et of Jet1 vs. Jet2 (accepted); x: Jet1/GeV ; Jet2/GeV",30,30);
643  hA[61]=new L2Histo(61,"diJet1 eta-phi (accepted); x: iEta [-1,+2] ; y: iPhi ~sector ",15,30);
644  hA[62]=new L2Histo(62,"diJet2 eta-phi (accepted); x: iEta [-1,+2] ; y: iPhi ~sector",15,30);
645 
646  hA[63]=new L2Histo(63,"diJet phi1-phi2 (accepted); x: iPhi1 ~sector ; y: iPhi2 ~sector ",30,30);
647 
648  hA[64]=new L2Histo(64,"diJet1 Et (accepted); x: Et (GeV)", 60);
649  hA[65]=new L2Histo(65,"diJet2 Et (accepted); x: Et (GeV)", 60);
650 
651  hA[66]=new L2Histo(66,"diJet1 eta (accepted); x: i Eta [-1,+2]", 15);
652  hA[67]=new L2Histo(67,"diJet2 eta (accepted); x: i Eta [-1,+2]", 15);
653  hA[68]=new L2Histo(68,"diJet1 phi (accepted); x: iPhi ~sector", 30);
654  hA[69]=new L2Histo(69,"diJet2 phi (accepted); x: iPhi ~sector", 30);
655  hA[70]=new L2Histo(70,"diJet delZeta (accepted); x: delta zeta (rad*10)", MxPhiRad10);
656  hA[71]=new L2Histo(71,"diJet delZeta vs. eta1 (accepted); x: iEta1 [-1,+2] ; y: delta zeta (rad*10)",15, MxPhiRad10);
657  hA[72]=new L2Histo(72,"diJet eta2 vs. eta1 (accepted); x: iEta1 [-1,+2] ;x: iEta2 [-1,+2] ",15,15);
658  hA[73]=new L2Histo(73,"diJet delZeta vs. avrPhi (accepted); x: (iphi1+iphi2)/2 (12 deg/bin); y: delta zeta (rad*10)",30, MxPhiRad10);
659  hA[74]=new L2Histo(74,"total Et diJet (accepted); x: Et (GeV)", 100);
660 
661  // Di-Jet_EE accepted ... all these histos #80-#94 added ... wwj 2/10
662  hA[80]=new L2Histo(80,"Et of Jet1 vs. Jet2 (accepted_EE); x: Jet1/GeV ; Jet2/GeV",30,30);
663  hA[81]=new L2Histo(81,"diJet1 eta-phi (accepted_EE); x: iEta [-1,+2] ; y: iPhi ~sector ",15,30);
664  hA[82]=new L2Histo(82,"diJet2 eta-phi (accepted_EE); x: iEta [-1,+2] ; y: iPhi ~sector",15,30);
665 
666  hA[83]=new L2Histo(83,"diJet phi1-phi2 (accepted_EE); x: iPhi1 ~sector ; y: iPhi2 ~sector ",30,30);
667 
668  hA[84]=new L2Histo(84,"diJet1 Et (accepted_EE); x: Et (GeV)", 60);
669  hA[85]=new L2Histo(85,"diJet2 Et (accepted_EE); x: Et (GeV)", 60);
670 
671  hA[86]=new L2Histo(86,"diJet1 eta (accepted_EE); x: i Eta [-1,+2]", 15);
672  hA[87]=new L2Histo(87,"diJet2 eta (accepted_EE); x: i Eta [-1,+2]", 15);
673  hA[88]=new L2Histo(88,"diJet1 phi (accepted_EE); x: iPhi ~sector", 30);
674  hA[89]=new L2Histo(89,"diJet2 phi (accepted_EE); x: iPhi ~sector", 30);
675  hA[90]=new L2Histo(90,"diJet delZeta (accepted_EE); x: delta zeta (rad*10)", MxPhiRad10);
676  hA[91]=new L2Histo(91,"diJet delZeta vs. eta1 (accepted_EE); x: iEta1 [-1,+2] ; y: delta zeta (rad*10)",15, MxPhiRad10);
677  hA[92]=new L2Histo(92,"diJet eta2 vs. eta1 (accepted_EE); x: iEta1 [-1,+2] ;x: iEta2 [-1,+2] ",15,15);
678  hA[93]=new L2Histo(93,"diJet delZeta vs. avrPhi (accepted_EE); x: (iphi1+iphi2)/2 (12 deg/bin); y: delta zeta (rad*10)",30, MxPhiRad10);
679  hA[94]=new L2Histo(94,"total Et diJet (accepted_EE); x: Et (GeV)", 60);
680 
681 // Di-Jet_EB accepted ... all these histos #100-#114 added ... wwj 2/10
682  hA[100]=new L2Histo(100,"Et of Jet1 vs. Jet2 (accepted_EB); x: Jet1/GeV ; Jet2/GeV",30,30);
683  hA[101]=new L2Histo(101,"diJet1 eta-phi (accepted_EB); x: iEta [-1,+2] ; y: iPhi ~sector ",15,30);
684  hA[102]=new L2Histo(102,"diJet2 eta-phi (accepted_EB); x: iEta [-1,+2] ; y: iPhi ~sector",15,30);
685 
686  hA[103]=new L2Histo(103,"diJet phi1-phi2 (accepted_EB); x: iPhi1 ~sector ; y: iPhi2 ~sector ",30,30);
687 
688  hA[104]=new L2Histo(104,"diJet1 Et (accepted_EB); x: Et (GeV)", 60);
689  hA[105]=new L2Histo(105,"diJet2 Et (accepted_EB); x: Et (GeV)", 60);
690 
691  hA[106]=new L2Histo(106,"diJet1 eta (accepted_EB); x: i Eta [-1,+2]", 15);
692  hA[107]=new L2Histo(107,"diJet2 eta (accepted_EB); x: i Eta [-1,+2]", 15);
693  hA[108]=new L2Histo(108,"diJet1 phi (accepted_EB); x: iPhi ~sector", 30);
694  hA[109]=new L2Histo(109,"diJet2 phi (accepted_EB); x: iPhi ~sector", 30);
695  hA[110]=new L2Histo(110,"diJet delZeta (accepted_EB); x: delta zeta (rad*10)", MxPhiRad10);
696  hA[111]=new L2Histo(111,"diJet delZeta vs. eta1 (accepted_EB); x: iEta1 [-1,+2] ; y: delta zeta (rad*10)",15, MxPhiRad10);
697  hA[112]=new L2Histo(112,"diJet eta2 vs. eta1 (accepted_EB); x: iEta1 [-1,+2] ;x: iEta2 [-1,+2] ",15,15);
698  hA[113]=new L2Histo(113,"diJet delZeta vs. avrPhi (accepted_EB); x: (iphi1+iphi2)/2 (12 deg/bin); y: delta zeta (rad*10)",30, MxPhiRad10);
699  hA[114]=new L2Histo(114,"total Et diJet (accepted_EB); x: Et (GeV)", 60);
700 
701 // Di-Jet_BB accepted ... all these histos #120-#134 added ... wwj 2/10
702  hA[120]=new L2Histo(120,"Et of Jet1 vs. Jet2 (accepted_BB); x: Jet1/GeV ; Jet2/GeV",30,30);
703  hA[121]=new L2Histo(121,"diJet1 eta-phi (accepted_BB); x: iEta [-1,+2] ; y: iPhi ~sector ",15,30);
704  hA[122]=new L2Histo(122,"diJet2 eta-phi (accepted_BB); x: iEta [-1,+2] ; y: iPhi ~sector",15,30);
705 
706  hA[123]=new L2Histo(123,"diJet phi1-phi2 (accepted_BB); x: iPhi1 ~sector ; y: iPhi2 ~sector ",30,30);
707 
708  hA[124]=new L2Histo(124,"diJet1 Et (accepted_BB); x: Et (GeV)", 60);
709  hA[125]=new L2Histo(125,"diJet2 Et (accepted_BB); x: Et (GeV)", 60);
710 
711  hA[126]=new L2Histo(126,"diJet1 eta (accepted_BB); x: i Eta [-1,+2]", 15);
712  hA[127]=new L2Histo(127,"diJet2 eta (accepted_BB); x: i Eta [-1,+2]", 15);
713  hA[128]=new L2Histo(128,"diJet1 phi (accepted_BB); x: iPhi ~sector", 30);
714  hA[129]=new L2Histo(129,"diJet2 phi (accepted_BB); x: iPhi ~sector", 30);
715  hA[130]=new L2Histo(130,"diJet delZeta (accepted_BB); x: delta zeta (rad*10)", MxPhiRad10);
716  hA[131]=new L2Histo(131,"diJet delZeta vs. eta1 (accepted_BB); x: iEta1 [-1,+2] ; y: delta zeta (rad*10)",15, MxPhiRad10);
717  hA[132]=new L2Histo(132,"diJet eta2 vs. eta1 (accepted_BB); x: iEta1 [-1,+2] ;x: iEta2 [-1,+2] ",15,15);
718  hA[133]=new L2Histo(133,"diJet delZeta vs. avrPhi (accepted_BB); x: (iphi1+iphi2)/2 (12 deg/bin); y: delta zeta (rad*10)",30, MxPhiRad10);
719  hA[134]=new L2Histo(134,"total Et diJet (accepted_BB); x: Et (GeV)", 100);
720 
721 }
722 
723 //=======================================
724 //=======================================
725 void
726 L2jetAlgo2012::clearEvent(){
727  mAccept=false;
728  memset(eve_patchEne,0,sizeof(eve_patchEne));
729  memset(eve_phiEne,0,sizeof(eve_phiEne));
730  eve_Jet[0]->clear();
731  eve_Jet[1]->clear();
732 }
733 
734 
735 
736 //=======================================
737 //=======================================
738 int
739 L2jetAlgo2012::projectAdc(const HitTower1 *hit,const int hitSize,
740  unsigned short *phiBin, unsigned short *patchBin,
741  L2Histo *hHot ){
742  int tmpNused=0; /* counts mapped & used ADC channels */
743  short rdo;
744  float low_noise_et;
745 
746  for(int i=0;i< hitSize;i++,hit++) {
747  rdo=hit->rdo;
748  low_noise_et=hit->low_noise_et;
749  eve_patchEne[patchBin[rdo]]+=low_noise_et;
750  eve_phiEne[phiBin[rdo]]+=low_noise_et;
751  tmpNused++;
752 
753  if(low_noise_et > par_hotTwEtThr) hHot->fill(rdo);//do we want rdo or i?
754 
755  }
756  return tmpNused;
757 }
758 
759 
760 //========================================
761 //========================================
762 float L2jetAlgo2012::true2Dscan(){
763  //implament a true 2-D scan to test the performance of the current algo
764  int iphiBinEdge = 0;
765  int ietaBinEdge = 0;
766  float maxPatchEt = 0.0;
767  float sumTot = 0.0;
768 
769  float *totalCaloEt=eve_phiEne;
770  float eneA[cl2jetMaxEtaBins];
771  float secondPatchArray[cl2jetMaxPhiBins][cl2jetMaxEtaBins-cl2jet_par_mxEtaBin+1];
772  memset(secondPatchArray,0,sizeof(secondPatchArray));
773 
774  int i;//loops over all phi bin. phi edge
775  int j;//collapses phi bins
776  int j1;//phi wrapped phi bin
777  int j2;//loops over eta bins and sums energy
778  int k;//loops over eta bins. eta edge
779 
780  for(i=0;i<cl2jetMaxPhiBins;i++)
781  {
782  sumTot+=totalCaloEt[i];
783  memset(eneA,0,sizeof(eneA));
784  for(j=0;j<cl2jet_par_mxPhiBin;j++)
785  {
786  j1 = (i+j)%cl2jetMaxPhiBins;//must wrap up phi
787  float *patchEneA=eve_patchEne+(j1*cl2jetMaxEtaBins);
788  for(j2=0;j2<cl2jetMaxEtaBins;j2++,patchEneA++)
789  {
790  eneA[j2]+=*patchEneA;
791  }
792  }
793  float *eneAp=eneA;
794  for(k=0;k<cl2jetMaxEtaBins-cl2jet_par_mxEtaBin+1;k++,eneAp++)
795  {
796  float sum=eneAp[0]+eneAp[1]+eneAp[2]+eneAp[3]+eneAp[4];
797  secondPatchArray[i][k]=sum;
798  if(maxPatchEt>sum) continue;
799  //float *eneSave=eneA;
800  maxPatchEt=sum;
801  iphiBinEdge=i;
802  ietaBinEdge=k;
803  }
804  }
805 
806  eve_Jet[0]->iphiBin=iphiBinEdge;
807  eve_Jet[0]->ietaBin=ietaBinEdge;
808  eve_Jet[0]->iene=maxPatchEt;
809 
810  //now find second highest patch
811  int iphiBinEdge2=-1;
812  int ietaBinEdge2=-1;
813  float maxPatchEt2=0.0;
814  char doWrap=0;
815  int a1=iphiBinEdge-par_minPhiBinDiff;
816  int a2=iphiBinEdge+par_minPhiBinDiff;
817  if (a1<0) { a1+=cl2jetMaxPhiBins; doWrap+=1; }
818  if (a2>=cl2jetMaxPhiBins) { a2-=cl2jetMaxPhiBins; doWrap+=2; }
819 
820  int b;
821  int b1;
822  if (!doWrap)//masking range does not wrap up
823  {
824  for(b=0;b<cl2jetMaxPhiBins;b++)
825  {
826  if(b>=a1 && b<=a2) continue;
827  for(b1=0;b1<cl2jetMaxEtaBins-cl2jet_par_mxEtaBin+1;b1++)
828  {
829  if(maxPatchEt2>secondPatchArray[b][b1]) continue;
830  maxPatchEt2=secondPatchArray[b][b1];
831  iphiBinEdge2=b;
832  ietaBinEdge2=b1;
833  }
834  }
835  }
836  else
837  {//masking range wraps up
838  for(b=a2;b<a1;b++)
839  {
840  for(b1=0;b1<cl2jetMaxEtaBins-cl2jet_par_mxEtaBin+1;b1++)
841  {
842  if(maxPatchEt2>secondPatchArray[b][b1]) continue;
843  maxPatchEt2=secondPatchArray[b][b1];
844  iphiBinEdge2=b;
845  ietaBinEdge2=b1;
846  }
847  }
848  }
849 
850  eve_Jet[1]->iphiBin=iphiBinEdge2;
851  eve_Jet[1]->ietaBin=ietaBinEdge2;
852  eve_Jet[1]->iene=maxPatchEt2;
853 
854  return sumTot;
855 }
856 
857 //========================================
858 //========================================
859 void L2jetAlgo2012::weightedEtaPhi(int iK){
860  L2Jet *J=eve_Jet[iK];
861 
862  // empty calo protection
863  if(J->iene<=0.1)
864  {
865  J->fphiBin=J->iphiBin+.333;
866  J->fetaBin=J->ietaBin+.333;
867  J->rmsPhiBin=0.0;
868  J->rmsEtaBin=0.0;
869  return;
870  }
871 
872  int iphi0=J->iphiBin;
873  int ieta0=J->ietaBin;
874  float iene0=J->iene;
875 
876  //variables for weighted phi
877  float sumY=0.0;
878  float sum1;//wrk variable
879 
880  //these variables for phi RMS:
881  double sumYY = 0.0;
882  int nEnePhi=0; // counts bins with energy for RMS computation
883 
884  //variable for weighted eta
885  float sumX=0.0;
886 
887  //these variables for eta RMS:
888  double sumXX = 0.0;
889  int nEneEta = 0;
890 
891  float etaEneA[cl2jet_par_mxEtaBin];
892  memset(etaEneA,0,sizeof(etaEneA));
893 
894  int ix,iy;
895  // pick 5x5 eta-phi patch and calculate weighted phi==Y
896  for(iy=iphi0;iy<iphi0+cl2jet_par_mxPhiBin;iy++)
897  {
898  int jy=iy % cl2jetMaxPhiBins;
899  float *patchEneA=eve_patchEne+(jy*cl2jetMaxEtaBins);// phi bins are consecutive
900 
901  sum1=0.0;
902  for(ix=ieta0;ix<ieta0+cl2jet_par_mxEtaBin;ix++)
903  {
904  sum1+=patchEneA[ix];
905  etaEneA[ix-ieta0]+=patchEneA[ix];
906  }
907  if(sum1>0) nEnePhi++;
908  sumY+=sum1*iy;
909  sumYY+=sum1*iy*iy;
910  }
911 
912  // pick 5x5 eta-phi patch and calculate weighted eta==X
913  int jx;
914  for(jx=ieta0;jx<ieta0+cl2jet_par_mxEtaBin;jx++)
915  {
916  if(etaEneA[jx-ieta0]>0) nEneEta++;
917  sumX+=jx*etaEneA[jx-ieta0];
918  sumXX+=(jx*jx*etaEneA[jx-ieta0]);
919  }
920 
921  //calculate phi mean and rms
922  float phiMean= 1.*sumY/iene0;//was phiSum
923  float fphiBinMax=0.5 +phiMean;
924  float rmsPhi=0.0; // for single bin RMS is 0
925  if(nEnePhi>1)
926  { //here I calculate RMS
927  rmsPhi=sqrt(nEnePhi/(nEnePhi-1.) *(sumYY/iene0 - phiMean*phiMean));
928  }
929  if(fphiBinMax>cl2jetMaxPhiBins) fphiBinMax-=cl2jetMaxPhiBins;
930 
931  //calculate eta mean and rms
932  float etaMean = 1.*sumX/iene0;
933  float fetaBin=0.5 + etaMean;
934  float rmsEta=0.0; // for single bin RMS is 0
935  if(nEneEta>1)
936  { //here I calculate RMS
937  rmsEta=sqrt(nEneEta/(nEneEta-1.) *(sumXX/iene0 - etaMean*etaMean));
938  }
939 
940  J->fphiBin=fphiBinMax;
941  J->fetaBin=fetaBin;
942  J->rmsPhiBin=rmsPhi;
943  J->rmsEtaBin=rmsEta;
944 
945 }
946 
947 //========================================
948 //========================================
949 void
950 L2jetAlgo2012:: dumpPatchEneA(){
951  // dump L2 array with energy
952  int ix,iy;
953  for(iy=0;iy<cl2jetMaxPhiBins;iy++) {
954  float *patchEneA=eve_patchEne+(iy*cl2jetMaxEtaBins);// phi bins are consecutive
955 
956  for(ix=0;ix<cl2jetMaxEtaBins;ix++,patchEneA++){
957  printf(" %6f",*patchEneA);
958  }
959  printf(" iPhi=%d\n",iy);
960  }
961 
962 }
963 
964 
965 //=======================================
966 //=======================================
967 void
968 L2jetAlgo2012::computeE(int token){
969 }
970 //=======================================
971 //=======================================
972 void
973 L2jetAlgo2012::finishRunHisto(){
974  // auxialiary operations on histograms at the end of the run
975 
976  const int *data20=hA[20]->getData();
977  const int *data30=hA[30]->getData();
978 
979  int bHotSum=1,bHotId=-1;
980  int eHotSum=1;
981 
982  const L2EmcDb2012::EmcCDbItem *xE=mDb->getByIndex(402), *xB=mDb->getByIndex(402);
983  int i;
984  for(i=0; i<EmcDbIndexMax; i++) {
985  const L2EmcDb2012::EmcCDbItem *x=mDb->getByIndex(i);
986  if(mDb->isEmpty(x)) continue;
987  if (mDb->isBTOW(x) ) {
988  int softId=atoi(x->tube+2);
989  int ieta= (x->eta-1);
990  int iphi= (x->sec-1)*10 + x->sub-'a' ;
991  hA[21]->fillW(softId,data20[x->rdo]);
992  hA[22]->fillW(ieta, iphi,data20[x->rdo]);
993  if(bHotSum<data20[x->rdo]) {
994  bHotSum=data20[x->rdo];
995  bHotId=softId;
996  xB=x;
997  }
998  }// end of BTOW
999  else if (mDb->isETOW(x) ) {
1000  int ihard=x->chan+(x->crate-1)*128;
1001  int ieta= 12-x->eta;
1002  int iphi= (x->sec-1)*5 + x->sub-'A' ;
1003  hA[31]->fillW(ihard,data30[x->rdo]);
1004  hA[32]->fillW(ieta, iphi,data30[x->rdo]);
1005  if(eHotSum<data30[x->rdo]) {
1006  eHotSum=data30[x->rdo];
1007  xE=x;
1008  }
1009 
1010  }// end of BTOW
1011  }
1012  if (mLogFile){
1013  fprintf(mLogFile,"L2jet::finishRun()\n");
1014  fprintf(mLogFile,"#BTOW_hot tower _candidate_ (bHotSum=%d) :, softID %d , crate %d , chan %d , name %s\n",bHotSum,bHotId,xB->crate,xB->chan,xB->name);
1015  fprintf(mLogFile,"#ETOW_hot tower _candidate_ (eHotSum=%d) :, name %s , crate %d , chan %d\n",eHotSum,xE->name,xE->crate,xE->chan);
1016  }
1017 
1018 }
1019 
1020 
1021 /**********************************************************************
1022  $Log: L2jetAlgo2012.cxx,v $
1023  Revision 1.4 2012/03/21 18:18:03 jml
1024  got rid of printfs from 2012 files
1025 
1026  Revision 1.3 2011/10/19 16:12:11 jml
1027  more 2012 stuff
1028 
1029  Revision 1.2 2011/10/19 15:39:44 jml
1030  2012
1031 
1032  Revision 1.1 2011/10/18 15:11:43 jml
1033  adding 2012 algorithms
1034 
1035  Revision 1.1 2010/04/17 16:42:14 pibero
1036  *** empty log message ***
1037 
1038  Revision 1.1 2007/11/19 22:18:27 balewski
1039  most L2algos provide triggerID's
1040 
1041  Revision 1.8 2007/11/14 03:58:14 balewski
1042  cleanup of common timing measurement
1043 
1044  Revision 1.7 2007/11/13 23:06:07 balewski
1045  toward more unified L2-algos
1046 
1047  Revision 1.6 2007/11/13 00:12:36 balewski
1048  added offline triggerID, take1
1049 
1050  Revision 1.5 2007/11/08 04:02:31 balewski
1051  run on l2ana as well
1052 
1053  Revision 1.4 2007/11/02 17:43:08 balewski
1054  cleanup & it started to work w/ L2upsilon
1055 
1056  Revision 1.3 2007/11/02 03:03:47 balewski
1057  modified L2VirtualAlgo
1058 
1059  Revision 1.2 2007/10/25 02:07:02 balewski
1060  added L2upsilon & binary event dump
1061 
1062  Revision 1.1 2007/10/11 00:33:19 balewski
1063  L2algo added
1064 
1065  Revision 1.5 2006/03/28 19:46:49 balewski
1066  ver16b, in l2new
1067 
1068  Revision 1.4 2006/03/11 17:08:33 balewski
1069  now CVS comments should work
1070 
1071 */
1072 
1073