LCOV - code coverage report
Current view: top level - ITS/ITSrec - AliITSOnlineSDDInjectors.cxx (source / functions) Hit Total Coverage
Test: coverage.info Lines: 1 489 0.2 %
Date: 2016-06-14 17:26:59 Functions: 1 35 2.9 %

          Line data    Source code
       1             : /**************************************************************************
       2             :  * Copyright(c) 2007-2009, ALICE Experiment at CERN, All rights reserved. *
       3             :  *                                                                        *
       4             :  * Author: The ALICE Off-line Project.                                    *
       5             :  * Contributors are mentioned in the code where appropriate.              *
       6             :  *                                                                        *
       7             :  * Permission to use, copy, modify and distribute this software and its   *
       8             :  * documentation strictly for non-commercial purposes is hereby granted   *
       9             :  * without fee, provided that the above copyright notice appears in all   *
      10             :  * copies and that both the copyright notice and this permission notice   *
      11             :  * appear in the supporting documentation. The authors make no claims     *
      12             :  * about the suitability of this software for any purpose. It is          *
      13             :  * provided "as is" without express or implied warranty.                  *
      14             :  **************************************************************************/
      15             : #include <TFile.h>
      16             : #include "AliITSOnlineSDDInjectors.h"
      17             : #include "AliLog.h"
      18             : #include <TH2F.h>
      19             : #include <TF1.h>
      20             : #include <TGraphErrors.h>
      21             : #include <TMath.h>
      22             : #include <TString.h>
      23             : 
      24             : /* $Id$ */
      25             : 
      26             : ///////////////////////////////////////////////////////////////////
      27             : //                                                               //
      28             : // Implementation of the class used for SDD injector analysis    //
      29             : // Origin: F.Prino, Torino, prino@to.infn.it                     //
      30             : //                                                               //
      31             : ///////////////////////////////////////////////////////////////////
      32             : 
      33         116 : ClassImp(AliITSOnlineSDDInjectors)
      34             : 
      35             : const Float_t AliITSOnlineSDDInjectors::fgkSaturation = 1008.;
      36             : const Float_t AliITSOnlineSDDInjectors::fgkDefaultLThreshold1 = 8.;
      37             : const Float_t AliITSOnlineSDDInjectors::fgkDefaultLThreshold = 15.;
      38             : const Float_t AliITSOnlineSDDInjectors::fgkDefaultHThreshold1 =15.;
      39             : const Float_t AliITSOnlineSDDInjectors::fgkDefaultHThreshold = 30.;
      40             : const Float_t AliITSOnlineSDDInjectors::fgkDefaultMinSpeed = 5.5;
      41             : const Float_t AliITSOnlineSDDInjectors::fgkDefaultMaxSpeed = 9.0;
      42             : const Float_t AliITSOnlineSDDInjectors::fgkDefaultMaxErr = 1.5;
      43             : const Int_t   AliITSOnlineSDDInjectors::fgkDefaultPolDegree = 3;
      44             : const Float_t AliITSOnlineSDDInjectors::fgkDefaultTimeStep = 50.;
      45             : const UShort_t AliITSOnlineSDDInjectors::fgkDefaultTbMin[kInjLines] = {10,50,100};
      46             : const UShort_t AliITSOnlineSDDInjectors::fgkDefaultTbMax[kInjLines] = {20,70,120};
      47             : 
      48             : //______________________________________________________________________
      49             : AliITSOnlineSDDInjectors::AliITSOnlineSDDInjectors():
      50           0 :   AliITSOnlineSDD(),
      51           0 :   fHisto(),
      52           0 :   fTbZero(0.),
      53           0 :   fRMSTbZero(0.),
      54           0 :   fNEvents(0),
      55           0 :   fParam(),
      56           0 :   fPolDegree(0),
      57           0 :   fActualPolDegree(0),
      58           0 :   fMinDriftSpeed(0.),
      59           0 :   fMaxDriftSpeed(0.),
      60           0 :   fMaxDriftSpeedErr(0.),
      61           0 :   fFirstPadForFit(0),
      62           0 :   fLastPadForFit(0),
      63           0 :   fPadStatusCutForFit(0),
      64           0 :   fTimeStep(0.),
      65           0 :   fUseTimeZeroSignal(kFALSE),
      66           0 :   fMaxCellsAboveThreshold(40)
      67           0 : {
      68             :   // default constructor
      69           0 :   SetPositions();
      70           0 :   SetDefaults();
      71           0 :   SetTimeStep(fgkDefaultTimeStep);
      72           0 :   for(Int_t i=0;i<kInjPads;i++){ 
      73           0 :     fSumDriftSpeed[i]=0.;
      74           0 :     fSumSqDriftSpeed[i]=0.;
      75           0 :     fSumPadStatus[i]=0;
      76           0 :     fSumPadStatusCut[i]=0;
      77           0 :     fNEventsInPad[i]=0;
      78             :   }
      79           0 :   Reset();
      80           0 : }
      81             : //______________________________________________________________________
      82             : AliITSOnlineSDDInjectors::AliITSOnlineSDDInjectors(Int_t nddl, Int_t ncarlos, Int_t sid):
      83           0 :   AliITSOnlineSDD(nddl,ncarlos,sid),
      84           0 :   fHisto(),
      85           0 :   fTbZero(0.),
      86           0 :   fRMSTbZero(0.),
      87           0 :   fNEvents(0),
      88           0 :   fParam(),
      89           0 :   fPolDegree(0),
      90           0 :   fActualPolDegree(0),
      91           0 :   fMinDriftSpeed(0.),
      92           0 :   fMaxDriftSpeed(0.),
      93           0 :   fMaxDriftSpeedErr(0.),
      94           0 :   fFirstPadForFit(0),
      95           0 :   fLastPadForFit(0),
      96           0 :   fPadStatusCutForFit(0),
      97           0 :   fTimeStep(0.),
      98           0 :   fUseTimeZeroSignal(kFALSE),
      99           0 :   fMaxCellsAboveThreshold(40)
     100           0 : { 
     101             : // standard constructor
     102           0 :   SetPositions();
     103           0 :   SetDefaults();
     104           0 :   SetTimeStep(fgkDefaultTimeStep);
     105           0 :   for(Int_t i=0;i<kInjPads;i++){ 
     106           0 :     fSumDriftSpeed[i]=0.;
     107           0 :     fSumSqDriftSpeed[i]=0.;
     108           0 :     fSumPadStatus[i]=0;
     109           0 :     fSumPadStatusCut[i]=0;
     110           0 :     fNEventsInPad[i]=0;
     111             :   }
     112           0 :   Reset();
     113           0 : }
     114             : //______________________________________________________________________
     115           0 : AliITSOnlineSDDInjectors::~AliITSOnlineSDDInjectors(){
     116             :   // Destructor
     117             :   // fHisto should not be deleted here because it points to an histo created 
     118             :   // by the external code which calls the method AnalyzeEvent
     119             :   // if(fHisto) delete fHisto;  
     120           0 :   if(fParam) delete [] fParam;
     121           0 : }
     122             : //______________________________________________________________________
     123             : void AliITSOnlineSDDInjectors::SetDefaults(){
     124             :   // Sets default values for parameters
     125           0 :   for(Int_t i=0;i<kInjLines;i++) {
     126           0 :     SetInjLineRange(i,fgkDefaultTbMin[i],fgkDefaultTbMax[i]);
     127           0 :     SetUseLine(i,kTRUE);
     128           0 :     SetThresholds(i,fgkDefaultLThreshold,fgkDefaultHThreshold);
     129             :   }
     130           0 :   SetThresholds(0,fgkDefaultLThreshold1,fgkDefaultHThreshold1);
     131           0 :   SetPolDegree(fgkDefaultPolDegree);
     132           0 :   SetMinDriftSpeed(fgkDefaultMinSpeed);
     133           0 :   SetMaxDriftSpeed(fgkDefaultMaxSpeed);
     134           0 :   SetMaxDriftSpeedErr(fgkDefaultMaxErr);
     135           0 :   SetFitLimits(1,kInjPads-2); // exclude first and last pad
     136           0 :   SetPadStatusCutForFit();
     137           0 : }
     138             : //______________________________________________________________________
     139             : void AliITSOnlineSDDInjectors::Set20MHzConfig(){
     140             :   // Sets specific parameters for 20 MHz running
     141           0 :   SetInjLineRange(0,10,20);
     142           0 :   SetInjLineRange(1,50,70);
     143           0 :   SetInjLineRange(2,100,120);
     144           0 :   SetTimeStep(50.);
     145           0 :   SetMaxNumberOfCellsPerAnode(40);
     146           0 : }
     147             : //______________________________________________________________________
     148             : void AliITSOnlineSDDInjectors::Set40MHzConfig(){
     149             :   // Sets specific parameters for 20 MHz running
     150           0 :   SetInjLineRange(0,20,50);
     151           0 :   SetInjLineRange(1,90,160);
     152           0 :   SetInjLineRange(2,170,240);
     153           0 :   SetTimeStep(25.);
     154           0 :   SetMaxNumberOfCellsPerAnode(80);
     155           0 : }
     156             : //______________________________________________________________________
     157             : void AliITSOnlineSDDInjectors::SetPositions(){
     158             :   // Sets drift distances for the 3 injector lines
     159           0 :   Double_t xLinFromCenterUm[kInjLines]={31860.,17460.,660.};
     160             :   Double_t xAnodeFromCenterUm=35085;
     161           0 :   for(Int_t i=0;i<kInjLines;i++){
     162           0 :     fPosition[i]=xAnodeFromCenterUm-xLinFromCenterUm[i];
     163           0 :     fPosition[i]/=10000.; // from microns to cm
     164             :   }
     165           0 : }
     166             : //______________________________________________________________________
     167             : void AliITSOnlineSDDInjectors::Reset(){
     168             :   // Resets all counters
     169           0 :   for(Int_t i=0;i<kInjPads;i++){ 
     170           0 :     fDriftSpeed[i]=0.;
     171           0 :     fDriftSpeedErr[i]=0.;
     172             :   }
     173           0 :   for(Int_t i=0;i<kInjPads;i++){
     174           0 :     for(Int_t j=0;j<kInjLines;j++){
     175           0 :       fGoodInj[i][j]=0;
     176           0 :       fCentroid[i][j]=0.;
     177           0 :       fRMSCentroid[i][j]=0.;
     178             :     }
     179             :   }
     180           0 : }
     181             : //______________________________________________________________________
     182             : void AliITSOnlineSDDInjectors::AnalyzeEvent(TH2F* his){
     183             :   // Analyze the current event
     184           0 :   AddEvent(his);
     185           0 :   FitDriftSpeedVsAnode();
     186           0 : }
     187             : //______________________________________________________________________
     188             : void AliITSOnlineSDDInjectors::AddEvent(TH2F* his){
     189             :   // Add the drift speed from current event to the average value
     190           0 :   if(fNEvents==0){
     191           0 :     for(Int_t i=0;i<kInjPads;i++){ 
     192           0 :       fSumDriftSpeed[i]=0.;
     193           0 :       fSumSqDriftSpeed[i]=0.;
     194           0 :       fSumPadStatus[i]=0;
     195           0 :       fSumPadStatusCut[i]=0;
     196           0 :       fNEventsInPad[i]=0;
     197             :     }
     198           0 :   }
     199           0 :   Reset();
     200           0 :   fHisto=his;
     201           0 :   FindGoodInjectors();
     202           0 :   FindCentroids();
     203           0 :   CalcTimeBinZero();
     204           0 :   for(Int_t j=0;j<kInjPads;j++){ 
     205           0 :     CalcDriftSpeed(j);
     206           0 :     Int_t padStatus=GetInjPadStatus(j);
     207           0 :     fSumPadStatus[j]+=padStatus;
     208           0 :     if(padStatus>fPadStatusCutForFit){
     209           0 :       fSumDriftSpeed[j]+=fDriftSpeed[j];
     210           0 :       fSumSqDriftSpeed[j]+=fDriftSpeed[j]*fDriftSpeed[j];
     211           0 :       fSumPadStatusCut[j]+=padStatus;
     212           0 :       fNEventsInPad[j]++;
     213           0 :     }
     214             :   }
     215           0 :   ++fNEvents;
     216           0 : }
     217             : //______________________________________________________________________
     218             : Double_t AliITSOnlineSDDInjectors::GetRMSDriftSpeed(Int_t ipad) const {
     219             :   // Compute RMS of drift speed distribution on one anode
     220           0 :   if(fNEventsInPad[ipad]<=1) return 0.;
     221           0 :   Double_t mean=fSumDriftSpeed[ipad]/(Double_t)fNEventsInPad[ipad];
     222           0 :   Double_t diff=fSumSqDriftSpeed[ipad]/(Double_t)fNEventsInPad[ipad]-mean*mean;
     223           0 :   if(diff<0.) diff=0.;
     224           0 :   return TMath::Sqrt(diff);
     225           0 : }
     226             : 
     227             : //______________________________________________________________________
     228             : void AliITSOnlineSDDInjectors::FitMeanDriftSpeedVsAnode(){
     229             :   // Fits the average drift speed vs.anode number
     230           0 :   if(fNEvents==0) return;
     231           0 :   for(Int_t i=0;i<kInjPads;i++){ 
     232           0 :     fDriftSpeed[i]=GetMeanDriftSpeed(i);
     233           0 :     Int_t padStatusCut=(Int_t)(GetMeanPadStatusCut(i)+0.5);
     234           0 :     for(Int_t ilin=0; ilin<kInjLines ; ilin++) fGoodInj[i][ilin]=(padStatusCut&1<<ilin)>>ilin;
     235           0 :     if(fNEventsInPad[i]>1){
     236           0 :       Double_t rms=GetRMSDriftSpeed(i);
     237           0 :       if(rms>0.) fDriftSpeedErr[i]=rms/TMath::Sqrt(fNEventsInPad[i]);
     238           0 :     }else{
     239           0 :       for(Int_t ilin=0; ilin<kInjLines ; ilin++) fGoodInj[i][ilin]=0;
     240             :     }
     241             :   }
     242           0 :   FitDriftSpeedVsAnode();
     243           0 :   for(Int_t i=0;i<kInjPads;i++){ 
     244           0 :     Int_t padStatus=(Int_t)(GetMeanPadStatusCut(i)+0.5);
     245           0 :     for(Int_t ilin=0; ilin<kInjLines ; ilin++) fGoodInj[i][ilin]=(padStatus&1<<ilin)>>ilin;
     246             :   }
     247           0 : }
     248             : //______________________________________________________________________
     249             : TGraphErrors* AliITSOnlineSDDInjectors::GetTimeVsDistGraph(Int_t jpad) const{
     250             :   // Builds the graph of drift time vs. drift distance
     251             :   const Int_t kPts=kInjLines+1;
     252           0 :   Float_t x[kPts],y[kPts],ex[kPts],ey[kPts];
     253           0 :   x[0]=0.;
     254           0 :   ex[0]=0.;
     255           0 :   y[0]=fTbZero;
     256           0 :   ey[0]=0.;
     257           0 :   for(Int_t i=0;i<kInjLines;i++){
     258           0 :     x[i+1]=fPosition[i];
     259           0 :     ex[i+1]=0.;
     260           0 :     y[i+1]=fCentroid[jpad][i];
     261           0 :     ey[i+1]=fRMSCentroid[jpad][i];
     262             :   }
     263           0 :   TGraphErrors *g=new TGraphErrors(4,x,y,ex,ey);
     264           0 :   return g;
     265           0 : }
     266             : 
     267             : //______________________________________________________________________
     268             : TGraphErrors* AliITSOnlineSDDInjectors::GetDriftSpeedGraph() const{
     269             :   // Builds the graph of drift speed vs. anode number
     270             :   Int_t ipt=0;
     271           0 :   TGraphErrors *g=new TGraphErrors(0);
     272           0 :   for(Int_t i=0;i<kInjPads;i++){
     273           0 :     if(fDriftSpeed[i]>0){ 
     274           0 :       g->SetPoint(ipt,GetAnodeNumber(i),fDriftSpeed[i]);
     275           0 :       g->SetPointError(ipt,0,fDriftSpeedErr[i]);
     276           0 :       ipt++;
     277           0 :     }
     278             :   }
     279           0 :   return g;
     280           0 : }
     281             : //______________________________________________________________________
     282             : TGraphErrors* AliITSOnlineSDDInjectors::GetSelectedDriftSpeedGraph(Int_t minAcceptStatus) const{
     283             :   // TGraphErrors with only pads with status of injector >= minAcceptStatus
     284             :   Int_t ipt=0;
     285           0 :   TGraphErrors *g=new TGraphErrors(0);
     286           0 :   for(Int_t i=0;i<kInjPads;i++){
     287           0 :     Int_t padStatus = GetInjPadStatus(i);
     288           0 :     if(fDriftSpeed[i]>0 && padStatus >= minAcceptStatus ){
     289           0 :       g->SetPoint(ipt,GetAnodeNumber(i),fDriftSpeed[i]);
     290           0 :       g->SetPointError(ipt,0,fDriftSpeedErr[i]);
     291           0 :       ipt++;
     292           0 :     }
     293             :   }
     294           0 :   return g;
     295           0 : }
     296             : //______________________________________________________________________
     297             : void AliITSOnlineSDDInjectors::CalcTimeBinZero(){
     298             :   // Get time zero from trigger signal
     299             :   Double_t tzero=0.,intCont=0.,rmsPeak=0.;
     300           0 :   Bool_t isTbUsed[256];
     301             :   Int_t nTbUsed=0;
     302           0 :   for(Int_t i=0;i<256;i++) isTbUsed[i]=0;
     303           0 :   for(Int_t ian=0;ian<fgkNAnodes;ian++){
     304           0 :     for(Int_t itb=1;itb<fTbMin[0];itb++){
     305           0 :       Double_t cont=fHisto->GetBinContent(itb,ian+1);
     306           0 :       Double_t contm1=fHisto->GetBinContent(itb+1,ian+1);
     307           0 :       Double_t contp1=fHisto->GetBinContent(itb-1,ian+1);
     308           0 :       if(cont>fLowThreshold[0]){
     309           0 :         if(cont>fHighThreshold[0] &&(contm1>fLowThreshold[0] || contp1>fLowThreshold[0])){
     310           0 :           tzero+=cont*float(itb);
     311           0 :           rmsPeak+=cont*float(itb)*float(itb);
     312           0 :           intCont+=cont;
     313           0 :           if(!isTbUsed[itb]){
     314           0 :             isTbUsed[itb]=1;
     315           0 :             ++nTbUsed;
     316           0 :           }
     317             :         }
     318             :       }
     319             :     }
     320             :   }
     321           0 :   if(intCont>0){ 
     322           0 :     fTbZero=tzero/intCont;
     323           0 :     fRMSTbZero=TMath::Sqrt(rmsPeak/intCont-fTbZero*fTbZero);
     324           0 :   }
     325           0 :   if(nTbUsed==1) fRMSTbZero=0.5; 
     326           0 : }
     327             : 
     328             : //______________________________________________________________________
     329             : void AliITSOnlineSDDInjectors::FitDriftSpeedVsAnode(){
     330             :   // fits the anode dependence of drift speed
     331             : 
     332             :   Float_t rangeForMax[2]={78.,178.};
     333           0 :   PolyFit(fPolDegree);
     334           0 :   fActualPolDegree=fPolDegree;
     335           0 :   if(fPolDegree==3){
     336           0 :     Double_t deltasq=fParam[2]*fParam[2]-3*fParam[1]*fParam[3];
     337             :     Double_t zero1=-999.;
     338             :     Double_t zero2=-999.;
     339           0 :     if(deltasq>=0. && TMath::Abs(fParam[3])>0.){
     340           0 :       Double_t delta=TMath::Sqrt(deltasq);
     341           0 :       zero1=(-fParam[2]+delta)/3./fParam[3];
     342           0 :       zero2=(-fParam[2]-delta)/3./fParam[3];
     343           0 :     }
     344             :     Bool_t twoZeroes=kFALSE;
     345             :     Bool_t oneZero=kFALSE;
     346           0 :     if(zero1>0. && zero1<256. && zero2>0. && zero2<256.) twoZeroes=kTRUE;
     347           0 :     if(zero1>rangeForMax[0] && zero1<rangeForMax[1]) oneZero=kTRUE;
     348           0 :     if(zero2>rangeForMax[0] && zero2<rangeForMax[1]) oneZero=kTRUE;
     349           0 :     if(!oneZero || twoZeroes){
     350           0 :       PolyFit(2);
     351             :       Double_t xmax=-999.;
     352           0 :       if(fParam[2]<0.) xmax=-fParam[1]/2./fParam[2];
     353           0 :       if(xmax>rangeForMax[0] && xmax<rangeForMax[1]){
     354           0 :         fActualPolDegree=2;
     355           0 :       }else{
     356             :         Double_t averSpeed=0.;
     357             :         Double_t sumWei=0.;
     358             :         Int_t nUsedPts=0;
     359           0 :         for(Int_t jpad=fFirstPadForFit; jpad<=fLastPadForFit; jpad++){
     360           0 :           if(fDriftSpeed[jpad]>0 && GetInjPadStatus(jpad)>fPadStatusCutForFit){
     361           0 :             Double_t wei=1./fDriftSpeedErr[jpad]/fDriftSpeedErr[jpad];
     362           0 :             averSpeed+=wei*fDriftSpeed[jpad];
     363           0 :             sumWei+=wei;
     364           0 :             nUsedPts++;
     365           0 :           }
     366             :         }
     367           0 :         if(sumWei>0.) averSpeed/=sumWei;
     368           0 :         if(nUsedPts<fPolDegree+1) averSpeed=0;
     369           0 :         fParam[0]=averSpeed;
     370           0 :         for(Int_t i=1; i < fPolDegree+1; i++) fParam[i]=0.;
     371           0 :         fActualPolDegree=0;
     372             :       }
     373           0 :     }
     374           0 :   }
     375           0 : }
     376             : //______________________________________________________________________
     377             : void AliITSOnlineSDDInjectors::PolyFit(Int_t degree){
     378             :   // fits the anode dependence of drift speed with a polynomial function
     379           0 :   const Int_t kNn=degree+1;
     380           0 :   const Int_t kDimens=fPolDegree+1;
     381             : 
     382           0 :   Double_t **mat = new Double_t*[kNn];
     383           0 :   for(Int_t i=0; i < kNn; i++) mat[i] = new Double_t[kNn];
     384           0 :   Double_t *vect = new Double_t[kNn];
     385             : 
     386           0 :   for(Int_t k1=0;k1<kNn;k1++){
     387           0 :     vect[k1]=0;
     388           0 :     for(Int_t k2=0;k2<kNn;k2++){
     389           0 :       mat[k1][k2]=0;
     390             :     }
     391             :   }
     392             :   Int_t npts = 0;
     393           0 :   for(Int_t k1=0;k1<kNn;k1++){
     394           0 :     for(Int_t jpad=fFirstPadForFit; jpad<=fLastPadForFit; jpad++){
     395           0 :       Double_t x=(Double_t)GetAnodeNumber(jpad);
     396           0 :       if(fDriftSpeed[jpad]>0 && GetInjPadStatus(jpad)>fPadStatusCutForFit){
     397           0 :         vect[k1]+=fDriftSpeed[jpad]*TMath::Power(x,k1)/TMath::Power(fDriftSpeedErr[jpad],2);    
     398           0 :         if(k1==0) npts++;
     399           0 :         for(Int_t k2=0;k2<kNn;k2++){
     400           0 :           mat[k1][k2]+=TMath::Power(x,k1+k2)/TMath::Power(fDriftSpeedErr[jpad],2);
     401             :         }
     402           0 :       }
     403             :     }
     404             :   }
     405           0 :   if(npts<fPolDegree+1){ 
     406           0 :     if(fParam) delete [] fParam;
     407           0 :     fParam=new Double_t[kDimens];
     408           0 :     for(Int_t i=0; i<kDimens;i++)fParam[i]=0;
     409           0 :   }else{
     410           0 :     Int_t *iPivot = new Int_t[kNn];
     411           0 :     Int_t *indxR = new Int_t[kNn];
     412           0 :     Int_t *indxC = new Int_t[kNn];
     413           0 :     for(Int_t i=0;i<kNn;i++) iPivot[i]=0;
     414             :     Int_t iCol=-1,iRow=-1;
     415           0 :     for(Int_t i=0;i<kNn;i++){
     416             :       Double_t big=0.;
     417           0 :       for(Int_t j=0;j<kNn;j++){
     418           0 :         if(iPivot[j]!=1){
     419           0 :           for(Int_t k=0;k<kNn;k++){
     420           0 :             if(iPivot[k]==0){
     421           0 :               if(TMath::Abs(mat[j][k])>=big){
     422           0 :                 big=TMath::Abs(mat[j][k]);
     423             :                 iRow=j;
     424             :                 iCol=k;
     425           0 :               }
     426             :             }
     427             :           }
     428           0 :         }
     429             :       }
     430           0 :       iPivot[iCol]++;
     431             :       Double_t aux;
     432           0 :       if(iRow!=iCol){
     433           0 :         for(Int_t l=0;l<kNn;l++){
     434           0 :           aux=mat[iRow][l];
     435           0 :           mat[iRow][l]=mat[iCol][l];
     436           0 :           mat[iCol][l]=aux;
     437             :         }
     438           0 :         aux=vect[iRow];
     439           0 :         vect[iRow]=vect[iCol];
     440           0 :         vect[iCol]=aux;
     441           0 :       }
     442           0 :       indxR[i]=iRow;
     443           0 :       indxC[i]=iCol;
     444           0 :       if(mat[iCol][iCol]==0) break;
     445           0 :       Double_t pivinv=1./mat[iCol][iCol];
     446           0 :       mat[iCol][iCol]=1;
     447           0 :       for(Int_t l=0;l<kNn;l++) mat[iCol][l]*=pivinv;
     448           0 :       vect[iCol]*=pivinv;
     449           0 :       for(Int_t m=0;m<kNn;m++){
     450           0 :         if(m!=iCol){
     451           0 :           aux=mat[m][iCol];
     452           0 :           mat[m][iCol]=0;
     453           0 :           for(Int_t n=0;n<kNn;n++) mat[m][n]-=mat[iCol][n]*aux;
     454           0 :           vect[m]-=vect[iCol]*aux;
     455           0 :         }
     456             :       }    
     457           0 :     }
     458           0 :     delete [] iPivot;
     459           0 :     delete [] indxR;
     460           0 :     delete [] indxC;
     461             :     
     462             :   
     463           0 :     if(fParam) delete [] fParam;
     464           0 :     fParam=new Double_t[kDimens];
     465           0 :     for(Int_t i=0; i<kNn;i++)fParam[i]=vect[i];
     466           0 :     if(degree<fPolDegree) for(Int_t i=kNn; i<kDimens;i++)fParam[i]=0.;
     467             :   }
     468             : 
     469           0 :   for(Int_t i=0; i < kNn; i++) delete [] mat[i];
     470           0 :   delete [] mat;
     471           0 :   delete [] vect;
     472           0 : }
     473             : //______________________________________________________________________
     474             : void AliITSOnlineSDDInjectors::CalcDriftSpeed(Int_t jpad){
     475             :   // Computes the drift speed from the fit to the 3 injector lines for each anode
     476             :   Double_t sumY=0,sumX=0,sumXX=0,sumYY=0.,sumXY=0,sumWEI=0.;
     477             :   Int_t npt=0;
     478           0 :   Double_t y[kInjLines],ey[kInjLines];
     479             :   Double_t tzero=0,erry=0;
     480           0 :   for(Int_t i=0;i<kInjLines;i++){ 
     481           0 :     y[i]=fCentroid[jpad][i];
     482           0 :     ey[i]=fRMSCentroid[jpad][i];
     483             :   }
     484           0 :   for(Int_t i=0;i<kInjLines;i++){
     485           0 :     if(!fUseLine[i]) continue;
     486           0 :     if(fGoodInj[jpad][i] && ey[i]!=0){
     487           0 :       sumY+=y[i]/ey[i]/ey[i];
     488           0 :       sumX+=fPosition[i]/ey[i]/ey[i];
     489           0 :       sumXX+=fPosition[i]*fPosition[i]/ey[i]/ey[i];
     490           0 :       sumYY+=y[i]*y[i]/ey[i]/ey[i];
     491           0 :       sumXY+=fPosition[i]*y[i]/ey[i]/ey[i];
     492           0 :       sumWEI+=1./ey[i]/ey[i];
     493           0 :       tzero=fTbZero/ey[i]/ey[i];
     494           0 :       erry=ey[i]/ey[i]/ey[i];
     495           0 :       npt++;
     496           0 :     }
     497             :   }
     498             :   Double_t slope=0.,eslope=0.;
     499           0 :   if(npt==1){
     500           0 :     slope=(sumY-tzero)/sumX;
     501           0 :     eslope=erry/sumX;
     502           0 :   }
     503           0 :   if(npt>1){ 
     504           0 :     if(fUseTimeZeroSignal){
     505           0 :       sumY+=fTbZero/fRMSTbZero/fRMSTbZero;
     506           0 :       sumX+=0.;
     507           0 :       sumXX+=0.;
     508           0 :       sumYY+=fTbZero*fTbZero/fRMSTbZero/fRMSTbZero;
     509           0 :       sumXY+=0.;
     510           0 :       sumWEI+=1./fRMSTbZero/fRMSTbZero;
     511           0 :     }
     512           0 :     slope=(sumWEI*sumXY-sumY*sumX)/(sumWEI*sumXX-sumX*sumX);
     513           0 :     eslope=TMath::Sqrt(sumWEI/(sumWEI*sumXX-sumX*sumX));
     514           0 :   }
     515             : 
     516             :   Double_t vel=0,evel=0;
     517           0 :   if(slope!=0. && fTimeStep>0.){
     518           0 :     vel=1./slope*10000./fTimeStep;// micron/ns
     519           0 :     evel=eslope/slope/slope*10000./fTimeStep;// micron/ns
     520           0 :   }
     521           0 :   if(vel>fMaxDriftSpeed||vel<fMinDriftSpeed || evel>fMaxDriftSpeedErr){ 
     522             :     vel=0.;
     523             :     evel=0.;
     524           0 :   }
     525           0 :   fDriftSpeed[jpad]=vel;
     526           0 :   fDriftSpeedErr[jpad]=evel;
     527           0 : }
     528             : //______________________________________________________________________
     529             : Int_t AliITSOnlineSDDInjectors::GetAnodeNumber(Int_t iInjPad) const{
     530             :   // Injectors location along anodes:
     531             :   // Side left  (UP)   - channel 0: injectors on anodes 0,7,15,...,247,255 
     532             :   // Side right (DOWN) - channel 1: injectors on anodes 0,8,16,...,248,255 
     533             :   Int_t ian=-1;
     534           0 :   if(iInjPad>=kInjPads) return ian;
     535           0 :   if(fSide==1){  // right side
     536           0 :     ian=iInjPad*8;
     537           0 :     if(iInjPad==32) ian--;
     538             :   }else{         // left side
     539           0 :     ian=iInjPad*8-1;
     540           0 :     if(iInjPad==0) ian=0;
     541             :   }
     542           0 :   return ian;
     543           0 : }
     544             : //______________________________________________________________________
     545             : Int_t AliITSOnlineSDDInjectors::GetInjPadNumberFromAnode(Int_t nAnode) const{
     546             :   // Converts anode number into injector pad index
     547             :   Int_t iInjPad=-1;
     548           0 :   if(fSide==1){  // right side
     549           0 :     if(nAnode%8==0) iInjPad=nAnode/8;
     550           0 :     if(nAnode==255) iInjPad=32;
     551             :   }else{         // left side
     552           0 :     if(nAnode%8==7) iInjPad=1+nAnode/8;
     553           0 :     if(nAnode==0) iInjPad=0;
     554             :   }
     555           0 :   if(nAnode>=256) iInjPad=-1;
     556           0 :   return iInjPad;
     557             : }
     558             : //______________________________________________________________________
     559             : Int_t AliITSOnlineSDDInjectors::GetInjPadStatus(Int_t jpad) const{
     560             :   // returns an integer value with status of injector lines for given pad/anode
     561             :   // status=7  -->  111  all injector are good
     562             :   // status=6  -->  110  1st line (close to anodes) is bad, other two are good
     563             :   // ....
     564             :   // status=1  -->  001  only 1st line (close to anodes) good
     565             :   // status=0  -->  000  all lines are bad
     566             :   Int_t istatus=0;
     567           0 :   if(jpad>=0 && jpad<kInjPads){
     568           0 :     for(Int_t jlin=0;jlin<kInjLines;jlin++) istatus+=fGoodInj[jpad][jlin]<<jlin;
     569           0 :   }
     570           0 :   return istatus;
     571             : }
     572             : //______________________________________________________________________
     573             : void AliITSOnlineSDDInjectors::FindGoodInjectors(){
     574             :   // Mark good injector pads
     575             :   // good = 1 cell above high threshold + 1 neighbour above low threshold
     576           0 :   for(Int_t jpad=0;jpad<kInjPads;jpad++){
     577           0 :     Int_t ian=GetAnodeNumber(jpad);
     578             :     Int_t countAbove=0;
     579           0 :     for(Int_t jjj=0; jjj<fHisto->GetNbinsX(); jjj++){
     580           0 :       Float_t c=fHisto->GetBinContent(jjj+1,ian+1);
     581           0 :       if(c>0.5) countAbove++;
     582             :     }
     583           0 :     if(countAbove>fMaxCellsAboveThreshold){
     584           0 :       for(Int_t jlin=0;jlin<kInjLines;jlin++) fGoodInj[jpad][jlin]=0;           
     585           0 :     }else{
     586           0 :       for(Int_t jlin=0;jlin<kInjLines;jlin++){
     587           0 :         for(Int_t jjj=fTbMin[jlin];jjj<fTbMax[jlin];jjj++){
     588           0 :           Float_t c1=fHisto->GetBinContent(jjj,ian+1);
     589           0 :           Float_t c2=fHisto->GetBinContent(jjj+1,ian+1);
     590             :           //    Float_t c3=fHisto->GetBinContent(jjj+2,ian+1);
     591           0 :           if(c1>fLowThreshold[jlin] && c2>fLowThreshold[jlin]){ 
     592           0 :             if(c1>fHighThreshold[jlin] || c2>fHighThreshold[jlin]){
     593           0 :               fGoodInj[jpad][jlin]=1;
     594           0 :               break;
     595             :             }
     596             :           }
     597           0 :         }
     598             :       }
     599             :     }
     600             :   }
     601           0 : }
     602             : //______________________________________________________________________
     603             : void AliITSOnlineSDDInjectors::FindCentroids(){
     604             :   // Computes the centroids (weighted mean) of teh injector pads
     605           0 :   for(Int_t jpad=0;jpad<kInjPads;jpad++){
     606           0 :     Int_t ian=GetAnodeNumber(jpad);
     607           0 :     for(Int_t jlin=0;jlin<kInjLines;jlin++){
     608           0 :       if(!fGoodInj[jpad][jlin]) continue;
     609             :       Double_t maxcont=0;
     610             :       Int_t ilmax=-1;
     611           0 :       for(Int_t jjj=fTbMin[jlin];jjj<fTbMax[jlin];jjj++){
     612           0 :         Double_t cont=fHisto->GetBinContent(jjj,ian+1);
     613           0 :         if(cont>maxcont){
     614             :           maxcont=cont;
     615             :           ilmax=jjj;
     616           0 :         }
     617             :       }
     618             :       Double_t intCont=0;
     619             :       Int_t jjj=ilmax;
     620           0 :       while(1){
     621           0 :         Double_t cont=fHisto->GetBinContent(jjj,ian+1);
     622           0 :         if(cont<fLowThreshold[jlin]) break;
     623           0 :         if(cont<fgkSaturation){
     624           0 :           fCentroid[jpad][jlin]+=cont*(Double_t)jjj;
     625           0 :           fRMSCentroid[jpad][jlin]+=cont*(Double_t)jjj*(Double_t)jjj;
     626           0 :           intCont+=cont;
     627           0 :         }
     628           0 :         jjj--;
     629           0 :       }
     630           0 :       jjj=ilmax+1;
     631           0 :       while(1){
     632           0 :         Double_t cont=fHisto->GetBinContent(jjj,ian+1);
     633           0 :         if(cont<fLowThreshold[jlin]) break;
     634           0 :         if(cont<fgkSaturation){
     635           0 :           fCentroid[jpad][jlin]+=cont*float(jjj);
     636           0 :           fRMSCentroid[jpad][jlin]+=cont*(Double_t)jjj*(Double_t)jjj;
     637           0 :           intCont+=cont;
     638           0 :         }
     639           0 :         jjj++;
     640           0 :       }
     641           0 :       if(intCont>0){ 
     642           0 :         fCentroid[jpad][jlin]/=intCont;
     643           0 :         fRMSCentroid[jpad][jlin]=TMath::Sqrt(fRMSCentroid[jpad][jlin]/intCont-fCentroid[jpad][jlin]*fCentroid[jpad][jlin])/TMath::Sqrt(intCont);
     644           0 :       }
     645             :       else{ 
     646           0 :         fCentroid[jpad][jlin]=0.;
     647           0 :         fRMSCentroid[jpad][jlin]=0.;
     648           0 :         fGoodInj[jpad][jlin]=0;
     649             :       }
     650           0 :       if(fRMSCentroid[jpad][jlin]==0) fGoodInj[jpad][jlin]=0;
     651           0 :     }
     652             :   }
     653           0 : }
     654             : //______________________________________________________________________
     655             : void AliITSOnlineSDDInjectors::PrintInjectorStatus(){
     656             :   // Dumps the status bit of injector pads
     657           0 :   for(Int_t jpad=0;jpad<kInjPads;jpad++){
     658           0 :     printf("Line%d-Anode%d: %d %d %d\n",jpad,GetAnodeNumber(jpad),fGoodInj[jpad][0],fGoodInj[jpad][1],fGoodInj[jpad][2]);
     659             :   }
     660           0 : }
     661             : //______________________________________________________________________
     662             : void AliITSOnlineSDDInjectors::PrintCentroids(){
     663             :   // Dumps the drift time centroids of injector pads
     664           0 :   for(Int_t jpad=0;jpad<kInjPads;jpad++){
     665           0 :     printf("Line%d-Anode%d: %f+-%f %f+-%f %f+-%f\n",jpad,GetAnodeNumber(jpad),fCentroid[jpad][0],fRMSCentroid[jpad][0],fCentroid[jpad][1],fRMSCentroid[jpad][1],fCentroid[jpad][2],fRMSCentroid[jpad][2]);
     666             :   }
     667           0 : }
     668             : //______________________________________________________________________
     669             : void AliITSOnlineSDDInjectors::WriteToASCII(Int_t evNumb, UInt_t timeStamp, Int_t optAppend){
     670             :   // writes drift speed vs. anode fit parameters into an ASCII file 
     671             :   // to be sent to FXS by the DA and processed by the SHUTTLE
     672             : 
     673           0 :   TString outfilnam;
     674           0 :   outfilnam.Form("SDDinj_ddl%02dc%02d_sid%d.data",fDDL,fCarlos,fSide);  
     675             :   FILE* outf;
     676           0 :   if(optAppend==0){ 
     677           0 :     outf=fopen(outfilnam.Data(),"w");
     678           0 :     fprintf(outf,"%d\n",fActualPolDegree);
     679             :   }
     680           0 :   else outf=fopen(outfilnam.Data(),"a");
     681           0 :   fprintf(outf,"%d   %d   ",evNumb,timeStamp);
     682           0 :   for(Int_t ic=0;ic<fPolDegree+1;ic++){
     683           0 :     fprintf(outf,"%G ",fParam[ic]);
     684             :   }
     685           0 :   fprintf(outf,"\n");
     686           0 :   fclose(outf);  
     687           0 : }
     688             : //______________________________________________________________________
     689             : TH1F* AliITSOnlineSDDInjectors::GetMeanDriftSpeedVsPadHisto() const{
     690             :   // Builds histogram of average drift speed vs. pad number
     691           0 :   TString hisnam;
     692           0 :   hisnam.Form("hdrsp%02dc%02ds%d",fDDL,fCarlos,fSide);
     693           0 :   TH1F* h=new TH1F(hisnam.Data(),"",kInjPads,-0.5,kInjPads-0.5);
     694           0 :   if(fNEvents>0){
     695           0 :     for(Int_t i=0;i<kInjPads;i++){ 
     696           0 :       h->SetBinContent(i+1,GetMeanDriftSpeed(i));    
     697           0 :       Double_t rms=GetRMSDriftSpeed(i);
     698             :       Double_t err=0.;
     699           0 :       if(rms>0.) err=rms/TMath::Sqrt(fNEventsInPad[i]);
     700           0 :       h->SetBinError(i+1,err);
     701             :     }
     702           0 :   }
     703             :   return h;
     704           0 : }
     705             : //______________________________________________________________________
     706             : Bool_t AliITSOnlineSDDInjectors::WriteToROOT(TFile *fil) const {
     707             :   // Writes the output histograms into a root file
     708           0 :   if(fil==0){ 
     709           0 :     AliWarning("Invalid pointer to ROOT file");
     710           0 :     return kFALSE;    
     711             :   }  
     712           0 :   TString hisnam;
     713           0 :   fil->cd();
     714           0 :   hisnam.Form("hdrsp%02dc%02ds%d",fDDL,fCarlos,fSide);
     715           0 :   TH1F hdsp(hisnam.Data(),"",kInjPads,-0.5,kInjPads-0.5);
     716           0 :   if(fNEvents==0){
     717           0 :     AliWarning("Zero analyzed events");
     718           0 :     return kFALSE;    
     719             :   }  
     720             :     
     721           0 :   for(Int_t i=0;i<kInjPads;i++){ 
     722           0 :     hdsp.SetBinContent(i+1,GetMeanDriftSpeed(i));    
     723           0 :     Double_t rms=GetRMSDriftSpeed(i);
     724             :     Double_t err=0.;
     725           0 :     if(rms>0.) err=rms/TMath::Sqrt(fNEventsInPad[i]);
     726           0 :     hdsp.SetBinError(i+1,err);
     727             :   }
     728           0 :   hdsp.Write();
     729           0 :   return kTRUE;    
     730           0 : }
     731             : //______________________________________________________________________
     732             : void AliITSOnlineSDDInjectors::WriteInjectorStatusToASCII(){
     733             :   // dump status of injectors encoded into UInt_t
     734             :   // 5 bits (value 0-31) to store number of pads with given status
     735           0 :   TString outfilnam;
     736           0 :   outfilnam.Form("SDDinj_ddl%02dc%02d_sid%d.data",fDDL,fCarlos,fSide);  
     737           0 :   FILE* outf=fopen(outfilnam.Data(),"a");
     738           0 :   Int_t n[8]={0,0,0,0,0,0,0,0};
     739           0 :   for(Int_t jpad=fFirstPadForFit; jpad<=fLastPadForFit; jpad++){
     740           0 :     Int_t statusPad=GetInjPadStatus(jpad);
     741           0 :     ++n[statusPad];
     742             :   }
     743             :   UInt_t statusInj=0;
     744           0 :   statusInj+=(n[7]&0x1F)<<25; // bits 25-29: n. of pads with status 7
     745           0 :   statusInj+=(n[6]&0x1F)<<20; // bits 20-24: n. of pads with status 6
     746           0 :   statusInj+=(n[5]&0x1F)<<15; // bits 15-19: n. of pads with status 5
     747           0 :   statusInj+=(n[4]&0x1F)<<10; // bits 10-14: n. of pads with status 4
     748           0 :   statusInj+=(n[3]&0x1F)<<5;  // bits  5- 9: n. of pads with status 3
     749           0 :   statusInj+=(n[2]&0x1F);     // bits  0- 4: n. of pads with status 2
     750             : 
     751           0 :   fprintf(outf,"-99 %u\n",statusInj); // -99 used in preprocessor to find line
     752             :                                       // with injector status info
     753           0 :   fclose(outf);  
     754             :   
     755           0 : }

Generated by: LCOV version 1.11