changes in the MagF constructor
[u/mrichter/AliRoot.git] / MUON / AliMUONChamberCalibrationTask.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 // $Id$
17
18 //-----------------------------------------------------------------------------
19 /// \file   AliMUONChamberCalibrationTask.cxx
20 /// \brief  Implementation of the AliMUONChamberCalibrationTask 
21 /// \author Andry Rakotozafindrabe CEA/IRFU/SPhN
22 //-----------------------------------------------------------------------------
23
24 #include <Riostream.h>
25
26 #include <TBranch.h>
27 #include <TChain.h>
28 #include <TClonesArray.h>
29 #include <TFile.h>
30 #include <TMath.h>
31 #include <TRandom.h>
32 #include <TROOT.h>
33 #include <TString.h>
34 #include <TTree.h>
35
36 #include "AliMUONChamberCalibrationTask.h"
37
38 // STEER includes
39 #include "AliCDBManager.h"
40 #include "AliESDEvent.h"
41 #include "AliESDInputHandler.h"
42 #include "AliESDtrack.h"
43 #include "AliESDMuonTrack.h"
44 #include "AliLog.h"
45 #include "AliMagF.h"
46 #include "AliRecoParam.h"
47 #include "AliTracker.h"
48
49 // ANALYSIS includes
50 #include "AliAnalysisDataSlot.h"
51 #include "AliAnalysisManager.h"
52
53 // MUON includes
54 #include "AliMpConstants.h"
55 #include "AliMpCDB.h"
56 #include "AliMpPad.h"
57 #include "AliMpSegmentation.h"
58 #include "AliMpVSegmentation.h"
59 #include "AliMUONCalibrationData.h"
60 #include "AliMUONClusterInfo.h"
61 #include "AliMUONESDInterface.h"
62 #include "AliMUONPadInfo.h"
63 #include "AliMUONRecoParam.h"
64 #include "AliMUONTrack.h"
65 #include "AliMUONTrackParam.h"
66 #include "AliMUONVCalibParam.h"
67 #include "AliMUONVCluster.h"
68 //#include "AliMUONVClusterStore.h"
69 #include "AliMUONVDigit.h"
70 #include "AliMUONVDigitStore.h"
71
72 /// \cond CLASSIMP
73 ClassImp( AliMUONChamberCalibrationTask )
74 /// \endcond CLASSIMP
75
76 //______________________________________________________________
77 AliMUONChamberCalibrationTask::AliMUONChamberCalibrationTask():
78   AliAnalysisTaskSE( "AliMUONChamberCalibrationTask" ),
79   fOCDBPath( "local://$ALICE_ROOT/OCDB" ),
80   fCalibChoice(NOGAIN),
81   fClusterInfoTree(0x0),
82   fMuonRecoParam(0x0),
83   fClusterInfo(0x0),
84   fCalibData(0x0),
85   fESDInterface(0x0),
86   fDigitStore(0x0),
87   fESDInputHandler(0x0),
88   fESDInputEvent(0x0)
89 {
90   //
91   /// Default constructor
92   //
93
94 }
95
96 //______________________________________________________________
97 AliMUONChamberCalibrationTask::AliMUONChamberCalibrationTask( const char* name,
98                                                               char* ocdbpath,
99                                                               const Int_t my_calib_option ):
100   AliAnalysisTaskSE( name ),
101   fOCDBPath( "local://$ALICE_ROOT/OCDB" ),
102   fCalibChoice(NOGAIN),
103   fClusterInfoTree(0x0),
104   fMuonRecoParam(0x0),
105   fClusterInfo(0x0),
106   fCalibData(0x0),
107   fESDInterface(0x0),
108   fDigitStore(0x0),
109   fESDInputHandler(0x0),
110   fESDInputEvent(0x0)
111 {
112   //
113   /// constructor
114   //
115
116   fOCDBPath = ocdbpath;
117   if ( (my_calib_option >= ((Int_t)NOGAIN)) && (my_calib_option <= ((Int_t)INJECTIONGAIN)) ) 
118     fCalibChoice = (Calibration_t)my_calib_option;
119   else {
120     AliWarning( Form("Wrong value of the calibration option %d not within [%d, %d] !!! Will use NOGAIN", 
121                      my_calib_option, (Int_t)NOGAIN, (Int_t)INJECTIONGAIN ) );
122     fCalibChoice = NOGAIN;
123   }
124 }
125
126 //______________________________________________________________
127 AliMUONChamberCalibrationTask::~AliMUONChamberCalibrationTask()
128 {
129   //
130   /// destructor
131   //
132
133   delete fMuonRecoParam;
134   delete fClusterInfo;
135   delete fESDInterface;
136
137 }
138 //______________________________________________________________
139 void AliMUONChamberCalibrationTask::CreateOutputObjects()
140 {
141   //
142   /// Creates the output TTree
143   //
144
145   AliDebug( 1, "" ); 
146
147   TFile* clusterInfoFile = OpenFile( 0, "RECREATE" );
148   if( clusterInfoFile ) clusterInfoFile->SetCompressionLevel(1);
149   else AliError( "no output file created !!!" );
150
151   if ( !fClusterInfoTree ) fClusterInfoTree = new TTree( "clusterInfoTree", "clusterInfoTree" ); 
152   fClusterInfoTree->Branch( "clusterInfo" , &fClusterInfo, 32000, 99);
153 }
154
155 //______________________________________________________________
156 void AliMUONChamberCalibrationTask::LocalInit()
157 {
158   //
159   /// Initialization
160   /// Initialize the cluster info and the ESD interface
161   /// Set the magnetic field, the mapping and the reconstruction parameters
162   //
163
164   AliDebug( 1, "" );
165
166   // initialize the cluster info and the ESD interface
167
168   fClusterInfo = new AliMUONClusterInfo();
169   fESDInterface = new AliMUONESDInterface();
170
171   gRandom->SetSeed(0);
172   
173   // set mag field
174
175   if ( !TGeoGlobalMagField::Instance()->GetField() ) {
176     AliInfo( "Loading field map..." );
177     AliMagF* field = new AliMagF( "Maps","Maps", 1., 1., AliMagF::k5kG );
178     TGeoGlobalMagField::Instance()->SetField( field );
179     TGeoGlobalMagField::Instance()->Lock();
180   }
181   
182   // Load mapping
183
184   AliCDBManager* man = AliCDBManager::Instance();
185   man->SetDefaultStorage( fOCDBPath );
186   man->SetSpecificStorage( "MUON/Calib/MappingData", fOCDBPath );
187   man->SetSpecificStorage( "MUON/Calib/MappingRunData", fOCDBPath ); // for the manu serial numbers
188   man->SetRun(0);
189   man->Print();
190   if ( ! AliMpCDB::LoadDDLStore() ) {
191     AliFatal( "Could not access mapping from OCDB !" );
192     exit(-1); 
193   }
194
195   // Set the reconstruction parameters for track refitting 
196   // (needed when applying any of the with-gain options)
197
198   fMuonRecoParam = AliMUONRecoParam::GetCosmicParam();
199
200   TString caliboption1 = "NOGAIN";
201   TString caliboption2 = "GAINCONSTANTCAPA";
202   TString caliboption3 = "GAIN";
203   TString caliboption4 = "INJECTIONGAIN";
204
205   TString caliboption = caliboption1;
206   if ( fCalibChoice == GAINCONSTANTCAPA ) caliboption = caliboption2;
207   if ( fCalibChoice == GAIN ) caliboption = caliboption3;  
208   if ( fCalibChoice == INJECTIONGAIN ) caliboption = caliboption4;
209   fMuonRecoParam->SetCalibrationMode(caliboption.Data());
210
211   for (Int_t iCh=0; iCh<10; iCh++) {
212     fMuonRecoParam->SetDefaultNonBendingReso( iCh, 0.152 ); // input ESD was aligned (default cosmic settings)
213     fMuonRecoParam->SetDefaultBendingReso( iCh, 0.027 );
214   }
215   fMuonRecoParam->SetMaxNonBendingDistanceToTrack(5.); // was at 1. in default cosmic settings
216   fMuonRecoParam->SetMaxBendingDistanceToTrack(5.);
217
218   fMuonRecoParam->RequestStation(1, kTRUE); // only St 4 and 5 enabled in default cosmic settings
219   fMuonRecoParam->ImproveTracks(kTRUE, 7.); // was 6. in default cosmic settings
220
221   AliInfo( "reconstruction parameters initialized as follows :" );
222   fMuonRecoParam->Print("FULL");
223
224   AliMUONESDInterface::ResetTracker(fMuonRecoParam);
225 }
226
227 //______________________________________________________________
228 void AliMUONChamberCalibrationTask::ConnectInputData( Option_t* /*option*/ )
229 {
230   //
231   /// Connect to ESD here
232   /// (called once)
233   //
234
235   AliDebug( 1, "" );
236
237   fESDInputHandler = dynamic_cast<AliESDInputHandler*>(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
238   TTree* tree = NULL;
239
240   if ( fESDInputHandler ) {
241
242     // The properly initialized ESD input handler reads ESD tree 
243     // and connect it to ESD event, so we only need to retrieve the later
244     fESDInputEvent = fESDInputHandler->GetEvent();
245     if ( !fESDInputEvent ) {
246
247       AliFatal( "Could not get input ESD event !!! ");
248
249     } 
250   } else {
251
252       AliError( "Could not get input ESD handler !!!" );
253       // If no input event handler we need to get the tree once
254       // from input slot 0 for the chain
255       tree = dynamic_cast<TTree*> (GetInputData(0));
256       if ( tree ) tree->GetReadEntry();
257       else AliFatal( "Could not read tree from input slot 0 !!!" );
258     }
259 }
260
261 //______________________________________________________________
262 void AliMUONChamberCalibrationTask::Exec( Option_t* /*option*/ )
263 {
264   //
265   /// Process the current event
266   /// (called for each event)
267   //
268
269   static Bool_t first = kTRUE;
270
271   if ( first ) { 
272     AliDebug( 1, "" );
273     first = kFALSE;
274   }
275
276   if ( !fESDInputEvent ) {
277     AliError( "Input ESD event not available !!! " );
278     return;
279   }
280
281   // load the current ESD event
282   fESDInterface->LoadEvent( *fESDInputEvent );
283
284   // get digit store
285   fDigitStore = fESDInterface->GetDigits();
286
287   // prepare access to calibration data 
288   if ( !fCalibData ) fCalibData = new AliMUONCalibrationData( fESDInputEvent->GetESDRun()->GetRunNumber() );
289
290   // --------------------------------------------------------------------
291   // fill cluster info from clusters attached to each track of this event
292   // --------------------------------------------------------------------
293
294   Int_t nTracks = (Int_t)fESDInputEvent->GetNumberOfMuonTracks(); 
295   if ( nTracks < 1 ) return; 
296
297   TIter nextTrack( fESDInterface->CreateTrackIterator() );
298   AliMUONTrack* track;
299
300   while ( (track = static_cast<AliMUONTrack*>(nextTrack())) ) { // loop over tracks
301
302     UInt_t muonClusterMap = BuildClusterMap( *track );
303     
304     AliMUONTrackParam* trackParam = 
305       static_cast<AliMUONTrackParam*>(track->GetTrackParamAtCluster()->First());
306
307     while ( trackParam ) { // loop over clusters
308
309       fClusterInfo->Clear("C");
310       
311       // fill cluster info
312       AliMUONVCluster* cluster = trackParam->GetClusterPtr();
313       fClusterInfo->SetRunId( fESDInputEvent->GetRunNumber() );
314       fClusterInfo->SetEventId( fESDInputEvent->GetEventNumberInFile() );
315       fClusterInfo->SetZ( cluster->GetZ() );
316       fClusterInfo->SetClusterId( cluster->GetUniqueID() );
317       fClusterInfo->SetClusterXY( cluster->GetX(), cluster->GetY() );
318       fClusterInfo->SetClusterXYErr( cluster->GetErrX(), cluster->GetErrY() );
319       fClusterInfo->SetClusterChi2( cluster->GetChi2() );
320       fClusterInfo->SetClusterCharge( cluster->GetCharge() );
321       
322       // fill track info
323       fClusterInfo->SetTrackId( track->GetUniqueID() );
324       fClusterInfo->SetTrackXY( trackParam->GetNonBendingCoor(), trackParam->GetBendingCoor() );
325       fClusterInfo->SetTrackThetaXY( TMath::ATan( trackParam->GetNonBendingSlope() ), 
326                                      TMath::ATan( trackParam->GetBendingSlope() ) );
327       fClusterInfo->SetTrackP( trackParam->P() );
328       const TMatrixD paramCov = trackParam->GetCovariances();
329       fClusterInfo->SetTrackXYErr( TMath::Sqrt( paramCov(0,0) ), 
330                                    TMath::Sqrt( paramCov(2,2) ) );
331       fClusterInfo->SetTrackChi2( track->GetNormalizedChi2() );
332       fClusterInfo->SetTrackCharge( (Short_t)trackParam->GetCharge() );
333       fClusterInfo->SetTrackNHits( track->GetNClusters() );
334       fClusterInfo->SetTrackChamberHitMap( muonClusterMap );
335       
336       // fill pad info if available       
337       for ( Int_t i=0; i<cluster->GetNDigits(); i++ ) {
338
339         AliMUONVDigit* digit = fDigitStore->FindObject( cluster->GetDigitId(i) );
340         if ( !digit ) continue;
341         
342         // pad location
343         const AliMpVSegmentation* seg = AliMpSegmentation::Instance()->
344           GetMpSegmentation( digit->DetElemId(), AliMp::GetCathodType( digit->Cathode() ) );
345         AliMpPad pad = seg->PadByIndices( digit->PadX(), digit->PadY() );
346         
347         // calibration parameters
348         AliMUONVCalibParam* ped = fCalibData->Pedestals( digit->DetElemId(), digit->ManuId() );
349         AliMUONVCalibParam* gain = fCalibData->Gains( digit->DetElemId(), digit->ManuId() );
350         Int_t manuChannel = digit->ManuChannel();
351         Int_t planeType = 0;
352         if ( digit->ManuId() & AliMpConstants::ManuMask(AliMp::kNonBendingPlane) ) {
353           planeType = 1;
354         }
355         
356         // fill pad info
357         AliMUONPadInfo padInfo;
358         padInfo.SetPadId( digit->GetUniqueID() );
359         padInfo.SetPadPlaneType( planeType );
360         padInfo.SetPadXY( pad.GetPositionX(), pad.GetPositionY() );
361         padInfo.SetPadDimXY( pad.GetDimensionX(), pad.GetDimensionY() );
362         padInfo.SetPadCharge( (Double_t)digit->Charge() );
363         padInfo.SetPadADC( digit->ADC() );
364         padInfo.SetSaturated( digit->IsSaturated() );
365         padInfo.SetCalibrated( digit->IsCalibrated() );
366         padInfo.SetPedestal( ped->ValueAsFloatFast(manuChannel,0), // mean
367                              ped->ValueAsFloatFast(manuChannel,1) ); // sigma
368         padInfo.SetGain( gain->ValueAsFloatFast(manuChannel,0), // a0
369                          gain->ValueAsFloatFast(manuChannel,1), // a1
370                          (Int_t)gain->ValueAsFloatFast(manuChannel,2), // threshold
371                          (Int_t)gain->ValueAsFloatFast(manuChannel,3) ); // fit quality
372         
373         fClusterInfo->AddPad( padInfo );
374       }
375             
376       // fill cluster info tree
377       fClusterInfoTree->Fill();
378       
379       trackParam = static_cast<AliMUONTrackParam*>(track->GetTrackParamAtCluster()->After(trackParam));
380     }
381     
382   }
383     // Added protection in case the derived task is not an AOD producer.
384     AliAnalysisDataSlot *out0 = GetOutputSlot(0);
385     if (out0 && out0->IsConnected()) PostData( 0, fClusterInfoTree );    
386
387 }
388
389
390 //______________________________________________________________
391 UInt_t AliMUONChamberCalibrationTask::BuildClusterMap( AliMUONTrack &track )
392 {
393   //
394   /// Build the map of clusters in tracking chambers
395   //
396
397   UInt_t muonClusterMap = 0;
398   
399   AliMUONTrackParam* trackParam = static_cast<AliMUONTrackParam*>(track.GetTrackParamAtCluster()->First());
400   while ( trackParam ) {
401     
402     muonClusterMap |= BIT(trackParam->GetClusterPtr()->GetChamberId());
403     
404     trackParam = static_cast<AliMUONTrackParam*>(track.GetTrackParamAtCluster()->After(trackParam));
405   }
406   
407   return muonClusterMap;
408 }
409
410 //______________________________________________________________
411 void AliMUONChamberCalibrationTask::Terminate( Option_t* /*option*/ )
412 {
413   //
414   /// Called once per task on the client machine at the end of the analysis.
415   //
416
417   AliDebug( 1, "" );
418 }