Using detector quality flag (taken from ALICE logbook) to decide whether to rpodcue...
[u/mrichter/AliRoot.git] / ZDC / AliZDCTowerCalibTask.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, 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
16 #include "TChain.h"
17 #include "TTree.h"
18 #include "TH1F.h"
19 #include "TCanvas.h"
20 #include "TDecompLU.h" 
21
22 #include "AliAnalysisTask.h"
23 #include "AliAnalysisManager.h"
24
25 #include "AliESDEvent.h"
26 #include "AliESDInputHandler.h"
27 #include "AliESDZDC.h"
28 #include "AliZDCTowerCalibTask.h"
29 #include "AliZDCTowerCalib.h"
30 #include "AliCDBManager.h"
31 #include "AliCDBId.h"
32 #include "AliCDBRunRange.h"
33 #include "AliCDBStorage.h"
34
35 ClassImp(AliZDCTowerCalibTask)
36
37 //________________________________________________________________________
38 AliZDCTowerCalibTask::AliZDCTowerCalibTask() 
39   : AliAnalysisTask(), 
40     fESD(0),
41     fAZNA(),
42     fAZNC(),
43     fBZNA(),
44     fBZNC(),
45     fADCMin(0)
46 {
47 }
48
49 //________________________________________________________________________
50 AliZDCTowerCalibTask::AliZDCTowerCalibTask(const char *name) 
51   : AliAnalysisTask(name, ""), 
52     fESD(0),
53     fAZNA(),
54     fAZNC(),
55     fBZNA(),
56     fBZNC(),
57     fADCMin(0)
58 {
59   // Default constructor
60   
61   // Define input and output slots here
62   // Input slot #0 works with a TChain
63   DefineInput(0, TChain::Class());
64 //   // Output slot #0 writes into a TH1 container
65 //   DefineOutput(0, TH1F::Class());
66   
67   fAZNA.ResizeTo(4,4); 
68   fAZNC.ResizeTo(4,4); 
69   fBZNA.ResizeTo(4); 
70   fBZNC.ResizeTo(4); 
71   
72   fAZNA.Zero(); 
73   fAZNC.Zero(); 
74   fBZNA.Zero(); 
75   fBZNC.Zero(); 
76   fADCMin = 0; 
77 }
78
79 //________________________________________________________________________
80 void AliZDCTowerCalibTask::ConnectInputData(Option_t *) 
81 {
82   // Connect ESD or AOD here
83   // Called once
84   
85   TTree* tree = dynamic_cast<TTree*> (GetInputData(0));
86   if (!tree) Printf("ERROR: Could not read chain from input slot 0");
87   else {
88     // Disable all branches and enable only the needed ones
89     // The next two lines are different when data produced as AliESDEvent is read
90     /*
91       tree->SetBranchStatus("*", kFALSE);
92       tree->SetBranchStatus("fTracks.*", kTRUE);
93     */
94     
95     AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
96     if (!esdH) Printf("ERROR: Could not get ESDInputHandler");
97     else fESD = esdH->GetEvent();
98   }
99 }
100
101 //________________________________________________________________________
102 void AliZDCTowerCalibTask::CreateOutputObjects(){
103   // Create histograms
104   // Called once
105 }
106
107 //________________________________________________________________________
108 void AliZDCTowerCalibTask::Exec(Option_t *) 
109 {
110   // Main loop
111   // Called for each event
112
113   if (!fESD) {
114     Printf("ERROR: fESD not available");
115     return;
116   }
117
118   AliESDZDC *zdcEvent = fESD->GetESDZDC();
119   Double_t *ezn1 = 0, *ezn2 = 0; 
120   if(fESD->GetEventType()==7){
121     ezn1 = (Double_t *) zdcEvent->GetZN1TowerEnergy(); // znc
122     ezn2 = (Double_t *) zdcEvent->GetZN2TowerEnergy(); // zna
123   }
124   Double_t ezn1sum = ezn1[1] + ezn1[2] + ezn1[3] + ezn1[4]; 
125   Double_t ezn2sum = ezn2[1] + ezn2[2] + ezn2[3] + ezn2[4]; 
126   // select event
127   Bool_t isZNAAccepted = kTRUE; 
128   if (ezn2[0] < fADCMin || ezn2sum < fADCMin) isZNAAccepted = kFALSE; 
129   Bool_t isZNCAccepted = kTRUE; 
130   if (ezn1[0] < fADCMin || ezn1sum < fADCMin) isZNCAccepted = kFALSE; 
131   // calculate coefficient matrix and known terms 
132   for (Int_t i=0; i<4; i++) { 
133     for (Int_t j=0; j<4; j++) { 
134       if (isZNAAccepted) fAZNA[i][j] += ezn2[ i + 1 ] * ezn2[ j + 1 ] / ezn2[0]; 
135       if (isZNCAccepted) fAZNC[i][j] += ezn1[ i + 1 ] * ezn1[ j + 1 ] / ezn1[0]; 
136     }
137     if (isZNAAccepted) fBZNA[i] += ezn2[ i + 1 ]; 
138     if (isZNCAccepted) fBZNC[i] += ezn1[ i + 1 ];     
139   }
140   // Post output data.
141   //  PostData(0, fHistPt);
142 }      
143
144 //________________________________________________________________________
145 void AliZDCTowerCalibTask::Terminate(Option_t *) 
146 {
147   // solve the system of linear equations giving the calibration coefficients 
148   // for zna and znc
149
150   TDecompLU luZNA(fAZNA);
151   Bool_t isZNAOk;
152   TVectorD coeffZNA = luZNA.Solve(fBZNA,isZNAOk);
153   if (isZNAOk) { 
154     printf ("ZNA: coefficient matrix:\n"); 
155     fAZNA.Print(); 
156     printf ("ZNA: known terms:\n"); 
157     fBZNA.Print(); 
158     printf ("ZNA: Fitted calibration coefficients:\n"); 
159     coeffZNA.Print(); 
160   }
161   else { 
162     printf ("Singular coefficient matrix for ZNA! Setting all coefficients to 1\n"); 
163     for (Int_t i=0; i<4; i++) coeffZNA[i] = 1; 
164   } 
165
166   TDecompLU luZNC(fAZNC);
167   Bool_t isZNCOk;
168   TVectorD coeffZNC = luZNC.Solve(fBZNC,isZNCOk);
169   if (isZNCOk) { 
170     printf ("ZNC: coefficient matrix:\n"); 
171     fAZNC.Print(); 
172     printf ("ZNC: known terms:\n"); 
173     fBZNC.Print(); 
174     printf ("ZNC: Fitted calibration coefficients:\n"); 
175     coeffZNC.Print(); 
176   }
177   else {
178     printf ("Singular coefficient matrix for ZNC! Setting all coefficients to 1\n"); 
179     for (Int_t i=0; i<4; i++) coeffZNC[i] = 1; 
180   } 
181   // write to OCDB
182
183
184
185   AliZDCTowerCalib *towerCalib = new AliZDCTowerCalib();
186   
187   towerCalib->SetZN1EqualCoeff(0, 1.);
188   towerCalib->SetZN2EqualCoeff(0, 1.);
189   
190   for(Int_t j=1; j<5; j++){  
191     towerCalib->SetZN1EqualCoeff(j, coeffZNC[j+1]);
192     towerCalib->SetZN2EqualCoeff(j, coeffZNA[j+1]);
193   }
194   towerCalib->Print("");
195   
196   AliCDBManager *manager = AliCDBManager::Instance();
197   manager->SetDefaultStorage("local://$ALICE_ROOT/OCDB");
198   
199   AliCDBMetaData *md = new AliCDBMetaData();
200   md->SetResponsible("Chiara Oppedisano");
201   md->SetComment("Calibration object for ZDC written by a macro");
202   md->SetObjectClassName("AliZDCTowerCalib");
203   md->SetAliRootVersion(gSystem->Getenv("ARVERSION"));
204   md->SetBeamPeriod(0);
205   
206   AliCDBId id("ZDC/Calib/TowerCalib",fESD->GetRunNumber(),AliCDBRunRange::Infinity());
207   manager->GetDefaultStorage()->Put(towerCalib,id, md);
208 }