StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
L2etowCalAlgo12.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: L2etowCalAlgo12.cxx,v 1.5 2012/03/21 18:18:03 jml Exp $
10  \author Jan Balewski, MIT, 2008
11  *****************************************************
12  Descripion:
13  calibrates Endcap towers, result is used by other L2-algos
14  *****************************************************/
15 
16 
17 #ifdef IS_REAL_L2 //in l2-ana environment
18  #include "../L2algoUtil/L2EmcDb2012.h"
19  #include "../L2algoUtil/L2Histo.h"
20 #else
21  #include "L2EmcDb2012.h"
22  #include "L2Histo.h"
23  #include "L2EmcGeom2012.h"
24 #endif
25 
26 #include "L2etowCalAlgo12.h"
27 
28 
29 //=================================================
30 //=================================================
31 L2etowCalAlgo12::L2etowCalAlgo12(const char* name, const char *uid, L2EmcDb2012* db, L2EmcGeom2012 *geoX, char* outDir, int resOff) : L2VirtualAlgo2012( name, uid, db, outDir, false, true, resOff) {
32  /* called once per days
33  all memory allocation must be done here
34  */
35 
36  mGeom=geoX;
37  if (!mGeom)
38  criticalError("L2etowCalAlgo is broken -- can't find geom.");
39 
40  setMaxHist(32);
41  createHisto();
42 
43  // initialize ETOW-Calibrated-data to zero
44  int k;
45  for(k=0;k<L2eventStream2012::mxToken;k++){
46  L2EtowCalibData12 & etowCalibData=globL2eventStream2012.etow[k];
47  etowCalibData.nInputBlock=0;
48  etowCalibData.hitSize=0;
49  }
50  }
51 
52 /* ========================================
53  ======================================== */
54 int
55 L2etowCalAlgo12::initRunUser( int runNo, int *rc_ints, float *rc_floats) {
56 
57  par_adcMask=(unsigned short)(-0x10); //not in .h yet
58  par_pedOff=0x10/2; //hex must match with above //not in .h yet
59 
60  // unpack params from run control GUI
61  par_dbg = rc_ints[0];
62  par_gainType = rc_ints[1];
63  par_nSigPed = rc_ints[2];
64 
65  par_twEneThres = rc_floats[0];
66  par_hotEtThres = rc_floats[1];;
67 
68  // verify consistency of input params
69  int kBad=0;
70  kBad+=0x00001 * (par_gainType<kGainZero || par_gainType>kGainOffline);
71  kBad+=0x00002 * (par_nSigPed<2 || par_nSigPed>5);
72  kBad+=0x00004 * (par_twEneThres<0.1 || par_twEneThres>1.5);
73 
74  if (mLogFile) {
75  fprintf(mLogFile,"L2%s algorithm initRun(R=%d), compiled: %s , %s\n params:\n",getName(),mRunNumber,__DATE__,__TIME__);
76  fprintf(mLogFile," - use ETOW=%d, gain Ideal=%d or Offline=%d, debug=%d\n",
77  par_gainType>=kGainIdeal, par_gainType==kGainIdeal, par_gainType==kGainOffline, par_dbg);
78  fprintf(mLogFile," - thresholds: ADC-ped> %d*sigPed .AND. energy>%.2f GeV \n", par_nSigPed, par_twEneThres);
79 
80  fprintf(mLogFile," - hot tower thresholds: ET/GeV=%.2f\n",par_hotEtThres);
81  fprintf(mLogFile,"initRun() params checked for consistency, Error flag=0x%04x\n",kBad);
82  }
83 
84  //initialize the etow calibrated data all zeroes.
85  int k;
86  for(k=0;k<L2eventStream2012::mxToken;k++){
87  L2EtowCalibData12 & etowCalibData=globL2eventStream2012.etow[k];
88  etowCalibData.nInputBlock=0;
89  etowCalibData.hitSize=0;
90  }
91 
92  if(kBad) return kBad;
93 
94  // clear content of all histograms
95  int i;
96  for (i=0; i<mxHA;i++) if(hA[i])hA[i]->reset();
97 
98  // upadate title of histos
99  char txt[1000];
100  sprintf(txt,"ETOW tower, E_T>%.2f GeV (input); x: ETOW RDO index=chan*6+fiber; y: counts",par_hotEtThres);
101  hA[10]->setTitle(txt);
102 
103  sprintf(txt,"ETOW tower, Et>%.2f GeV (input); x: ETOW softID=i#phi+60*i#eta",par_hotEtThres);
104  hA[11]->setTitle(txt);
105  sprintf(txt,"ETOW tower, Et>%.2f GeV (input); x: eta bin, [+2,+1]; y: phi bin ~ TPC sector",par_hotEtThres);
106  hA[12] ->setTitle(txt);
107 
108  sprintf(txt,"#ETOW towers / event , Et>%.2f GeV; x: # ETOW towers; y: counts",par_hotEtThres);
109  hA[14] ->setTitle(txt);
110 
111  // re-caluclate geometry properties
112  mGeom->etow.clear();
113  int nT=0; /* counts # of unmasekd towers */
114  int nTg=0; /* counts # of reasonable calibrated towers */
115  int nEneThr=0, nPedThr=0; //ETOW count # of towers above & below threshold
116  if(par_gainType>=kGainIdeal) // this disables the whole loop below
117  for(i=0; i<EmcDbIndexMax; i++) {
118  const L2EmcDb2012::EmcCDbItem *x=mDb->getByIndex(i);
119  if(mDb->isEmpty(x)) continue; /* dropped not mapped channels */
120  /*....... E N D C A P .................*/
121  if (!mDb->isETOW(x) ) continue; /* drop if not ETOW */
122  if(x->fail) continue; /* dropped masked channels */
123  if(x->gain<=0) continue; /* dropped uncalibrated towers , tmp */
124  nT++;
125 
126  float adcThres=x->ped+par_nSigPed* fabs(x->sigPed);
127  float otherThr=x->ped+par_twEneThres*x->gain;
128 
129  if(adcThres<otherThr) { //use energy threshold if higher
130  adcThres=otherThr;
131  nEneThr++;
132  } else {
133  nPedThr++;
134  }
135 
136  /* use rdo index to match RDO order in the ADC data banks */
137  if(x->eta<=0 || x->eta>EtowGeom::mxEtaBin) return -90;
138  int ietaTw= (x->eta-1); /* correct */
139 
140  // use ideal gains for now, hardcoded
141  if (par_gainType!=kGainIdeal) return -102;
142  mGeom->etow.gain2Ene_rdo[x->rdo]=mGeom->etow.idealGain2Ene[ietaTw];
143  mGeom->etow.gain2ET_rdo[x->rdo]=mGeom->getIdealAdc2ET();
144 
145  mGeom->etow.thr_rdo[x->rdo]=(int) (adcThres);
146  mGeom->etow.ped_rdo[x->rdo]=(int) (x->ped);
147  mGeom->etow.ped_shifted_rdo[x->rdo]=(unsigned short)(par_pedOff - x->ped);
148  nTg++;
149  }
150 
151  if (mLogFile) {
152  fprintf(mLogFile," found towers working=%d calibrated=%d, based on ASCII DB\n",nT,nTg);
153  fprintf(mLogFile," thresh defined by energy=%d or NsigPed=%d \n",nEneThr, nPedThr);
154  }
155 
156  return 0; //OK
157 
158 
159 }
160 
161 /* ========================================
162  ======================================== */
163 void
164 L2etowCalAlgo12::calibrateEtow(int token, int eemcIn, unsigned short *rawAdc){
165  // Etow calibration is a special case, must have one exit at the end
166 
167  computeStart();
168  token&=L2eventStream2012::tokenMask; // only to protect against a bad token, Gerard's trick
169 
170  //...... now token is valid ........
171  L2EtowCalibData12 & etowCalibData=globL2eventStream2012.etow[token];
172  // clear data for this token from previous event
173  etowCalibData.nInputBlock++;
174  etowCalibData.hitSize=0;
175 
176  int nTower=0; /* counts mapped & used ADC channels */
177  int nHotTower=0;
178  if(eemcIn && par_gainType>kGainZero) { // EVEVEVEVEVE
179  // ............process this event ...............
180  short rdo;
181  int adc; // pedestal subtracted
182  int low_noise_adc;
183  float et;
184  float low_noise_et; // bit-chopped
185  unsigned short *thr=mGeom->etow.thr_rdo;
186  unsigned short *ped=mGeom->etow.ped_rdo;
187  unsigned short *ped_shifted=mGeom->etow.ped_shifted_rdo;
188  float *gain2ET=mGeom->etow.gain2ET_rdo;
189  float *gain2Ene=mGeom->etow.gain2Ene_rdo;
190  HitTower1 *hit=etowCalibData.hit;
191  for(rdo=0; rdo<EtowGeom::mxRdo; rdo++){
192  if(rawAdc[rdo]<thr[rdo])continue;
193  if(nTower>=L2EtowCalibData12::mxListSize) break; // overflow protection
194  adc=rawAdc[rdo]-ped[rdo]; //did NOT correct for common pedestal noise - bad for the jet finder
195  et=adc/gain2ET[rdo];
196  low_noise_adc=(rawAdc[rdo]+ped_shifted[rdo]) & par_adcMask;//note that ped_shifted is defined as pedOffset-ped.
197  low_noise_et=low_noise_adc/gain2ET[rdo];
198  hit->rdo=rdo;
199  hit->adc=adc;
200  hit->et=et;
201  hit->low_noise_et=low_noise_et;
202  hit->ene=adc/gain2Ene[rdo];
203  hit++;
204  nTower++;
205  // only monitoring
206  if(et >par_hotEtThres) {
207  hA[10]->fill(rdo);
208  nHotTower++;
209  }
210  }
211  etowCalibData.hitSize=nTower;
212 
213  // QA histos
214  hA[13]->fill(nTower);
215  hA[14]->fill(nHotTower);
216  if(nTower>=L2EtowCalibData12::mxListSize) mhN->fill(5); // was overflow
217  } // EVEVEVEVEVE
218 
219  // debugging should be off for any time critical computation
220  if(par_dbg>0){
221  LOG(DBG,"L2-%s-compute: set adcL size=%d\n",getName(),nTower);
222  LOG(DBG,"dbg=%s: found nTw=%d\n",getName(),nTower);
223  if(par_dbg>0) print0();
224  printCalibratedData(token);
225  }
226 
227  computeStop(token);
228 
229 }
230 /* ========================================
231  ======================================== */
232 
233 void
234 L2etowCalAlgo12::clear(int token){
235  token&=L2eventStream2012::tokenMask; // only protect against bad token, Gerard's trick
236  //...... token is valid no w ........
237  L2EtowCalibData12 & etowCalibData=globL2eventStream2012.etow[token];
238  etowCalibData.hitSize=0;
239  return;
240 }
241 
242 
243 /* ========================================
244  ======================================== */
245 void
246 L2etowCalAlgo12::computeUser(int token ){
247 
248  LOG(CRIT,"computeUser-%s FATAL CRASH\n If you see this message it means l2new is very badly misconfigured \n and L2-etow-calib algo was not executed properly\n before calling other individual L2-algos. \n\n l2new will aborted now - fix the code, Jan B.\n",getName());
249  criticalError("L2etowCalAlgo12::computeUser has been called and should not have been. Serious problem in L2");
250 
251 }
252 
253 
254 /* ========================================
255  ======================================== */
256 void
257 L2etowCalAlgo12::finishRunUser() {
258  /* called once at the end of the run
259  write do whatever you want, log-file & histo-file are still open
260  Here it seraches for hot tower, re-project histos vs. other representations
261  */
262 
263  //report mean and RMS of nonzero towers per event:
264  int mean=0, RMS=0;
265  hA[13]->findMean(&mean,&RMS);
266  if (mLogFile){
267  fprintf(mLogFile,"#ETOW_nonzero_towers_per_event: mean=%d, rms=%d\n",mean,RMS);
268  }
269 
270 
271  int eHotSum=1,eHotId=-1;
272  const int *data20=hA[10]->getData();
273  const L2EmcDb2012::EmcCDbItem *xE=0; //mDb->getByIndex(502); // some wired default?
274 
275  int i;
276  for(i=0; i<EmcDbIndexMax; i++) {
277  const L2EmcDb2012::EmcCDbItem *x=mDb->getByIndex(i);
278  if(mDb->isEmpty(x)) continue;
279  if (!mDb->isETOW(x) ) continue;
280  int ieta= (x->eta-1);
281  int iphi= (x->sec-1)*EtowGeom::mxSubs + x->sub-'A' ;
282  int softId= iphi+EtowGeom::mxPhiBin*ieta;
283  hA[11]->fillW(softId,data20[x->rdo]);
284  hA[12]->fillW(ieta, iphi,data20[x->rdo]);
285  if(eHotSum<data20[x->rdo]) {
286  eHotSum=data20[x->rdo];
287  eHotId=softId;
288  xE=x;
289  }
290  }
291 
292  int par_nHotThresh=20;
293  if (mLogFile && eHotSum>par_nHotThresh){
294  fprintf(mLogFile,"#ETOW_hot tower _candidate_ (eHotSum=%d of %d eve) :, softID %d , crate %d , chan %d , name %s\n",eHotSum,mEventsInRun,eHotId,xE->crate,xE->chan,xE->name);
295  }
296 
297  //...... QA tokens .....
298  int tkn1=99999, tkn2=0; // min/max token
299  int nTkn=0;
300  int tkn3=-1, nTkn3=-1; // most often used token
301 
302  int k;
303  for(k=0;k<L2eventStream2012::mxToken;k++){
304  L2EtowCalibData12 & etowCalibData=globL2eventStream2012.etow[k];
305  if(etowCalibData.nInputBlock==0) continue;
306  hA[1]->fillW(k,etowCalibData.nInputBlock);
307  if(nTkn3<etowCalibData.nInputBlock){
308  nTkn3=etowCalibData.nInputBlock;
309  tkn3=k;
310  }
311 
312  nTkn++;
313  if(tkn1>k) tkn1=k;
314  if(tkn2<k) tkn2=k;
315  }
316  if (mLogFile){
317  fprintf(mLogFile,"#ETOW_token_QA: _candidate_ hot token=%d used %d for %d events, token range [%d, %d], used %d tokens\n",tkn3,nTkn3,mEventsInRun,tkn1,tkn2,nTkn);
318  }
319 
320 }
321 
322 
323 //=======================================
324 //=======================================
325 void
326 L2etowCalAlgo12::createHisto() {
327  memset(hA,0,mxHA*sizeof(L2Histo*));
328  //token related spectra
329  hA[1]=new L2Histo(1,"L2-etow-calib: seen tokens; x: token value; y: events ",L2eventStream2012::mxToken);
330 
331  // ETOW raw spectra (zz 4 lines)
332  hA[10]=new L2Histo(10,"etow hot tower 1", EtowGeom::mxRdo); // title upadted in initRun
333  hA[11]=new L2Histo(11,"etow hot tower 2", EtowGeom::mxRdo); // title upadted in initRun
334  hA[12]=new L2Histo(12,"etow hot tower 3", EtowGeom::mxEtaBin,EtowGeom::mxPhiBin); // title upadted in initRun
335  hA[13]=new L2Histo(13,"ETOW #tower w/ energy /event; x: # ETOW towers; y: counts", 100);
336  hA[14]=new L2Histo(14,"# hot towers/event", 30);
337 
338 }
339 
340 
341 
342 /* ========================================
343  ======================================== */
344 void
345 L2etowCalAlgo12::print0(){ // full raw input ADC array
346  // empty
347  }
348 
349 
350 /**********************************************************************
351  $Log: L2etowCalAlgo12.cxx,v $
352  Revision 1.5 2012/03/21 18:18:03 jml
353  got rid of printfs from 2012 files
354 
355  Revision 1.4 2011/10/19 16:12:10 jml
356  more 2012 stuff
357 
358  Revision 1.3 2011/10/19 15:58:06 jml
359  more compile offline
360 
361  Revision 1.2 2011/10/19 15:39:42 jml
362  2012
363 
364  Revision 1.1 2011/10/18 15:11:41 jml
365  adding 2012 algorithms
366 
367  Revision 1.2 2010/04/18 02:53:35 pibero
368  Fixed memset bug.
369 
370  Revision 1.1 2010/04/17 17:27:31 pibero
371  *** empty log message ***
372 
373  Revision 1.4 2008/11/12 00:16:40 rcorliss
374  add low_noise_et to BTOW/ETOW calibration
375 
376  Revision 1.3 2008/02/01 00:16:40 balewski
377  add mxListSize to BTOW/ETOW calibration
378 
379  Revision 1.2 2008/01/30 21:56:40 balewski
380  E+B high-enery-filter L2-algo fuly functional
381 
382  Revision 1.1 2008/01/30 00:47:16 balewski
383  Added L2-Etow-calib
384 
385 
386 
387 */
388 
389