TRD nSigma OADB related new codes and modifications and OADB root file -- Xianguo Lu
authorxlu <xianguo.lu@cern.ch>
Thu, 31 Jul 2014 21:02:45 +0000 (23:02 +0200)
committermorsch <andreas.morsch@cern.ch>
Tue, 30 Sep 2014 13:17:39 +0000 (15:17 +0200)
OADB/COMMON/PID/data/TRDdEdxParams.root [new file with mode: 0644]
STEER/CMakelibSTEERBase.pkg
STEER/STEERBase/AliPIDResponse.cxx
STEER/STEERBase/AliPIDResponse.h
STEER/STEERBase/AliTRDPIDResponse.cxx
STEER/STEERBase/AliTRDPIDResponse.h
STEER/STEERBase/AliTRDdEdxParams.cxx [new file with mode: 0644]
STEER/STEERBase/AliTRDdEdxParams.h [new file with mode: 0644]
STEER/STEERBaseLinkDef.h
TRD/makeTRDdEdxOADB.C [new file with mode: 0644]

diff --git a/OADB/COMMON/PID/data/TRDdEdxParams.root b/OADB/COMMON/PID/data/TRDdEdxParams.root
new file mode 100644 (file)
index 0000000..474bd8b
Binary files /dev/null and b/OADB/COMMON/PID/data/TRDdEdxParams.root differ
index 3731a82..ea697eb 100644 (file)
@@ -65,6 +65,7 @@ set ( SRCS
     STEERBase/AliTRDPIDResponseObject.cxx
     STEERBase/AliTRDNDFast.cxx
     STEERBase/AliTRDTKDInterpolator.cxx
+    STEERBase/AliTRDdEdxParams.cxx
     STEERBase/AliPIDResponse.cxx 
     STEERBase/AliITSPIDResponse.cxx 
     STEERBase/AliTPCPIDResponse.cxx 
index f590a3f..3c45a64 100644 (file)
@@ -42,6 +42,7 @@
 #include <AliPID.h>
 #include <AliOADBContainer.h>
 #include <AliTRDPIDResponseObject.h>
+#include <AliTRDdEdxParams.h>
 #include <AliTOFPIDParams.h>
 #include <AliHMPIDPIDParams.h>
 
@@ -89,6 +90,7 @@ fOADBvoltageMaps(NULL),
 fUseTPCEtaCorrection(kFALSE),
 fUseTPCMultiplicityCorrection(kFALSE),
 fTRDPIDResponseObject(NULL),
+fTRDdEdxParams(NULL),
 fTOFtail(0.9),
 fTOFPIDParams(NULL),
 fHMPIDPIDParams(NULL),
@@ -114,6 +116,7 @@ AliPIDResponse::~AliPIDResponse()
   //
   delete fArrPidResponseMaster;
   delete fTRDPIDResponseObject;
+  delete fTRDdEdxParams;
   delete fTOFPIDParams;
 }
 
@@ -153,6 +156,7 @@ fOADBvoltageMaps(NULL),
 fUseTPCEtaCorrection(other.fUseTPCEtaCorrection),
 fUseTPCMultiplicityCorrection(other.fUseTPCMultiplicityCorrection),
 fTRDPIDResponseObject(NULL),
+fTRDdEdxParams(NULL),
 fTOFtail(0.9),
 fTOFPIDParams(NULL),
 fHMPIDPIDParams(NULL),
@@ -209,6 +213,7 @@ AliPIDResponse& AliPIDResponse::operator=(const AliPIDResponse &other)
     fUseTPCEtaCorrection=other.fUseTPCEtaCorrection;
     fUseTPCMultiplicityCorrection=other.fUseTPCMultiplicityCorrection;
     fTRDPIDResponseObject=NULL;
+    fTRDdEdxParams=NULL;
     fEMCALPIDParams=NULL;
     fTOFtail=0.9;
     fTOFPIDParams=NULL;
@@ -623,6 +628,8 @@ void AliPIDResponse::ExecNewRun()
   SetTPCEtaMaps();
 
   SetTRDPidResponseMaster(); 
+  //has to precede InitializeTRDResponse(), otherwise the read-out fTRDdEdxParams is not pased in TRDResponse
+  SetTRDdEdxParams();
   InitializeTRDResponse();
 
   SetEMCALPidResponseMaster(); 
@@ -1565,6 +1572,7 @@ void AliPIDResponse::InitializeTRDResponse(){
   // Set PID Params and references to the TRD PID response
   // 
     fTRDResponse.SetPIDResponseObject(fTRDPIDResponseObject);
+    fTRDResponse.SetdEdxParams(fTRDdEdxParams);
 }
 
 //______________________________________________________________________________
@@ -1588,6 +1596,32 @@ void AliPIDResponse::SetTRDSlices(UInt_t TRDslicesForPID[2],AliTRDPIDResponse::E
     AliDebug(1,Form("Slice Range set to %d - %d",TRDslicesForPID[0],TRDslicesForPID[1]));
 }
 
+//______________________________________________________________________________ 
+void AliPIDResponse::SetTRDdEdxParams()
+{
+  if(fTRDdEdxParams) return;
+
+  const TString containerName = "TRDdEdxParamsContainer";
+  AliOADBContainer cont(containerName.Data()); 
+  
+  const TString filePathNamePackage=Form("%s/COMMON/PID/data/TRDdEdxParams.root", fOADBPath.Data());
+
+  const Int_t statusCont = cont.InitFromFile(filePathNamePackage.Data(), cont.GetName());
+  if (statusCont){
+    AliFatal("Failed initializing settings from OADB"); 
+  }
+  else{
+    AliInfo(Form("AliPIDResponse::SetTRDdEdxParams loading %s from %s\n", cont.GetName(), filePathNamePackage.Data()));
+
+    fTRDdEdxParams = (AliTRDdEdxParams*)(cont.GetObject(fRun, "default"));
+    //fTRDdEdxParams->Print();
+
+    if(!fTRDdEdxParams){
+      AliError(Form("TRD dEdx Params default not found"));
+    }
+  }
+}
+
 //______________________________________________________________________________
 void AliPIDResponse::SetTOFPidResponseMaster()
 {
index 280c68e..106ee67 100644 (file)
@@ -30,6 +30,7 @@ class TLinearFitter;
 
 class AliVEvent;
 class AliTRDPIDResponseObject;
+class AliTRDdEdxParams;
 class AliTOFPIDParams;
 class AliHMPIDPIDParams;
 class AliOADBContainer;
@@ -237,6 +238,7 @@ private:
   Bool_t fUseTPCMultiplicityCorrection; // Use TPC multiplicity correction
 
   AliTRDPIDResponseObject *fTRDPIDResponseObject; //! TRD PID Response Object
+  AliTRDdEdxParams * fTRDdEdxParams; //! TRD dEdx Response for truncated mean signal
 
   Float_t fTOFtail;                    //! TOF tail effect used in TOF probability
   AliTOFPIDParams *fTOFPIDParams;      //! TOF PID Params - period depending (OADB loaded)
@@ -275,6 +277,7 @@ private:
   void SetTRDPidResponseMaster();
   void InitializeTRDResponse();
   void SetTRDSlices(UInt_t TRDslicesForPID[2],AliTRDPIDResponse::ETRDPIDMethod method) const;
+  void SetTRDdEdxParams();
 
   //TOF
   void SetTOFPidResponseMaster();
index 180339b..221ebe0 100644 (file)
@@ -42,6 +42,7 @@
 #include "AliTRDPIDResponse.h"
 //#include "AliTRDTKDInterpolator.h"
 #include "AliTRDNDFast.h"
+#include "AliTRDdEdxParams.h"
 
 ClassImp(AliTRDPIDResponse)
 
@@ -49,6 +50,7 @@ ClassImp(AliTRDPIDResponse)
 AliTRDPIDResponse::AliTRDPIDResponse():
   TObject()
   ,fkPIDResponseObject(NULL)
+  ,fkTRDdEdxParams(NULL)
   ,fGainNormalisationFactor(1.)
 {
   //
@@ -60,6 +62,7 @@ AliTRDPIDResponse::AliTRDPIDResponse():
 AliTRDPIDResponse::AliTRDPIDResponse(const AliTRDPIDResponse &ref):
   TObject(ref)
   ,fkPIDResponseObject(NULL)
+  ,fkTRDdEdxParams(NULL)
   ,fGainNormalisationFactor(ref.fGainNormalisationFactor)
 {
   //
@@ -78,7 +81,8 @@ AliTRDPIDResponse &AliTRDPIDResponse::operator=(const AliTRDPIDResponse &ref){
   TObject::operator=(ref);
   fGainNormalisationFactor = ref.fGainNormalisationFactor;
   fkPIDResponseObject = ref.fkPIDResponseObject;
-  
+  fkTRDdEdxParams = ref.fkTRDdEdxParams;
+
   return *this;
 }
 
@@ -87,7 +91,10 @@ AliTRDPIDResponse::~AliTRDPIDResponse(){
   //
   // Destructor
   //
-    if(IsOwner()) delete fkPIDResponseObject;
+  if(IsOwner()) {
+    delete fkPIDResponseObject;
+    delete fkTRDdEdxParams;
+  }
 }
 
 //____________________________________________________________
@@ -148,42 +155,13 @@ Double_t AliTRDPIDResponse::GetSignalDelta( const AliVTrack* track, AliPID::EPar
   //output other information in info[]
   //
 
-  //============================================================================================>>>>>>>>>>>>
-  //Data production dependent, in the end should go to OADB
-  //============================================================================================
-  //~/.task/TRDdEdx/2013output/resolution$ r bgmean.root
-  // hfit;1  1.521e+00 , 7.294e-01 , 8.066e+00 , 2.146e-02 , 3.137e+01 , 7.160e-15 , 2.838e+00 , 1.100e+01 , (14290.45/67) conv_runlist_LHC13b.pass3.root trdncls/trdnch>=18 && trdnch>=6 && tpcncls>120 && pidV0>=1          
-  // hscn0_wid;1     -4.122e-01 , 1.019e+00 , -1.319e-01 , (15071.70/120) conv_runlist_LHC13b.pass3.root tpcncls>120 && pidV0>=1
-
-  //mean signal parametrization
-  Double_t meanpars[]={1.521e+00 , 7.294e-01 , 8.066e+00 , 2.146e-02 , 3.137e+01 , 7.160e-15 , 2.838e+00 , 1.100e+01};
-  const Int_t nmp = sizeof(meanpars)/sizeof(Double_t);
-
-  if(type==AliPID::kProton){
-    const Double_t tmpproton[]={2.026e+00 , -1.462e-04 , 1.202e+00 , 1.162e-01 , 2.092e+00 , -3.018e-02 , 3.665e+00 , 3.487e-07};
-    for(Int_t ii=0; ii<nmp; ii++){
-      meanpars[ii] = tmpproton[ii];
-    }
-  }
-  else if(type==AliPID::kPion || type==AliPID::kElectron){
-    const Double_t tmppe[]={1.508e+00 , 7.356e-01 , 8.002e+00 , 1.932e-02 , 3.058e+01 , 5.114e-18 , 3.561e+00 , 1.408e+01};
-    for(Int_t ii=0; ii<nmp; ii++){
-      meanpars[ii] = tmppe[ii];
-    }
-  }
-
-  //resolution parametrization
-  const Double_t respars[]={-4.122e-01 , 1.019e+00 , -1.319e-01};
+  const Double_t badval = -9999;
 
   //cut on number of chambers
   const Double_t cutch = 6;
 
   //cut on mean number of clusters per chamber
   const Double_t cutclsperch = 18;
-  //============================================================================================<<<<<<<<<<<<<
-
-
-  const Double_t badval = -9999;
 
   const Double_t nch = track->GetTRDNchamberdEdx();
   if(nch < cutch){
@@ -205,12 +183,32 @@ Double_t AliTRDPIDResponse::GetSignalDelta( const AliVTrack* track, AliPID::EPar
     return badval;
   }
 
+  if(!fkTRDdEdxParams){
+    AliError("AliTRDPIDResponse::GetSignalDelta fkTRDdEdxParams null");
+    return -99999;
+  }
+
+  TVectorF meanvec = fkTRDdEdxParams->GetMeanParameter(type);
+  TVectorF resvec  = fkTRDdEdxParams->GetSigmaParameter(type);
+
+  Double_t meanpar[meanvec.GetNrows()];
+  Double_t respar[resvec.GetNrows()];
+
+  for(Int_t ii=0; ii<meanvec.GetNrows(); ii++){
+    meanpar[ii]=meanvec[ii];
+  }
+  for(Int_t ii=0; ii<resvec.GetNrows(); ii++){
+    respar[ii]=resvec[ii];
+  }
+
+  //============================================================================================<<<<<<<<<<<<<
+
   const Double_t bg = pTRD/AliPID::ParticleMass(type);
-  const Double_t expsig = MeandEdxTR(&bg, meanpars);
+  const Double_t expsig = MeandEdxTR(&bg, meanpar);
 
   if(info){
     info[0]= expsig;
-    info[1]= ResolutiondEdxTR(&ncls, respars);
+    info[1]= ResolutiondEdxTR(&ncls, respar);
   }
 
   const Double_t eps = 1e-10;
@@ -553,3 +551,11 @@ Bool_t AliTRDPIDResponse::SetPIDResponseObject(const AliTRDPIDResponseObject * o
     if((AliLog::GetDebugLevel("",IsA()->GetName()))>0)fkPIDResponseObject->Print("");
     return kTRUE;
 }
+
+
+//____________________________________________________________    
+Bool_t AliTRDPIDResponse::SetdEdxParams(const AliTRDdEdxParams * par)
+{
+  fkTRDdEdxParams = par;
+  return kTRUE;
+}
index d8c156f..d0dd4e0 100644 (file)
@@ -31,6 +31,7 @@
 class TObjArray;
 class AliVTrack;
 class AliTRDPIDResponseObject;
+class AliTRDdEdxParams;
 
 class AliTRDPIDResponse : public TObject {
   public:
@@ -75,6 +76,7 @@ class AliTRDPIDResponse : public TObject {
     void      SetGainNormalisationFactor(Double_t gainFactor) { fGainNormalisationFactor = gainFactor; }
 
     Bool_t SetPIDResponseObject(const AliTRDPIDResponseObject * obj);
+    Bool_t SetdEdxParams(const AliTRDdEdxParams * par);
     
     Bool_t    Load(const Char_t *filename = NULL);
   
@@ -85,8 +87,8 @@ class AliTRDPIDResponse : public TObject {
     Double_t  GetProbabilitySingleLayer(Int_t species, Double_t plocal, Double_t *dEdx,ETRDPIDMethod PIDmethod=kLQ1D) const;
     
     const AliTRDPIDResponseObject *fkPIDResponseObject;   // PID References and Params
+    const AliTRDdEdxParams * fkTRDdEdxParams; //parametrisation for truncated mean
     Double_t  fGainNormalisationFactor;         // Gain normalisation factor
-      
   
   ClassDef(AliTRDPIDResponse, 3)    // Tool for TRD PID
 };
diff --git a/STEER/STEERBase/AliTRDdEdxParams.cxx b/STEER/STEERBase/AliTRDdEdxParams.cxx
new file mode 100644 (file)
index 0000000..4e15df9
--- /dev/null
@@ -0,0 +1,88 @@
+/**************************************************************************
+ * 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.                  *
+ **************************************************************************/
+//
+//  Xianguo Lu <lu@physi.uni-heidelberg.de>
+//
+
+#include "AliLog.h"
+#include "AliTRDdEdxParams.h"
+
+ClassImp(AliTRDdEdxParams);
+
+AliTRDdEdxParams::AliTRDdEdxParams(const TString name, const TString title): TNamed(name,title)
+{
+  //
+  //constructor
+  //
+}
+
+void AliTRDdEdxParams::CheckType(const Int_t itype, const TString tag) const 
+{
+  //
+  //check if itype is in range
+  //
+
+  if(itype<0 || itype>=MAXNPAR){
+    AliError(Form("AliTRDdEdxParams::CheckType %s itype out of range %d\n", tag.Data(), itype));
+  }
+}
+
+const TVectorF AliTRDdEdxParams::GetParameter(const TVectorF par[], const Int_t itype)const
+{
+  //
+  //return parameter for particle itype from par[]
+  //
+
+  CheckType(itype, "GetParameter");
+
+  return par[itype];
+}
+
+void AliTRDdEdxParams::SetParameter(TVectorF par[], const Int_t itype, const Int_t npar, const Float_t vals[])
+{
+  //
+  //set parameter, vals of dimension npar, for particle itype
+  //
+
+  CheckType(itype, "SetParameter");
+
+  TVectorF p2(npar, vals);
+
+  par[itype].ResizeTo(p2);
+  par[itype] = p2;
+}
+
+void AliTRDdEdxParams::Print(Option_t* option) const
+{
+  //
+  //print all members
+  //
+
+  TObject::Print(option);
+
+  printf("\n======================= Mean ========================\n");
+  for(Int_t ii=0; ii<MAXNPAR; ii++){
+    printf("%d: Nrows() %d\n",ii, fMeanPar[ii].GetNrows());
+    if(fMeanPar[ii].GetNrows()) fMeanPar[ii].Print();
+  }
+
+  printf("\n======================= Sigma ========================\n");
+
+  for(Int_t ii=0; ii<MAXNPAR; ii++){
+    printf("%d: Nrows() %d\n",ii, fSigmaPar[ii].GetNrows());
+    if(fSigmaPar[ii].GetNrows()) fSigmaPar[ii].Print();
+  }
+  printf("AliTRDdEdxParams::Print done.\n\n");
+}
diff --git a/STEER/STEERBase/AliTRDdEdxParams.h b/STEER/STEERBase/AliTRDdEdxParams.h
new file mode 100644 (file)
index 0000000..7fc9ff8
--- /dev/null
@@ -0,0 +1,53 @@
+/**************************************************************************
+ * 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.                  *
+ **************************************************************************/
+//
+//  Xianguo Lu <lu@physi.uni-heidelberg.de>
+//
+
+#ifndef ALITRDDEDXPARAMS_H
+#define ALITRDDEDXPARAMS_H
+
+#include "TNamed.h"
+#include "TVectorT.h"
+#include "TString.h"
+
+//maximum number of paticle types
+#define MAXNPAR 10
+
+class AliTRDdEdxParams: public TNamed
+{
+ public:
+  AliTRDdEdxParams(const TString name="name", const TString title="title");
+  void Print(Option_t* option = "") const;
+
+  const TVectorF GetMeanParameter(const Int_t itype) const { return GetParameter(fMeanPar, itype);}
+  const TVectorF GetSigmaParameter(const Int_t itype) const { return GetParameter(fSigmaPar, itype);}
+
+  void SetMeanParameter(const Int_t itype, const Int_t npar, const Float_t vals[]){ SetParameter(fMeanPar, itype, npar, vals); }
+  void SetSigmaParameter(const Int_t itype, const Int_t npar, const Float_t vals[]){ SetParameter(fSigmaPar, itype, npar, vals); } 
+
+ private:
+  const TVectorF GetParameter(const TVectorF par[], const Int_t itype) const;
+  void SetParameter(TVectorF par[], const Int_t itype, const Int_t npar, const Float_t vals[]);
+
+  TVectorF fMeanPar[MAXNPAR];
+  TVectorF fSigmaPar[MAXNPAR];
+
+  void CheckType(const Int_t itype, const TString tag) const;
+
+  ClassDef(AliTRDdEdxParams,1);
+};
+
+#endif
index 34c917e..a8a0d21 100644 (file)
@@ -92,6 +92,7 @@
 #pragma link C++ class AliTRDPIDParams::AliTRDPIDThresholds+;
 #pragma link C++ class AliTRDPIDParams::AliTRDPIDCentrality+;
 /* #endif */
+#pragma link C++ class AliTRDdEdxParams+;
 #pragma link C++ class AliTRDPIDResponseObject+;
 #pragma link C++ class AliTRDNDFast+;
 #pragma link C++ class AliTRDTKDInterpolator+;
diff --git a/TRD/makeTRDdEdxOADB.C b/TRD/makeTRDdEdxOADB.C
new file mode 100644 (file)
index 0000000..532590c
--- /dev/null
@@ -0,0 +1,253 @@
+/*
+g++ makeTRDdEdxOADB.C  -g -O3 -Wall -Werror -I$ROOTSYS/include -I$ALICE_BUILD/include -I$ALICE_ROOT/TRD -I$ALICE_ROOT/TRD/Cal  -L$ROOTSYS/lib -lCore -lCint -lRIO -lNet -lHist -lGraf -lGraf3d -lGpad -lTree -lRint -lPostscript -lMatrix -lPhysics -lMathCore -lThread -pthread -lm -ldl -rdynamic -lMinuit -lEG -lGeom -lVMC -lProof -lProofPlayer -lXMLParser -lXMLIO -lSpectrum -lTreePlayer -lMLP -lGui -L$ALICE_BUILD/lib/tgt_linuxx8664gcc -lCDB -lSTEER -lRAWDatarec -lESD -lSTEERBase -lANALYSIS -lRAWDatabase  -lANALYSISalice -lAOD -lTPCrec -lTPCbase -lTRDbase -lTRDrec -lSTAT   -o makeTRDdEdxOADB
+
+#generate OADB
+makeTRDdEdxOADB 1
+
+#read and print OADB content
+makeTRDdEdxOADB 0
+ */
+
+#include <math.h>
+#include <stdio.h>
+#include <fstream>
+#include <string>
+#include <iostream>
+using namespace std;
+
+#include "Math/Functor.h"
+#include "Math/Factory.h"
+#include "Math/Minimizer.h"
+
+#include "TASImage.h"
+#include "TAxis.h"
+#include "TColor.h"
+#include "TCut.h"
+#include "TCanvas.h"
+#include "TChain.h"
+#include "TDatabasePDG.h"
+#include "TDecompLU.h"
+#include "TDecompSVD.h"
+#include "TDirectory.h"
+#include "TEventList.h"
+#include "TF1.h"
+#include "TF2.h"
+#include "TFile.h"
+#include "TGaxis.h"
+#include "TGeoManager.h"
+#include "TGeoGlobalMagField.h"
+#include "TRandom3.h"
+#include "TGraph.h"
+#include "TGraphAsymmErrors.h"
+#include "TGraphErrors.h"
+#include "TGraphPolar.h"
+#include "TGrid.h"
+#include "TH1F.h"
+#include "TH2F.h"
+#include "TH1D.h"
+#include "TH2D.h"
+#include "TH3D.h"
+#include "THnSparse.h"
+#include "TLatex.h"
+#include "TLegend.h"
+#include "TLegendEntry.h"
+#include "TLinearFitter.h"
+#include "TMarker.h"
+#include "TMath.h"
+#include "TMatrixD.h"
+#include "TMinuit.h"
+#include "TPaletteAxis.h"
+#include "TPaveText.h"
+#include "TPolyMarker.h"
+#include "TProfile.h"
+#include "TROOT.h"
+#include "TStopwatch.h"
+#include "TString.h"
+#include "TStyle.h"
+#include "TSystem.h"
+#include "TSystemDirectory.h"
+#include "TTree.h"
+#include "TTimeStamp.h"
+#include "TUUID.h"
+#include "TVector3.h"
+#include "TVectorD.h"
+
+#include "AliAnalysisManager.h"
+#include "AliAnalysisTask.h"
+#include "AliAODInputHandler.h"
+#include "AliCDBEntry.h"
+#include "AliCDBManager.h"
+#include "AliCentralitySelectionTask.h"
+#include "AliCTPRawStream.h"
+#include "AliESDEvent.h"
+#include "AliESDfriend.h"
+#include "AliESDHeader.h"
+#include "AliESDInputHandler.h"
+#include "AliESDtrackCuts.h"
+#include "AliESDtrack.h"
+#include "AliESDVertex.h"
+#include "AliExternalTrackParam.h"
+#include "AliGeomManager.h"
+#include "AliGRPManager.h"
+#include "AliGRPObject.h"
+#include "AliLog.h"
+#include "AliMagF.h"
+#include "AliMathBase.h"
+#include "AliMCEventHandler.h"
+#include "AliPhysicsSelectionTask.h"
+#include "AliPID.h"
+#include "AliRawReaderRoot.h"
+#include "AliReconstruction.h"
+//#include "AliTPCcalibDB.h"
+#include "AliTPCclusterMI.h"
+//#include "AliTPCParam.h"
+//#include "AliTPCRecoParam.h"
+//#include "AliTPCseed.h"
+//#include "AliTPCTransform.h"
+#include "AliTrackerBase.h"
+#include "AliTracker.h"
+#include "AliTrackPointArray.h"
+#include "AliTRDdEdxBaseUtils.h"
+#include "AliTRDdigitsManager.h"
+#include "AliTRDrawStream.h"
+#include "AliTRDseedV1.h"
+#include "AliTriggerAnalysis.h"
+#include "AliSysInfo.h"
+
+
+#include "AliTRDdEdxParams.h"
+#include "AliOADBContainer.h"
+#include "AliLog.h"
+
+AliTRDdEdxParams * getDefault()
+{
+  //
+  //version 31/07/2014 by Xianguo Lu based on pPb small sample only for nch=6 ncls/nch>=18, no binning in eta, no binning on multiplicity
+  //
+
+  AliTRDdEdxParams *defobj=new AliTRDdEdxParams("default","default");
+
+  const Float_t meanpars[]={1.521e+00 , 7.294e-01 , 8.066e+00 , 2.146e-02 , 3.137e+01 , 7.160e-15 , 2.838e+00 , 1.100e+01};
+  const Float_t respars[]={-4.122e-01 , 1.019e+00 , -1.319e-01};
+
+  for(Int_t ipar=0; ipar< MAXNPAR; ipar++){
+    if(ipar == AliPID::kProton){
+      const Float_t tmpproton[]={2.026e+00 , -1.462e-04 , 1.202e+00 , 1.162e-01 , 2.092e+00 , -3.018e-02 , 3.665e+00 , 3.487e-07}; 
+      defobj->SetMeanParameter(ipar, sizeof(tmpproton)/sizeof(Float_t), tmpproton);
+    }
+    else if(ipar == AliPID::kPion || ipar ==AliPID::kElectron){       
+      const Float_t tmppe[]={1.508e+00 , 7.356e-01 , 8.002e+00 , 1.932e-02 , 3.058e+01 , 5.114e-18 , 3.561e+00 , 1.408e+01};  
+      defobj->SetMeanParameter(ipar, sizeof(tmppe)/sizeof(Float_t), tmppe);
+    }
+    else{
+      defobj->SetMeanParameter(ipar, sizeof(meanpars)/sizeof(Float_t), meanpars);
+    }
+
+    defobj->SetSigmaParameter(ipar,  sizeof(respars)/sizeof(Float_t), respars);
+  }
+
+  defobj->Print();
+
+  return defobj;
+}
+
+void makeTRDdEdxOADB(const Bool_t kMC=0)
+{
+  //
+  //make OADB object
+  //currently only default
+  //non-default has to be tested after finer inspection on dependence on eta, multiplicity, etc. 
+  //
+
+  AliTRDdEdxParams * defobj = getDefault();
+
+  TString containerName = "TRDdEdxParamsContainer";
+  AliOADBContainer* cont = new AliOADBContainer(containerName.Data());
+
+  //current convention as other detetor is only to provide for data
+  //no MC is seen in current (31/07/2014)  AliROOT except for HMPID
+
+  const TString filePathNamePackage=Form("$ALICE_ROOT/OADB/COMMON/PID/%s/TRDdEdxParams.root", kMC?"MC":"data");
+  Int_t statusCont = cont->InitFromFile(filePathNamePackage.Data(), cont->GetName());
+  if (statusCont) {
+    printf("No OADBContainer for the current settings found - creating a new one...\n");
+  }
+
+  //non-default params
+  /*
+  const Double_t runLow = 1;
+  const Double_t runUp =  10;
+  const Int_t index = cont->GetIndexForRun((runUp + runLow) / 2.0);
+  if (index < 0) {
+    printf("Creating new object for run range %d - %d...\n", runLow, runUp);
+    cont->AppendObject(hh, runLow, runUp);
+  }
+  else {
+    printf("Updating existing object for run range %d - %d...\n", runLow, runUp);
+    cont->UpdateObject(index, hh, runLow, runUp);
+  }
+  */
+
+  //default params
+  cont->CleanDefaultList(); // Remove old default objects at first
+  cont->AddDefaultObject(defobj);
+
+
+  TFile* f = new TFile(filePathNamePackage.Data(), "update");
+  f->Delete(cont->GetName());
+  cont->Write(0, TObject::kOverwrite);
+  f->Purge();
+  f->Close();
+
+  printf("makeTRDdEdxOADB done\n");
+}
+
+void readTRDdEdxOADB(const Bool_t kMC=0)
+{
+  //
+  //only read in default
+  //non-default has to be uncommented if needed
+  //
+
+  const TString containerName = "TRDdEdxParamsContainer";
+  AliOADBContainer cont(containerName.Data()); 
+  
+  const TString filePathNamePackage=Form("$ALICE_ROOT/OADB/COMMON/PID/%s/TRDdEdxParams.root", kMC?"MC":"data");
+
+  const Int_t statusCont = cont.InitFromFile(filePathNamePackage.Data(), cont.GetName());
+  if (statusCont){
+    printf("Failed initializing settings from OADB\n"); exit(1);
+  }
+
+  /*
+  Int_t fRun = 5;
+  AliTRDdEdxParams *jj=(AliTRDdEdxParams*)(cont.GetObject(fRun, defaultObj.Data()));
+  jj->Print();
+  */
+
+  TString defaultObj = "default";
+  AliTRDdEdxParams *jj=(AliTRDdEdxParams*)(cont.GetObject(1000, defaultObj.Data()));
+  jj->Print();
+
+}
+
+int main(int argc, char * argv[])
+{
+  for(Int_t ii=0; ii<argc; ii++){
+    printf("%d: %s\n", ii, argv[ii]);
+  }
+
+  if(argc!=2){
+    printf("argc!=2\n");
+    return 1; 
+  }
+
+  if(atoi(argv[1])){
+    makeTRDdEdxOADB();
+  }
+  else{
+    readTRDdEdxOADB();
+  }
+
+  return 0;
+}