StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
L2btowCalAlgo08.cxx
1 #include <stdio.h>
2 #include <assert.h>
3 #include <string.h>
4 #include <stdlib.h>
5 #include <time.h>
6 #include <math.h>
7 
8 /*********************************************************
9  $Id: L2btowCalAlgo08.cxx,v 1.5 2008/02/01 00:16:40 balewski Exp $
10  \author Jan Balewski, MIT, 2008
11  *****************************************************
12  Descripion:
13  calibrates Barrel towers, result is used by other L2-algo
14  *****************************************************/
15 
16 
17 #ifdef IS_REAL_L2 //in l2-ana environment
18  #include "../L2algoUtil/L2EmcDb.h"
19  #include "../L2algoUtil/L2Histo.h"
20 #else
21  #include "L2EmcDb.h"
22  #include "L2Histo.h"
23  #include "L2EmcGeom.h"
24 #endif
25 
26 #include "L2btowCalAlgo08.h"
27 
28 /************************
29 the _only_ copy of L2event for all L2-algos
30 **************************/
31 L2eventStream2008 globL2eventStream2008;
32 /*************************/
33 
34 //=================================================
35 //=================================================
36 L2btowCalAlgo08::L2btowCalAlgo08(const char* name, L2EmcDb* db, L2EmcGeom *geoX, char* outDir) : L2VirtualAlgo2008( name, db, outDir) {
37  /* called one per days
38  all memory allocation must be done here
39  */
40 
41  geom=geoX; assert(geom);
42 
43  setMaxHist(32);
44  createHisto();
45 
46  // initilalize BTOW-Calibrated-data
47  int k;
48  for(k=0;k<L2eventStream2008::mxToken;k++){
49  L2BtowCalibData08 & btowCalibData=globL2eventStream2008.btow[k];
50  btowCalibData.nInputBlock=0;
51  }
52  }
53 
54 /* ========================================
55  ======================================== */
56 int
57 L2btowCalAlgo08::initRunUser( int runNo, int *rc_ints, float *rc_floats) {
58 
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 BTOW=%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  if(kBad) return kBad;
85 
86  // clear content of all histograms
87  int i;
88  for (i=0; i<mxHA;i++) if(hA[i])hA[i]->reset();
89 
90  // upadate title of histos
91  char txt[1000];
92  sprintf(txt,"BTOW tower, E_T>%.2f GeV (input); x: BTOW RDO index=chan*30+fiber; y: counts",par_hotEtThres);
93  hA[10]->setTitle(txt);
94 
95  sprintf(txt,"BTOW tower, Et>%.2f GeV (input); x: BTOW softID",par_hotEtThres);
96  hA[11]->setTitle(txt);
97  sprintf(txt,"BTOW tower, Et>%.2f GeV (input); x: eta bin, [-1,+1]; y: phi bin ~ TPC sector",par_hotEtThres);
98  hA[12] ->setTitle(txt);
99 
100  sprintf(txt,"#BTOW towers / event , Et>%.2f GeV; x: # BTOW towers; y: counts",par_hotEtThres);
101  hA[14] ->setTitle(txt);
102 
103  // re-caluclate geometry properties
104  geom->btow.clear();
105  int nB=0; /* counts # of unmasekd towers */
106  int nBg=0; /* counts # of reasonable calibrated towers */
107  int nEneThr=0, nPedThr=0; //BTOW count # of towers above & below threshold
108  if(par_gainType>=kGainIdeal) // this disables the whole loop below
109  for(i=0; i<EmcDbIndexMax; i++) {
110  const L2EmcDb::EmcCDbItem *x=mDb->getByIndex(i);
111  if(mDb->isEmpty(x)) continue; /* dropped not mapped channels */
112  /*....... B A R R E L .................*/
113  if (!mDb->isBTOW(x) ) continue; /* drop if not BTOW */
114  if(x->fail) continue; /* dropped masked channels */
115  if(x->gain<=0) continue; /* dropped uncalibrated towers , tmp */
116  nB++;
117 
118  float adcThres=x->ped+par_nSigPed* fabs(x->sigPed);
119  float otherThr=x->ped+par_twEneThres*x->gain;
120  // if(strstr("01tg34",x->name)) printf("%s g=%f adcThrEne=%.1f adcThrAbovePed=%.1f rdo=%d\n",x->name,x->gain,adcThres,otherThr,x->rdo);
121 
122  if(adcThres<otherThr) { //use energy threshold if higher
123  adcThres=otherThr;
124  nEneThr++;
125  } else {
126  nPedThr++;
127  }
128 
129  /* use rdo index to match RDO order in the ADC data banks */
130  if(x->eta<=0 || x->eta>BtowGeom::mxEtaBin) return -90;
131  int ietaTw= (x->eta-1); /* correct */
132 
133  // use ideal gains for now, hardcoded
134  assert(par_gainType==kGainIdeal); // offline gains not implemented - should be changed here, Jan
135  geom->btow.gain2Ene_rdo[x->rdo]=geom->btow.idealGain2Ene[ietaTw];
136  geom->btow.gain2ET_rdo[x->rdo]=geom->getIdealAdc2ET();
137 
138  geom->btow.thr_rdo[x->rdo]=(int) (adcThres);
139  geom->btow.ped_rdo[x->rdo]=(int) (x->ped);
140  nBg++;
141  }
142 
143  if (mLogFile) {
144  fprintf(mLogFile," found towers working=%d calibrated=%d, based on ASCII DB\n",nB,nBg);
145  fprintf(mLogFile," thresh defined by energy=%d or NsigPed=%d \n",nEneThr, nPedThr);
146  }
147 
148  return 0; //OK
149 
150 
151 }
152 
153 /* ========================================
154  ======================================== */
155 void
156 L2btowCalAlgo08::calibrateBtow(int token, int bemcIn, unsigned short *rawAdc){
157  // Btow calibration is a special case, must have one exit at the end
158 
159  computeStart();
160  token&=L2eventStream2008::tokenMask; // only protect against bad token, Gerard's trick
161 
162  //...... token is valid no w ........
163  L2BtowCalibData08 & btowCalibData=globL2eventStream2008.btow[token];
164  btowCalibData.nInputBlock++;
165 
166  // clear data for this token from previous event
167  btowCalibData.hitSize=0;
168 
169  int nTower=0; /* counts mapped & used ADC channels */
170  int nHotTower=0;
171  if(bemcIn && par_gainType>kGainZero) { // EVEVEVEVEVE
172  // ............process this event ...............
173  short rdo;
174  int adc; // pedestal subtracted
175  float et;
176  unsigned short *thr=geom->btow.thr_rdo;
177  unsigned short *ped=geom->btow.ped_rdo;
178  float *gain2ET=geom->btow.gain2ET_rdo;
179  float *gain2Ene=geom->btow.gain2Ene_rdo;
180  HitTower1 *hit=btowCalibData.hit;
181  for(rdo=0; rdo<BtowGeom::mxRdo; rdo++){
182  if(rawAdc[rdo]<thr[rdo])continue;
183  if(nTower>=L2BtowCalibData08::mxListSize) break; // overflow protection
184  adc=rawAdc[rdo]-ped[rdo]; //do NOT correct for common pedestal noise
185  et=adc/gain2ET[rdo];
186  hit->rdo=rdo;
187  hit->adc=adc;
188  hit->et=et;
189  hit->ene=adc/gain2Ene[rdo];
190  hit++;
191  nTower++; //printf("nBt=%d\n",nTower);
192  // only monitoring
193  if(et >par_hotEtThres) {
194  hA[10]->fill(rdo);
195  nHotTower++;
196  }
197  }
198  btowCalibData.hitSize=nTower;
199 
200  // QA histos
201  hA[13]->fill(nTower);
202  hA[14]->fill(nHotTower);
203  if(nTower>=L2BtowCalibData08::mxListSize) mhN->fill(5); // was overflow
204 
205  } // EVEVEVEVEVE
206 
207  // debugging should be off for any time critical computation
208  if(par_dbg>0){
209  printf("L2-%s-compute: set adcL size=%d, get=%d\n",getName(),nTower,999); //tmp
210  printf("dbg=%s: found nTw=%d\n",getName(),nTower);
211  if(par_dbg>0) print0();
212  printCalibratedData(token);
213  }
214 
215  computeStop(token);
216 
217 }
218 
219 /* ========================================
220  ======================================== */
221 void
222 L2btowCalAlgo08::computeUser(int token ){
223 
224  printf("computeUser-%s FATAL CRASH\n If you see this message it means l2new is very badly misconfigured \n and L2-btow-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());
225  assert(1==2);
226 }
227 
228 
229 /* ========================================
230  ======================================== */
231 void
232 L2btowCalAlgo08::finishRunUser() { /* called once at the end of the run */
233  // do whatever you want, log-file & histo-file are still open
234 
235  // serach for hot tower, re-project histos vs. other representations
236  int bHotSum=1,bHotId=-1;
237  const int *data20=hA[10]->getData();
238  const L2EmcDb::EmcCDbItem *xB=mDb->getByIndex(402); // some wired default?
239 
240  int i;
241  for(i=0; i<EmcDbIndexMax; i++) {
242  const L2EmcDb::EmcCDbItem *x=mDb->getByIndex(i);
243  if(mDb->isEmpty(x)) continue;
244  if (!mDb->isBTOW(x) ) continue;
245  int softId=atoi(x->tube+2);
246  int ieta= (x->eta-1);
247  int iphi= (x->sec-1)*BtowGeom::mxSubs + x->sub-'a' ;
248  hA[11]->fillW(softId,data20[x->rdo]);
249  hA[12]->fillW(ieta, iphi,data20[x->rdo]);
250  if(bHotSum<data20[x->rdo]) {
251  bHotSum=data20[x->rdo];
252  bHotId=softId;
253  xB=x;
254  }
255  }
256 
257  if (mLogFile){
258  fprintf(mLogFile,"#BTOW_hot tower _candidate_ (bHotSum=%d of %d eve) :, softID %d , crate %d , chan %d , name %s\n",bHotSum,mEventsInRun,bHotId,xB->crate,xB->chan,xB->name);
259  }
260 
261  //...... QA tokens .....
262  int tkn1=99999, tkn2=0; // min/max token
263  int nTkn=0;
264  int tkn3=-1, nTkn3=-1; // most often used token
265 
266  int k;
267  for(k=0;k<L2eventStream2008::mxToken;k++){
268  L2BtowCalibData08 & btowCalibData=globL2eventStream2008.btow[k];
269  if(btowCalibData.nInputBlock==0) continue;
270  hA[1]->fillW(k,btowCalibData.nInputBlock);
271  if(nTkn3<btowCalibData.nInputBlock){
272  nTkn3=btowCalibData.nInputBlock;
273  tkn3=k;
274  }
275 
276  nTkn++;
277  if(tkn1>k) tkn1=k;
278  if(tkn2<k) tkn2=k;
279  }
280  if (mLogFile){
281  fprintf(mLogFile,"#BTOW_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);
282  }
283 
284 }
285 
286 
287 //=======================================
288 //=======================================
289 void
290 L2btowCalAlgo08::createHisto() {
291  memset(hA,0,sizeof(hA));
292  //token related spectra
293  hA[1]=new L2Histo(1,"L2-btow-calib: seen tokens; x: token value; y: events ",L2eventStream2008::mxToken);
294 
295  // BTOW raw spectra
296  hA[10]=new L2Histo(10,"btow hot tower 1", BtowGeom::mxRdo); // title upadted in initRun
297  hA[11]=new L2Histo(11,"btow hot tower 2", BtowGeom::mxRdo); // title upadted in initRun
298  hA[12]=new L2Histo(12,"btow hot tower 3", BtowGeom::mxEtaBin,BtowGeom::mxPhiBin); // title upadted in initRun
299  hA[13]=new L2Histo(13,"BTOW #tower w/ energy /event; x: # BTOW towers; y: counts", 100);
300  hA[14]=new L2Histo(14,"# hot towers/event", 100);
301 
302 }
303 
304 
305 
306 /* ========================================
307  ======================================== */
308 void
309 L2btowCalAlgo08::print0(){ // full raw input ADC array
310  // empty
311  }
312 
313 
314 /****************************************************
315  $Log: L2btowCalAlgo08.cxx,v $
316  Revision 1.5 2008/02/01 00:16:40 balewski
317  add mxListSize to BTOW/ETOW calibration
318 
319  Revision 1.4 2008/01/30 00:47:15 balewski
320  Added L2-Etow-calib
321 
322  Revision 1.3 2008/01/18 23:29:13 balewski
323  now L2result is exported
324 
325  Revision 1.2 2008/01/16 23:32:34 balewski
326  toward token dependent compute()
327 
328  Revision 1.1 2007/12/19 02:30:18 balewski
329  new L2-btow-calib-2008
330 
331  Revision 1.1 2007/11/19 22:18:25 balewski
332  most L2algos provide triggerID's
333 
334 
335 */
336 
337