Task for tower inter-calibration (A. De Falco)

diff --git a/ZDC/AliZDCTowerCalibTask.cxx b/ZDC/AliZDCTowerCalibTask.cxx
new file mode 100644 (file)
index 0000000..bcb2a7e
--- /dev/null
+++ b/ZDC/AliZDCTowerCalibTask.cxx
@@ -0,0 +1,208 @@
+/**************************************************************************
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ *                                                                        *
+ * Author: The ALICE Off-line Project.                                    *
+ * Contributors are mentioned in the code where appropriate.              *
+ *                                                                        *
+ * Permission to use, copy, modify and distribute this software and its   *
+ * documentation strictly for non-commercial purposes is hereby granted   *
+ * without fee, provided that the above copyright notice appears in all   *
+ * copies and that both the copyright notice and this permission notice   *
+ * appear in the supporting documentation. The authors make no claims     *
+ * about the suitability of this software for any purpose. It is          *
+ * provided "as is" without express or implied warranty.                  *
+ **************************************************************************/
+
+#include "TChain.h"
+#include "TTree.h"
+#include "TH1F.h"
+#include "TCanvas.h"
+#include "TDecompLU.h" 
+
+#include "AliAnalysisTask.h"
+#include "AliAnalysisManager.h"
+
+#include "AliESDEvent.h"
+#include "AliESDInputHandler.h"
+#include "AliESDZDC.h"
+#include "AliZDCTowerCalibTask.h"
+#include "AliZDCTowerCalib.h"
+#include "AliCDBManager.h"
+#include "AliCDBId.h"
+#include "AliCDBRunRange.h"
+#include "AliCDBStorage.h"
+
+ClassImp(AliZDCTowerCalibTask)
+
+//________________________________________________________________________
+AliZDCTowerCalibTask::AliZDCTowerCalibTask() 
+  : AliAnalysisTask(), 
+    fESD(0),
+    fAZNA(),
+    fAZNC(),
+    fBZNA(),
+    fBZNC(),
+    fADCMin(0)
+{
+}
+
+//________________________________________________________________________
+AliZDCTowerCalibTask::AliZDCTowerCalibTask(const char *name) 
+  : AliAnalysisTask(name, ""), 
+    fESD(0),
+    fAZNA(),
+    fAZNC(),
+    fBZNA(),
+    fBZNC(),
+    fADCMin(0)
+{
+  // Default constructor
+  
+  // Define input and output slots here
+  // Input slot #0 works with a TChain
+  DefineInput(0, TChain::Class());
+//   // Output slot #0 writes into a TH1 container
+//   DefineOutput(0, TH1F::Class());
+  
+  fAZNA.ResizeTo(4,4); 
+  fAZNC.ResizeTo(4,4); 
+  fBZNA.ResizeTo(4); 
+  fBZNC.ResizeTo(4); 
+  
+  fAZNA.Zero(); 
+  fAZNC.Zero(); 
+  fBZNA.Zero(); 
+  fBZNC.Zero(); 
+  fADCMin = 0; 
+}
+
+//________________________________________________________________________
+void AliZDCTowerCalibTask::ConnectInputData(Option_t *) 
+{
+  // Connect ESD or AOD here
+  // Called once
+  
+  TTree* tree = dynamic_cast<TTree*> (GetInputData(0));
+  if (!tree) Printf("ERROR: Could not read chain from input slot 0");
+  else {
+    // Disable all branches and enable only the needed ones
+    // The next two lines are different when data produced as AliESDEvent is read
+    /*
+      tree->SetBranchStatus("*", kFALSE);
+      tree->SetBranchStatus("fTracks.*", kTRUE);
+    */
+    
+    AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
+    if (!esdH) Printf("ERROR: Could not get ESDInputHandler");
+    else fESD = esdH->GetEvent();
+  }
+}
+
+//________________________________________________________________________
+void AliZDCTowerCalibTask::CreateOutputObjects(){
+  // Create histograms
+  // Called once
+}
+
+//________________________________________________________________________
+void AliZDCTowerCalibTask::Exec(Option_t *) 
+{
+  // Main loop
+  // Called for each event
+
+  if (!fESD) {
+    Printf("ERROR: fESD not available");
+    return;
+  }
+
+  AliESDZDC *zdcEvent = fESD->GetESDZDC();
+  Double_t *ezn1 = 0, *ezn2 = 0; 
+  if(fESD->GetEventType()==7){
+    ezn1 = (Double_t *) zdcEvent->GetZN1TowerEnergy(); // znc
+    ezn2 = (Double_t *) zdcEvent->GetZN2TowerEnergy(); // zna
+  }
+  Double_t ezn1sum = ezn1[1] + ezn1[2] + ezn1[3] + ezn1[4]; 
+  Double_t ezn2sum = ezn2[1] + ezn2[2] + ezn2[3] + ezn2[4]; 
+  // select event
+  Bool_t isZNAAccepted = kTRUE; 
+  if (ezn2[0] < fADCMin || ezn2sum < fADCMin) isZNAAccepted = kFALSE; 
+  Bool_t isZNCAccepted = kTRUE; 
+  if (ezn1[0] < fADCMin || ezn1sum < fADCMin) isZNCAccepted = kFALSE; 
+  // calculate coefficient matrix and known terms 
+  for (Int_t i=0; i<4; i++) { 
+    for (Int_t j=0; j<4; j++) { 
+      if (isZNAAccepted) fAZNA[i][j] += ezn2[ i + 1 ] * ezn2[ j + 1 ] / ezn2[0]; 
+      if (isZNCAccepted) fAZNC[i][j] += ezn1[ i + 1 ] * ezn1[ j + 1 ] / ezn1[0]; 
+    }
+    if (isZNAAccepted) fBZNA[i] += ezn2[ i + 1 ]; 
+    if (isZNCAccepted) fBZNC[i] += ezn1[ i + 1 ];     
+  }
+  // Post output data.
+  //  PostData(0, fHistPt);
+}      
+
+//________________________________________________________________________
+void AliZDCTowerCalibTask::Terminate(Option_t *) 
+{
+  // solve the system of linear equations giving the calibration coefficients 
+  // for zna and znc
+
+  TDecompLU luZNA(fAZNA);
+  Bool_t isZNAOk;
+  TVectorD coeffZNA = luZNA.Solve(fBZNA,isZNAOk);
+  if (isZNAOk) { 
+    printf ("ZNA: coefficient matrix:\n"); 
+    fAZNA.Print(); 
+    printf ("ZNA: known terms:\n"); 
+    fBZNA.Print(); 
+    printf ("ZNA: Fitted calibration coefficients:\n"); 
+    coeffZNA.Print(); 
+  }
+  else { 
+    printf ("Singular coefficient matrix for ZNA! Setting all coefficients to 1\n"); 
+    for (Int_t i=0; i<4; i++) coeffZNA[i] = 1; 
+  } 
+
+  TDecompLU luZNC(fAZNC);
+  Bool_t isZNCOk;
+  TVectorD coeffZNC = luZNC.Solve(fBZNC,isZNCOk);
+  if (isZNCOk) { 
+    printf ("ZNC: coefficient matrix:\n"); 
+    fAZNC.Print(); 
+    printf ("ZNC: known terms:\n"); 
+    fBZNC.Print(); 
+    printf ("ZNC: Fitted calibration coefficients:\n"); 
+    coeffZNC.Print(); 
+  }
+  else {
+    printf ("Singular coefficient matrix for ZNC! Setting all coefficients to 1\n"); 
+    for (Int_t i=0; i<4; i++) coeffZNC[i] = 1; 
+  } 
+  // write to OCDB
+
+
+
+  AliZDCTowerCalib *towerCalib = new AliZDCTowerCalib();
+  
+  towerCalib->SetZN1EqualCoeff(0, 1.);
+  towerCalib->SetZN2EqualCoeff(0, 1.);
+  
+  for(Int_t j=1; j<5; j++){  
+    towerCalib->SetZN1EqualCoeff(j, coeffZNC[j+1]);
+    towerCalib->SetZN2EqualCoeff(j, coeffZNA[j+1]);
+  }
+  towerCalib->Print("");
+  
+  AliCDBManager *manager = AliCDBManager::Instance();
+  manager->SetDefaultStorage("local://$ALICE_ROOT/OCDB");
+  
+  AliCDBMetaData *md = new AliCDBMetaData();
+  md->SetResponsible("Chiara Oppedisano");
+  md->SetComment("Calibration object for ZDC written by a macro");
+  md->SetObjectClassName("AliZDCTowerCalib");
+  md->SetAliRootVersion(gSystem->Getenv("ARVERSION"));
+  md->SetBeamPeriod(0);
+  
+  AliCDBId id("ZDC/Calib/TowerCalib",fESD->GetRunNumber(),AliCDBRunRange::Infinity());
+  manager->GetDefaultStorage()->Put(towerCalib,id, md);
+}

diff --git a/ZDC/AliZDCTowerCalibTask.h b/ZDC/AliZDCTowerCalibTask.h
new file mode 100644 (file)
index 0000000..cc90b7c
--- /dev/null
+++ b/ZDC/AliZDCTowerCalibTask.h
@@ -0,0 +1,41 @@
+#ifndef AliZDCTowerCalibTask_h
+#define AliZDCTowerCalibTask_h
+
+// analysis task for ZN towers intercalibration
+// Author: Alessandro De Falco, INFN Cagliari
+
+#include "TMatrixD.h"
+#include "TVectorD.h"
+
+class TH1F;
+class AliESDEvent;
+
+#include "AliAnalysisTask.h"
+
+class AliZDCTowerCalibTask : public AliAnalysisTask {
+ public:
+  AliZDCTowerCalibTask();
+  AliZDCTowerCalibTask(const char *name);
+  virtual ~AliZDCTowerCalibTask() {}
+  
+  virtual void   ConnectInputData(Option_t *);
+  virtual void   CreateOutputObjects();
+  virtual void   Exec(Option_t *option);
+  virtual void   Terminate(Option_t *);
+  virtual void   SetADCMin(Double_t adcmin) { fADCMin = adcmin; } 
+
+  
+ private:
+  AliESDEvent *fESD;   //ESD object
+  TMatrixD fAZNA;      // coefficient matrix for ZNA calorimeter
+  TMatrixD fAZNC;      // coefficient matrix for ZNC calorimeter
+  TVectorD fBZNA;      // vector of known terms for ZNA calorimeter
+  TVectorD fBZNC;      // vector of known terms for ZNC calorimeter
+  Double_t fADCMin; 
+  AliZDCTowerCalibTask(const AliZDCTowerCalibTask&); // not implemented
+  AliZDCTowerCalibTask& operator=(const AliZDCTowerCalibTask&); // not implemented
+  
+  ClassDef(AliZDCTowerCalibTask, 1); 
+};
+
+#endif