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