StRoot: StTriggerUtilities/L2Emulator/L2algoUtil/L2etowCalAlgo08.cxx Source File

00001 #include <stdio.h>
00002 #include <assert.h>
00003 #include <string.h>
00004 #include <stdlib.h>
00005 #include <time.h>
00006 #include <math.h>
00007 
00008 /*********************************************************
00009   $Id: L2etowCalAlgo08.cxx,v 1.3 2008/02/01 00:16:40 balewski Exp $
00010   \author Jan Balewski, MIT, 2008 
00011  *****************************************************
00012   Descripion:
00013   calibrates Endcap towers, result is used by other L2-algos
00014  *****************************************************/
00015 
00016 
00017 #ifdef  IS_REAL_L2  //in l2-ana  environment
00018   #include "../L2algoUtil/L2EmcDb.h"
00019   #include "../L2algoUtil/L2Histo.h"
00020 #else
00021   #include "L2EmcDb.h"
00022   #include "L2Histo.h"
00023   #include "L2EmcGeom.h"
00024 #endif
00025 
00026 #include "L2etowCalAlgo08.h"
00027 
00028 
00029 //=================================================
00030 //=================================================
00031 L2etowCalAlgo08::L2etowCalAlgo08(const char* name, L2EmcDb* db, L2EmcGeom *geoX, char* outDir)  :  L2VirtualAlgo2008( name,  db,  outDir) { 
00032   /* called once per days
00033      all memory allocation must be done here
00034   */
00035 
00036   mGeom=geoX; assert(mGeom);
00037 
00038   setMaxHist(32);
00039   createHisto();
00040 
00041   // initilalize ETOW-Calibrated-data to zero
00042   int k;
00043   for(k=0;k<L2eventStream2008::mxToken;k++){
00044     L2EtowCalibData08 & etowCalibData=globL2eventStream2008.etow[k];  
00045     etowCalibData.nInputBlock=0;
00046     etowCalibData.hitSize=0;
00047   }
00048  }
00049 
00050 /* ========================================
00051   ======================================== */
00052 int 
00053 L2etowCalAlgo08::initRunUser( int runNo, int *rc_ints, float *rc_floats) {
00054 
00055  
00056   // unpack params from run control GUI
00057   par_dbg       =  rc_ints[0];
00058   par_gainType  =  rc_ints[1];
00059   par_nSigPed   =  rc_ints[2];
00060 
00061   par_twEneThres = rc_floats[0];
00062   par_hotEtThres = rc_floats[1];;
00063 
00064   // verify consistency of input params
00065   int kBad=0;
00066   kBad+=0x00001 * (par_gainType<kGainZero || par_gainType>kGainOffline);
00067   kBad+=0x00002 * (par_nSigPed<2 || par_nSigPed>5);
00068   kBad+=0x00004 * (par_twEneThres<0.1 ||  par_twEneThres>1.5);
00069 
00070   if (mLogFile) { 
00071     fprintf(mLogFile,"L2%s algorithm initRun(R=%d), compiled: %s , %s\n params:\n",getName(),mRunNumber,__DATE__,__TIME__);
00072     fprintf(mLogFile," - use ETOW=%d,  gain Ideal=%d or Offline=%d, debug=%d\n",
00073             par_gainType>=kGainIdeal, par_gainType==kGainIdeal, par_gainType==kGainOffline, par_dbg);
00074     fprintf(mLogFile," - thresholds: ADC-ped> %d*sigPed .AND. energy>%.2f GeV \n", par_nSigPed, par_twEneThres);
00075 
00076     fprintf(mLogFile," - hot tower thresholds:  ET/GeV=%.2f\n",par_hotEtThres);
00077     fprintf(mLogFile,"initRun() params checked for consistency, Error flag=0x%04x\n",kBad);
00078   }
00079   
00080   if(kBad) return kBad;
00081 
00082   // clear content of all histograms
00083   int i;
00084   for (i=0; i<mxHA;i++) if(hA[i])hA[i]->reset();
00085 
00086   // upadate title of histos
00087   char txt[1000];
00088   sprintf(txt,"ETOW tower, E_T>%.2f GeV (input); x: ETOW RDO index=chan*6+fiber; y: counts",par_hotEtThres);
00089   hA[10]->setTitle(txt);
00090   
00091   sprintf(txt,"ETOW tower, Et>%.2f GeV (input); x: ETOW softID=i#phi+60*i#eta",par_hotEtThres);
00092   hA[11]->setTitle(txt);
00093   sprintf(txt,"ETOW tower, Et>%.2f GeV (input); x: eta bin, [-1,+1];  y: phi bin ~ TPC sector",par_hotEtThres);
00094   hA[12] ->setTitle(txt);
00095   
00096   sprintf(txt,"#ETOW towers / event , Et>%.2f GeV; x: # ETOW towers; y: counts",par_hotEtThres);
00097   hA[14] ->setTitle(txt);
00098   
00099   // re-caluclate geometry properties
00100   mGeom->etow.clear(); 
00101   int nT=0; /* counts # of unmasekd towers */ 
00102   int nTg=0; /* counts # of reasonable calibrated towers */ 
00103   int nEneThr=0, nPedThr=0; //ETOW count # of towers above & below threshold
00104   if(par_gainType>=kGainIdeal)  // this disables the whole loop below
00105   for(i=0; i<EmcDbIndexMax; i++) {
00106     const L2EmcDb::EmcCDbItem *x=mDb->getByIndex(i);
00107     if(mDb->isEmpty(x)) continue;  /* dropped not mapped  channels */
00108     /*....... E N D C A P  .................*/
00109     if (!mDb->isETOW(x) ) continue; /* drop if not ETOW */
00110     if(x->fail) continue;          /* dropped masked channels */
00111     if(x->gain<=0) continue;       /* dropped uncalibrated towers , tmp */
00112     nT++;   
00113 
00114     float adcThres=x->ped+par_nSigPed* fabs(x->sigPed);
00115     float otherThr=x->ped+par_twEneThres*x->gain;
00116  
00117     if(adcThres<otherThr) { //use energy threshold if higher
00118       adcThres=otherThr;
00119       nEneThr++;
00120     } else {
00121       nPedThr++;
00122     }
00123     
00124     /* use rdo index to match RDO order in the ADC data banks */    
00125     if(x->eta<=0 || x->eta>EtowGeom::mxEtaBin) return -90;
00126     int ietaTw= (x->eta-1); /* correct */
00127 
00128     // use ideal gains for now, hardcoded
00129     assert(par_gainType==kGainIdeal); // offline gains not implemented - should be changed here, Jan
00130     mGeom->etow.gain2Ene_rdo[x->rdo]=mGeom->etow.idealGain2Ene[ietaTw];
00131     mGeom->etow.gain2ET_rdo[x->rdo]=mGeom->getIdealAdc2ET();
00132     
00133     mGeom->etow.thr_rdo[x->rdo]=(int) (adcThres);
00134     mGeom->etow.ped_rdo[x->rdo]=(int) (x->ped);
00135     nTg++;      
00136   }
00137   
00138   if (mLogFile) {
00139     fprintf(mLogFile,"  found  towers working=%d calibrated=%d, based on ASCII DB\n",nT,nTg);
00140     fprintf(mLogFile,"  thresh defined by energy=%d  or NsigPed=%d \n",nEneThr, nPedThr);
00141   }
00142 
00143   return 0; //OK
00144 
00145 
00146 }               
00147 
00148 /* ========================================
00149   ======================================== */
00150 void 
00151 L2etowCalAlgo08::calibrateEtow(int token, int eemcIn, ushort *rawAdc){
00152   // Etow calibration is a special case, must have one exit  at the end
00153 
00154   computeStart();
00155   token&=L2eventStream2008::tokenMask; // only to protect against a bad token, Gerard's trick
00156  
00157   //...... now token is valid  ........
00158   L2EtowCalibData08 & etowCalibData=globL2eventStream2008.etow[token];  
00159   // clear data for this token from previous event
00160   etowCalibData.nInputBlock++;
00161   etowCalibData.hitSize=0;
00162 
00163   int nTower=0; /* counts mapped & used ADC channels */
00164   int nHotTower=0;
00165   if(eemcIn && par_gainType>kGainZero) { // EVEVEVEVEVE
00166     // ............process this event ...............
00167     short rdo;
00168     int adc; // pedestal subtracted 
00169     float et;
00170     ushort *thr=mGeom->etow.thr_rdo;
00171     ushort *ped=mGeom->etow.ped_rdo;
00172     float *gain2ET=mGeom->etow.gain2ET_rdo;
00173     float *gain2Ene=mGeom->etow.gain2Ene_rdo;
00174     HitTower1 *hit=etowCalibData.hit;
00175     for(rdo=0; rdo<EtowGeom::mxRdo; rdo++){
00176       if(rawAdc[rdo]<thr[rdo])continue;
00177       if(nTower>=L2EtowCalibData08::mxListSize) break; // overflow protection
00178       adc=rawAdc[rdo]-ped[rdo];  //did NOT correct for common pedestal noise - bad for the jet finder
00179       et=adc/gain2ET[rdo]; 
00180       hit->rdo=rdo;
00181       hit->adc=adc;
00182       hit->et=et;
00183       hit->ene=adc/gain2Ene[rdo]; 
00184       hit++;
00185       nTower++; 
00186       // only monitoring
00187       if(et >par_hotEtThres) {
00188         hA[10]->fill(rdo);
00189         nHotTower++;
00190       }
00191     }
00192     etowCalibData.hitSize=nTower;
00193     
00194     // QA histos
00195     hA[13]->fill(nTower);
00196     hA[14]->fill(nHotTower);
00197     if(nTower>=L2EtowCalibData08::mxListSize) mhN->fill(5); // was overflow 
00198   } // EVEVEVEVEVE
00199 
00200   // debugging should be off for any time critical computation
00201   if(par_dbg>0){
00202     printf("L2-%s-compute: set adcL size=%d\n",getName(),nTower); 
00203     printf("dbg=%s: found  nTw=%d\n",getName(),nTower);
00204     if(par_dbg>0)   print0();
00205     printCalibratedData(token);
00206   } 
00207   
00208   computeStop(token);
00209 
00210 } 
00211 
00212 /* ========================================
00213   ======================================== */
00214 void 
00215 L2etowCalAlgo08::computeUser(int token ){
00216 
00217   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());
00218   assert(1==2);
00219 }
00220 
00221 
00222 /* ========================================
00223   ======================================== */
00224 void
00225 L2etowCalAlgo08::finishRunUser() {  
00226   /* called once at the end of the run
00227      write do whatever you want, log-file & histo-file are still open
00228      Here it seraches for hot tower, re-project histos vs. other representations
00229   */
00230   
00231   int eHotSum=1,eHotId=-1;
00232   const int *data20=hA[10]->getData();
00233   const L2EmcDb::EmcCDbItem *xE=0; //mDb->getByIndex(502); // some wired default?
00234   
00235   int i;
00236   for(i=0; i<EmcDbIndexMax; i++) {
00237     const L2EmcDb::EmcCDbItem *x=mDb->getByIndex(i);
00238     if(mDb->isEmpty(x)) continue;
00239     if (!mDb->isETOW(x) ) continue;
00240     int ieta= (x->eta-1);
00241     int iphi= (x->sec-1)*EtowGeom::mxSubs + x->sub-'A' ;
00242     int softId= iphi+EtowGeom::mxPhiBin*ieta;
00243     hA[11]->fillW(softId,data20[x->rdo]);
00244     hA[12]->fillW(ieta, iphi,data20[x->rdo]);
00245     if(eHotSum<data20[x->rdo]) {
00246       eHotSum=data20[x->rdo];
00247       eHotId=softId;
00248       xE=x;
00249     }
00250   }
00251   
00252   int par_nHotThresh=20;
00253   if (mLogFile && eHotSum>par_nHotThresh){
00254     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);
00255   }
00256   
00257   //...... QA tokens .....
00258   int tkn1=99999, tkn2=0; // min/max token
00259   int nTkn=0;
00260   int tkn3=-1, nTkn3=-1; // most often used token
00261   
00262   int k;
00263   for(k=0;k<L2eventStream2008::mxToken;k++){
00264     L2EtowCalibData08 & etowCalibData=globL2eventStream2008.etow[k];  
00265     if(etowCalibData.nInputBlock==0) continue;
00266     hA[1]->fillW(k,etowCalibData.nInputBlock);
00267     if(nTkn3<etowCalibData.nInputBlock){
00268       nTkn3=etowCalibData.nInputBlock;
00269       tkn3=k;
00270     }
00271 
00272     nTkn++;
00273     if(tkn1>k) tkn1=k;
00274     if(tkn2<k) tkn2=k;
00275   }
00276   if (mLogFile){
00277     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);
00278   }
00279 
00280 }
00281 
00282 
00283 //=======================================
00284 //=======================================
00285 void 
00286 L2etowCalAlgo08::createHisto() {
00287   memset(hA,0,sizeof(hA));
00288   //token related spectra
00289   hA[1]=new  L2Histo(1,"L2-etow-calib: seen tokens;  x:  token value; y: events ",L2eventStream2008::mxToken);
00290   
00291   // ETOW  raw spectra (zz 4 lines)
00292   hA[10]=new L2Histo(10,"etow hot tower 1", EtowGeom::mxRdo); // title upadted in initRun
00293   hA[11]=new L2Histo(11,"etow hot tower 2", EtowGeom::mxRdo); // title upadted in initRun
00294   hA[12]=new L2Histo(12,"etow hot tower 3", EtowGeom::mxEtaBin,EtowGeom::mxPhiBin); // title upadted in initRun  
00295   hA[13]=new L2Histo(13,"ETOW #tower w/ energy /event; x: # ETOW towers; y: counts", 30); 
00296   hA[14]=new L2Histo(14,"# hot towers/event", 30); 
00297 
00298 }
00299 
00300 
00301 
00302 /* ========================================
00303   ======================================== */
00304 void 
00305 L2etowCalAlgo08::print0(){ // full raw input  ADC array
00306   // empty
00307  }
00308 
00309 
00310 /**********************************************************************
00311   $Log: L2etowCalAlgo08.cxx,v $
00312   Revision 1.3  2008/02/01 00:16:40  balewski
00313   add mxListSize to BTOW/ETOW calibration
00314 
00315   Revision 1.2  2008/01/30 21:56:40  balewski
00316   E+B high-enery-filter L2-algo fuly functional
00317 
00318   Revision 1.1  2008/01/30 00:47:16  balewski
00319   Added L2-Etow-calib
00320 
00321  
00322  
00323 */
00324 
00325