//-----------------------------------------------------------------------------
/// \class AliMUONAlignment
-/// Alignment class for the ALICE DiMuon spectrometer
+/// Alignment class for the ALICE DiMuon spectrometer
///
-/// MUON specific alignment class which interface to AliMillepede.
+/// MUON specific alignment class which interface to AliMillepede.
/// For each track ProcessTrack calculates the local and global derivatives
/// at each cluster and fill the corresponding local equations. Provide methods
-/// for fixing or constraining detection elements for best results.
+/// for fixing or constraining detection elements for best results.
///
/// \author Bruce Becker, Javier Castillo
//-----------------------------------------------------------------------------
ClassImp(AliMUONAlignment)
/// \endcond
- Int_t AliMUONAlignment::fgNDetElem = 4*2+4*2+18*2+26*2+26*2;
- Int_t AliMUONAlignment::fgNDetElemCh[10] = {4,4,4,4,18,18,26,26,26,26};
- Int_t AliMUONAlignment::fgSNDetElemCh[10] = {4,8,12,16,34,52,78,104,130,156};
- Int_t AliMUONAlignment::fgNParCh = 4;
- Int_t AliMUONAlignment::fgNTrkMod = 16;
- Int_t AliMUONAlignment::fgNCh = 10;
- Int_t AliMUONAlignment::fgNSt = 5;
-
-AliMUONAlignment::AliMUONAlignment()
+//_____________________________________________________________________
+// static variables
+Int_t AliMUONAlignment::fgNDetElem = 4*2+4*2+18*2+26*2+26*2;
+Int_t AliMUONAlignment::fgNDetElemCh[10] = {4,4,4,4,18,18,26,26,26,26};
+Int_t AliMUONAlignment::fgSNDetElemCh[10] = {4,8,12,16,34,52,78,104,130,156};
+Int_t AliMUONAlignment::fgNParCh = 4;
+Int_t AliMUONAlignment::fgNTrkMod = 16;
+Int_t AliMUONAlignment::fgNCh = 10;
+Int_t AliMUONAlignment::fgNSt = 5;
+
+//_____________________________________________________________________
+AliMUONAlignment::AliMUONAlignment()
: TObject(),
fBFieldOn(kTRUE),
- fStartFac(256.),
- fResCutInitial(100.),
+ fStartFac(256.),
+ fResCutInitial(100.),
fResCut(100.),
fMillepede(0),
- fTrackParamAtCluster(0),
fTrack(0),
fCluster(0),
fTrackParam(0),
fPhi(0.),
fCosPhi(1.),
fSinPhi(0.),
+ fTrackRecord(),
fTransform(0)
{
- /// Default constructor, setting define alignment parameters
- fSigma[0] = 1.0e-1;
+
+ fSigma[0] = 1.5e-1;
fSigma[1] = 1.0e-2;
AliInfo(Form("fSigma[0]: %f\t fSigma[1]: %f",fSigma[0],fSigma[1]));
fDoF[0] = kTRUE; fDoF[1] = kTRUE; fDoF[2] = kTRUE; fDoF[3] = kTRUE;
fAllowVar[0] = 0.05; fAllowVar[1] = 0.05; fAllowVar[2] = 0.001; fAllowVar[3] = 0.5;
-
+
AliInfo(Form("fAllowVar[0]: %f\t fAllowVar[1]: %f\t fPhi: %f\t fgNDetElem: %i\t fNGlobal: %i\t fNLocal: %i",fAllowVar[0],fAllowVar[1],fPhi,fgNDetElem,fNGlobal,fNLocal));
fMillepede = new AliMillepede();
}
-AliMUONAlignment::~AliMUONAlignment() {
- /// Destructor
-}
+//_____________________________________________________________________
+AliMUONAlignment::~AliMUONAlignment()
+{}
-void AliMUONAlignment::Init(Int_t nGlobal, /* number of global paramers */
- Int_t nLocal, /* number of local parameters */
- Int_t nStdDev /* std dev cut */ )
+//_____________________________________________________________________
+void AliMUONAlignment::Init(
+ Int_t nGlobal, /* number of global paramers */
+ Int_t nLocal, /* number of local parameters */
+ Int_t nStdDev /* std dev cut */ )
{
/// Initialization of AliMillepede. Fix parameters, define constraints ...
fMillepede->InitMille(nGlobal,nLocal,nStdDev,fResCut,fResCutInitial);
-// Bool_t bStOnOff[5] = {kTRUE,kTRUE,kTRUE,kTRUE,kTRUE};
-// Bool_t bChOnOff[10] = {kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE};
-// Bool_t bSpecLROnOff[2] = {kTRUE,kTRUE};
+ // Bool_t bStOnOff[5] = {kTRUE,kTRUE,kTRUE,kTRUE,kTRUE};
+ // Bool_t bChOnOff[10] = {kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE};
+ // Bool_t bSpecLROnOff[2] = {kTRUE,kTRUE};
-// AllowVariations(bChOnOff);
+ // AllowVariations(bChOnOff);
// Fix parameters or add constraints here
-// for (Int_t iSt=0; iSt<5; iSt++)
-// if (!bStOnOff[iSt]) FixStation(iSt+1);
-// for (Int_t iCh=0; iCh<10; iCh++)
-// if (!bChOnOff[iCh]) FixChamber(iCh+1);
+ // for (Int_t iSt=0; iSt<5; iSt++)
+ // { if (!bStOnOff[iSt]) FixStation(iSt+1); }
+
+ // for (Int_t iCh=0; iCh<10; iCh++)
+ // { if (!bChOnOff[iCh]) FixChamber(iCh+1); }
// FixHalfSpectrometer(bChOnOff,bSpecLROnOff);
ResetConstraints();
-
+
// Define global constrains to be applied
// X, Y, P, XvsZ, YvsZ, PvsZ, XvsY, YvsY, PvsY
Bool_t bVarXYT[9] = {kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE};
// Other possible way to add constrains
bVarXYT[0] = kFALSE; bVarXYT[1] = kFALSE; bVarXYT[2] = kTRUE;
bDetTLBR[0] = kFALSE; bDetTLBR[1] = kTRUE; bDetTLBR[2] = kFALSE; bDetTLBR[3] = kFALSE;
-// AddConstraints(bStOnOff,bVarXYT,bDetTLBR);
+ // AddConstraints(bStOnOff,bVarXYT,bDetTLBR);
bVarXYT[0] = kTRUE; bVarXYT[1] = kTRUE; bVarXYT[2] = kFALSE;
// AddConstraints(bStOnOff,bVarXYT);
-
+
// Set iterations
- if (fStartFac>1) fMillepede->SetIterations(fStartFac);
+ if (fStartFac>1) fMillepede->SetIterations(fStartFac);
}
-void AliMUONAlignment::FixStation(Int_t iSt){
+//_____________________________________________________________________
+void AliMUONAlignment::FixStation(Int_t iSt)
+{
/// Fix all detection elements of station iSt
- Int_t iDetElemFirst = (iSt>1) ? fgSNDetElemCh[2*(iSt-1)-1] : 0;
- Int_t iDetElemLast = fgSNDetElemCh[2*(iSt)-1];
- for (Int_t i = iDetElemFirst; i < iDetElemLast; i++){
+ Int_t iDetElemFirst = (iSt>1) ? fgSNDetElemCh[2*(iSt-1)-1] : 0;
+ Int_t iDetElemLast = fgSNDetElemCh[2*(iSt)-1];
+ for (Int_t i = iDetElemFirst; i < iDetElemLast; i++)
+ {
FixParameter(i*fgNParCh+0, 0.0);
FixParameter(i*fgNParCh+1, 0.0);
FixParameter(i*fgNParCh+2, 0.0);
FixParameter(i*fgNParCh+3, 0.0);
}
+
}
-void AliMUONAlignment::FixChamber(Int_t iCh){
+//_____________________________________________________________________
+void AliMUONAlignment::FixChamber(Int_t iCh)
+{
/// Fix all detection elements of chamber iCh
- Int_t iDetElemFirst = (iCh>1) ? fgSNDetElemCh[iCh-2] : 0;
- Int_t iDetElemLast = fgSNDetElemCh[iCh-1];
- for (Int_t i = iDetElemFirst; i < iDetElemLast; i++){
+ Int_t iDetElemFirst = (iCh>1) ? fgSNDetElemCh[iCh-2] : 0;
+ Int_t iDetElemLast = fgSNDetElemCh[iCh-1];
+ for (Int_t i = iDetElemFirst; i < iDetElemLast; i++)
+ {
FixParameter(i*fgNParCh+0, 0.0);
FixParameter(i*fgNParCh+1, 0.0);
FixParameter(i*fgNParCh+2, 0.0);
}
}
-void AliMUONAlignment::FixDetElem(Int_t iDetElemId, TString sVarXYT){
+//_____________________________________________________________________
+void AliMUONAlignment::FixDetElem(Int_t iDetElemId, TString sVarXYT)
+{
+
/// Fix a given detection element
Int_t iDetElemNumber = iDetElemId%100;
- for (int iCh=0; iCh<iDetElemId/100-1; iCh++){
+ for (int iCh=0; iCh<iDetElemId/100-1; iCh++)
+ {
iDetElemNumber += fgNDetElemCh[iCh];
}
- if (sVarXYT.Contains("X")) { // X constraint
+
+ if (sVarXYT.Contains("X"))
+ {
+ // X constraint
FixParameter(iDetElemNumber*fgNParCh+0, 0.0);
}
- if (sVarXYT.Contains("Y")) { // Y constraint
+
+ if (sVarXYT.Contains("Y"))
+ {
+ // Y constraint
FixParameter(iDetElemNumber*fgNParCh+1, 0.0);
}
- if (sVarXYT.Contains("T")) { // T constraint
+
+ if (sVarXYT.Contains("T"))
+ {
+ // T constraint
FixParameter(iDetElemNumber*fgNParCh+2, 0.0);
}
- if (sVarXYT.Contains("Z")) { // T constraint
+
+ if (sVarXYT.Contains("Z"))
+ {
+ // T constraint
FixParameter(iDetElemNumber*fgNParCh+3, 0.0);
}
+
}
-void AliMUONAlignment::FixHalfSpectrometer(const Bool_t *lChOnOff, const Bool_t *lSpecLROnOff){
+//_____________________________________________________________________
+void AliMUONAlignment::FixHalfSpectrometer(const Bool_t *lChOnOff, const Bool_t *lSpecLROnOff)
+{
+
/// Fix left or right detector
- for (Int_t i = 0; i < fgNDetElem; i++){
+ for (Int_t i = 0; i < fgNDetElem; i++)
+ {
+
Int_t iCh=0;
- for (iCh=1; iCh<=fgNCh; iCh++){
- if (i<fgSNDetElemCh[iCh-1]) break;
- }
- if (lChOnOff[iCh-1]){
+ for (iCh=1; iCh<=fgNCh; iCh++)
+ { if (i<fgSNDetElemCh[iCh-1]) break; }
+
+ if (lChOnOff[iCh-1])
+ {
+
Int_t lDetElemNumber = (iCh==1) ? i : i-fgSNDetElemCh[iCh-2];
- if (iCh>=1 && iCh<=4){
- if ((lDetElemNumber==1 || lDetElemNumber==2) && !lSpecLROnOff[0]){ // From track crossings
- FixParameter(i*fgNParCh+0, 0.0);
- FixParameter(i*fgNParCh+1, 0.0);
- FixParameter(i*fgNParCh+2, 0.0);
- FixParameter(i*fgNParCh+3, 0.0);
- }
- if ((lDetElemNumber==0 || lDetElemNumber==3) && !lSpecLROnOff[1]){ // From track crossings
- FixParameter(i*fgNParCh+0, 0.0);
- FixParameter(i*fgNParCh+1, 0.0);
- FixParameter(i*fgNParCh+2, 0.0);
- FixParameter(i*fgNParCh+3, 0.0);
- }
+ if (iCh>=1 && iCh<=4)
+ {
+
+ if ((lDetElemNumber==1 || lDetElemNumber==2) && !lSpecLROnOff[0])
+ {
+ // From track crossings
+ FixParameter(i*fgNParCh+0, 0.0);
+ FixParameter(i*fgNParCh+1, 0.0);
+ FixParameter(i*fgNParCh+2, 0.0);
+ FixParameter(i*fgNParCh+3, 0.0);
+ }
+
+ if ((lDetElemNumber==0 || lDetElemNumber==3) && !lSpecLROnOff[1])
+ {
+ // From track crossings
+ FixParameter(i*fgNParCh+0, 0.0);
+ FixParameter(i*fgNParCh+1, 0.0);
+ FixParameter(i*fgNParCh+2, 0.0);
+ FixParameter(i*fgNParCh+3, 0.0);
+ }
+
}
- if (iCh>=5 && iCh<=6){
- if ((lDetElemNumber>=5&&lDetElemNumber<=13) && !lSpecLROnOff[0]){
- FixParameter(i*fgNParCh+0, 0.0);
- FixParameter(i*fgNParCh+1, 0.0);
- FixParameter(i*fgNParCh+2, 0.0);
- FixParameter(i*fgNParCh+3, 0.0);
- }
- if (((lDetElemNumber>=0&&lDetElemNumber<=4) ||
- (lDetElemNumber>=14&&lDetElemNumber<=17)) && !lSpecLROnOff[1]){
- FixParameter(i*fgNParCh+0, 0.0);
- FixParameter(i*fgNParCh+1, 0.0);
- FixParameter(i*fgNParCh+2, 0.0);
- FixParameter(i*fgNParCh+3, 0.0);
- }
+
+ if (iCh>=5 && iCh<=6)
+ {
+
+ if ((lDetElemNumber>=5&&lDetElemNumber<=13) && !lSpecLROnOff[0])
+ {
+ FixParameter(i*fgNParCh+0, 0.0);
+ FixParameter(i*fgNParCh+1, 0.0);
+ FixParameter(i*fgNParCh+2, 0.0);
+ FixParameter(i*fgNParCh+3, 0.0);
+ }
+
+ if (((lDetElemNumber>=0&&lDetElemNumber<=4) ||
+ (lDetElemNumber>=14&&lDetElemNumber<=17)) && !lSpecLROnOff[1])
+ {
+
+ FixParameter(i*fgNParCh+0, 0.0);
+ FixParameter(i*fgNParCh+1, 0.0);
+ FixParameter(i*fgNParCh+2, 0.0);
+ FixParameter(i*fgNParCh+3, 0.0);
+ }
+
}
- if (iCh>=7 && iCh<=10){
- if ((lDetElemNumber>=7&&lDetElemNumber<=19) && !lSpecLROnOff[0]){
- FixParameter(i*fgNParCh+0, 0.0);
- FixParameter(i*fgNParCh+1, 0.0);
- FixParameter(i*fgNParCh+2, 0.0);
- FixParameter(i*fgNParCh+3, 0.0);
- }
- if (((lDetElemNumber>=0&&lDetElemNumber<=6) ||
- (lDetElemNumber>=20&&lDetElemNumber<=25)) && !lSpecLROnOff[1]){
- FixParameter(i*fgNParCh+0, 0.0);
- FixParameter(i*fgNParCh+1, 0.0);
- FixParameter(i*fgNParCh+2, 0.0);
- FixParameter(i*fgNParCh+3, 0.0);
- }
+
+ if (iCh>=7 && iCh<=10)
+ {
+
+ if ((lDetElemNumber>=7&&lDetElemNumber<=19) && !lSpecLROnOff[0])
+ {
+ FixParameter(i*fgNParCh+0, 0.0);
+ FixParameter(i*fgNParCh+1, 0.0);
+ FixParameter(i*fgNParCh+2, 0.0);
+ FixParameter(i*fgNParCh+3, 0.0);
+ }
+
+ if (((lDetElemNumber>=0&&lDetElemNumber<=6) ||
+ (lDetElemNumber>=20&&lDetElemNumber<=25)) && !lSpecLROnOff[1])
+ {
+ FixParameter(i*fgNParCh+0, 0.0);
+ FixParameter(i*fgNParCh+1, 0.0);
+ FixParameter(i*fgNParCh+2, 0.0);
+ FixParameter(i*fgNParCh+3, 0.0);
+ }
+
}
+
}
+
}
+
}
-void AliMUONAlignment::SetNonLinear(const Bool_t *lChOnOff, const Bool_t *lVarXYT){
+//______________________________________________________________________
+void AliMUONAlignment::SetNonLinear(const Bool_t *lChOnOff, const Bool_t *lVarXYT)
+{
+
/// Set non linear parameter flag selected chambers and degrees of freedom
- for (Int_t i = 0; i < fgNDetElem; i++){
+ for (Int_t i = 0; i < fgNDetElem; i++)
+ {
+
Int_t iCh=0;
- for (iCh=1; iCh<=fgNCh; iCh++){
- if (i<fgSNDetElemCh[iCh-1]) break;
- }
- if (lChOnOff[iCh-1]){
- if (lVarXYT[0]) { // X constraint
- SetNonLinear(i*fgNParCh+0);
+ for (iCh=1; iCh<=fgNCh; iCh++)
+ { if (i<fgSNDetElemCh[iCh-1]) break; }
+
+ if (lChOnOff[iCh-1])
+ {
+
+ if (lVarXYT[0])
+ {
+ // X constraint
+ SetNonLinear(i*fgNParCh+0);
}
- if (lVarXYT[1]) { // Y constraint
- SetNonLinear(i*fgNParCh+1);
+
+ if (lVarXYT[1])
+ {
+ // Y constraint
+ SetNonLinear(i*fgNParCh+1);
}
- if (lVarXYT[2]) { // T constraint
- SetNonLinear(i*fgNParCh+2);
+
+ if (lVarXYT[2])
+ {
+ // T constraint
+ SetNonLinear(i*fgNParCh+2);
}
- if (lVarXYT[3]) { // Z constraint
- SetNonLinear(i*fgNParCh+3);
+
+ if (lVarXYT[3])
+ {
+ // Z constraint
+ SetNonLinear(i*fgNParCh+3);
}
+
}
+
}
+
}
-void AliMUONAlignment::AddConstraints(const Bool_t *lChOnOff, const Bool_t *lVarXYT){
- /// Add constraint equations for selected chambers and degrees of freedom
- for (Int_t i = 0; i < fgNDetElem; i++){
+//______________________________________________________________________
+void AliMUONAlignment::AddConstraints(const Bool_t *lChOnOff, const Bool_t *lVarXYT)
+{
+
+ /// Add constraint equations for selected chambers and degrees of freedom
+ for (Int_t i = 0; i < fgNDetElem; i++)
+ {
+
Int_t iCh=0;
- for (iCh=1; iCh<=fgNCh; iCh++){
+ for (iCh=1; iCh<=fgNCh; iCh++)
+ {
if (i<fgSNDetElemCh[iCh-1]) break;
}
- if (lChOnOff[iCh-1]){
- if (lVarXYT[0]) { // X constraint
- fConstraintX[i*fgNParCh+0]=1.0;
+
+ if (lChOnOff[iCh-1])
+ {
+ if (lVarXYT[0])
+ {
+ // X constraint
+ fConstraintX[i*fgNParCh+0]=1.0;
}
- if (lVarXYT[1]) { // Y constraint
- fConstraintY[i*fgNParCh+1]=1.0;
+
+ if (lVarXYT[1])
+ {
+ // Y constraint
+ fConstraintY[i*fgNParCh+1]=1.0;
}
- if (lVarXYT[2]) { // T constraint
- fConstraintP[i*fgNParCh+2]=1.0;
+
+ if (lVarXYT[2])
+ {
+ // T constraint
+ fConstraintP[i*fgNParCh+2]=1.0;
}
-// if (lVarXYT[3]) { // Z constraint
-// fConstraintP[i*fgNParCh+3]=1.0;
-// }
+ // if (lVarXYT[3]) { // Z constraint
+ // fConstraintP[i*fgNParCh+3]=1.0;
+ // }
+
}
}
- if (lVarXYT[0]) { // X constraint
+
+ if (lVarXYT[0])
+ {
+ // X constraint
AddConstraint(fConstraintX,0.0);
}
- if (lVarXYT[1]) { // Y constraint
+
+ if (lVarXYT[1])
+ {
+ // Y constraint
AddConstraint(fConstraintY,0.0);
}
- if (lVarXYT[2]) { // T constraint
+
+ if (lVarXYT[2])
+ {
+ // T constraint
AddConstraint(fConstraintP,0.0);
}
+
// if (lVarXYT[3]) { // Z constraint
// AddConstraint(fConstraintP,0.0);
// }
}
-void AliMUONAlignment::AddConstraints(const Bool_t *lChOnOff, const Bool_t *lVarXYT, const Bool_t *lDetTLBR, const Bool_t *lSpecLROnOff){
- /// Add constraint equations for selected chambers, degrees of freedom and detector half
+//______________________________________________________________________
+void AliMUONAlignment::AddConstraints(const Bool_t *lChOnOff, const Bool_t *lVarXYT, const Bool_t *lDetTLBR, const Bool_t *lSpecLROnOff)
+{
+ /// Add constraint equations for selected chambers, degrees of freedom and detector half
Double_t lDetElemLocX = 0.;
Double_t lDetElemLocY = 0.;
Double_t lDetElemLocZ = 0.;
Double_t lMeanZ = 0.;
Double_t lSigmaZ = 0.;
Int_t lNDetElem = 0;
- for (Int_t i = 0; i < fgNDetElem; i++){
+ for (Int_t i = 0; i < fgNDetElem; i++)
+ {
+
Int_t iCh=0;
for (iCh=1; iCh<=fgNCh; iCh++){
if (i<fgSNDetElemCh[iCh-1]) break;
}
- if (lChOnOff[iCh-1]){
+ if (lChOnOff[iCh-1]){
Int_t lDetElemNumber = (iCh==1) ? i : i-fgSNDetElemCh[iCh-2];
Int_t lDetElemId = iCh*100+lDetElemNumber;
fTransform->Local2Global(lDetElemId,lDetElemLocX,lDetElemLocY,lDetElemLocZ,
- lDetElemGloX,lDetElemGloY,lDetElemGloZ);
+ lDetElemGloX,lDetElemGloY,lDetElemGloZ);
if (iCh>=1 && iCh<=4){
- if ((lDetElemNumber==1 || lDetElemNumber==2) && lSpecLROnOff[0]){ // From track crossings
- lMeanY += lDetElemGloY;
- lSigmaY += lDetElemGloY*lDetElemGloY;
- lMeanZ += lDetElemGloZ;
- lSigmaZ += lDetElemGloZ*lDetElemGloZ;
- lNDetElem++;
- }
- if ((lDetElemNumber==0 || lDetElemNumber==3) && lSpecLROnOff[1]){ // From track crossings
- lMeanY += lDetElemGloY;
- lSigmaY += lDetElemGloY*lDetElemGloY;
- lMeanZ += lDetElemGloZ;
- lSigmaZ += lDetElemGloZ*lDetElemGloZ;
- lNDetElem++;
- }
+ if ((lDetElemNumber==1 || lDetElemNumber==2) && lSpecLROnOff[0]){ // From track crossings
+ lMeanY += lDetElemGloY;
+ lSigmaY += lDetElemGloY*lDetElemGloY;
+ lMeanZ += lDetElemGloZ;
+ lSigmaZ += lDetElemGloZ*lDetElemGloZ;
+ lNDetElem++;
+ }
+ if ((lDetElemNumber==0 || lDetElemNumber==3) && lSpecLROnOff[1]){ // From track crossings
+ lMeanY += lDetElemGloY;
+ lSigmaY += lDetElemGloY*lDetElemGloY;
+ lMeanZ += lDetElemGloZ;
+ lSigmaZ += lDetElemGloZ*lDetElemGloZ;
+ lNDetElem++;
+ }
}
if (iCh>=5 && iCh<=6){
- if ((lDetElemNumber>=5&&lDetElemNumber<=13) && lSpecLROnOff[0]){
- lMeanY += lDetElemGloY;
- lSigmaY += lDetElemGloY*lDetElemGloY;
- lMeanZ += lDetElemGloZ;
- lSigmaZ += lDetElemGloZ*lDetElemGloZ;
- lNDetElem++;
- }
- if (((lDetElemNumber>=0&&lDetElemNumber<=4) ||
- (lDetElemNumber>=14&&lDetElemNumber<=17)) && lSpecLROnOff[1]){
- lMeanY += lDetElemGloY;
- lSigmaY += lDetElemGloY*lDetElemGloY;
- lMeanZ += lDetElemGloZ;
- lSigmaZ += lDetElemGloZ*lDetElemGloZ;
- lNDetElem++;
- }
+ if ((lDetElemNumber>=5&&lDetElemNumber<=13) && lSpecLROnOff[0]){
+ lMeanY += lDetElemGloY;
+ lSigmaY += lDetElemGloY*lDetElemGloY;
+ lMeanZ += lDetElemGloZ;
+ lSigmaZ += lDetElemGloZ*lDetElemGloZ;
+ lNDetElem++;
+ }
+ if (((lDetElemNumber>=0&&lDetElemNumber<=4) ||
+ (lDetElemNumber>=14&&lDetElemNumber<=17)) && lSpecLROnOff[1]){
+ lMeanY += lDetElemGloY;
+ lSigmaY += lDetElemGloY*lDetElemGloY;
+ lMeanZ += lDetElemGloZ;
+ lSigmaZ += lDetElemGloZ*lDetElemGloZ;
+ lNDetElem++;
+ }
}
if (iCh>=7 && iCh<=10){
- if ((lDetElemNumber>=7&&lDetElemNumber<=19) && lSpecLROnOff[0]){
- lMeanY += lDetElemGloY;
- lSigmaY += lDetElemGloY*lDetElemGloY;
- lMeanZ += lDetElemGloZ;
- lSigmaZ += lDetElemGloZ*lDetElemGloZ;
- lNDetElem++;
- }
- if (((lDetElemNumber>=0&&lDetElemNumber<=6) ||
- (lDetElemNumber>=20&&lDetElemNumber<=25)) && lSpecLROnOff[1]){
- lMeanY += lDetElemGloY;
- lSigmaY += lDetElemGloY*lDetElemGloY;
- lMeanZ += lDetElemGloZ;
- lSigmaZ += lDetElemGloZ*lDetElemGloZ;
- lNDetElem++;
- }
+ if ((lDetElemNumber>=7&&lDetElemNumber<=19) && lSpecLROnOff[0]){
+ lMeanY += lDetElemGloY;
+ lSigmaY += lDetElemGloY*lDetElemGloY;
+ lMeanZ += lDetElemGloZ;
+ lSigmaZ += lDetElemGloZ*lDetElemGloZ;
+ lNDetElem++;
+ }
+ if (((lDetElemNumber>=0&&lDetElemNumber<=6) ||
+ (lDetElemNumber>=20&&lDetElemNumber<=25)) && lSpecLROnOff[1]){
+ lMeanY += lDetElemGloY;
+ lSigmaY += lDetElemGloY*lDetElemGloY;
+ lMeanZ += lDetElemGloZ;
+ lSigmaZ += lDetElemGloZ*lDetElemGloZ;
+ lNDetElem++;
+ }
}
}
}
lMeanZ /= lNDetElem;
lSigmaZ /= lNDetElem;
lSigmaZ = TMath::Sqrt(lSigmaZ-lMeanZ*lMeanZ);
- AliInfo(Form("Used %i DetElem, MeanZ= %f , SigmaZ= %f", lNDetElem,lMeanZ,lSigmaZ));
+ AliInfo(Form("Used %i DetElem, MeanZ= %f , SigmaZ= %f", lNDetElem,lMeanZ,lSigmaZ));
- for (Int_t i = 0; i < fgNDetElem; i++){
+ for (Int_t i = 0; i < fgNDetElem; i++){
Int_t iCh=0;
for (iCh=1; iCh<=fgNCh; iCh++){
if (i<fgSNDetElemCh[iCh-1]) break;
Int_t lDetElemNumber = (iCh==1) ? i : i-fgSNDetElemCh[iCh-2];
Int_t lDetElemId = iCh*100+lDetElemNumber;
fTransform->Local2Global(lDetElemId,lDetElemLocX,lDetElemLocY,lDetElemLocZ,
- lDetElemGloX,lDetElemGloY,lDetElemGloZ);
+ lDetElemGloX,lDetElemGloY,lDetElemGloZ);
if (lVarXYT[0]) { // X constraint
- if (lDetTLBR[0]) ConstrainT(i,iCh,fConstraintXT,0); // Top half
- if (lDetTLBR[1]) ConstrainL(i,iCh,fConstraintXL,0); // Left half
- if (lDetTLBR[2]) ConstrainB(i,iCh,fConstraintXB,0); // Bottom half
- if (lDetTLBR[3]) ConstrainR(i,iCh,fConstraintXR,0); // Right half
+ if (lDetTLBR[0]) ConstrainT(i,iCh,fConstraintXT,0); // Top half
+ if (lDetTLBR[1]) ConstrainL(i,iCh,fConstraintXL,0); // Left half
+ if (lDetTLBR[2]) ConstrainB(i,iCh,fConstraintXB,0); // Bottom half
+ if (lDetTLBR[3]) ConstrainR(i,iCh,fConstraintXR,0); // Right half
}
if (lVarXYT[1]) { // Y constraint
- if (lDetTLBR[0]) ConstrainT(i,iCh,fConstraintYT,1); // Top half
- if (lDetTLBR[1]) ConstrainL(i,iCh,fConstraintYL,1); // Left half
- if (lDetTLBR[2]) ConstrainB(i,iCh,fConstraintYB,1); // Bottom half
- if (lDetTLBR[3]) ConstrainR(i,iCh,fConstraintYR,1); // Right half
+ if (lDetTLBR[0]) ConstrainT(i,iCh,fConstraintYT,1); // Top half
+ if (lDetTLBR[1]) ConstrainL(i,iCh,fConstraintYL,1); // Left half
+ if (lDetTLBR[2]) ConstrainB(i,iCh,fConstraintYB,1); // Bottom half
+ if (lDetTLBR[3]) ConstrainR(i,iCh,fConstraintYR,1); // Right half
}
if (lVarXYT[2]) { // P constraint
- if (lDetTLBR[0]) ConstrainT(i,iCh,fConstraintPT,2); // Top half
- if (lDetTLBR[1]) ConstrainL(i,iCh,fConstraintPL,2); // Left half
- if (lDetTLBR[2]) ConstrainB(i,iCh,fConstraintPB,2); // Bottom half
- if (lDetTLBR[3]) ConstrainR(i,iCh,fConstraintPR,2); // Right half
+ if (lDetTLBR[0]) ConstrainT(i,iCh,fConstraintPT,2); // Top half
+ if (lDetTLBR[1]) ConstrainL(i,iCh,fConstraintPL,2); // Left half
+ if (lDetTLBR[2]) ConstrainB(i,iCh,fConstraintPB,2); // Bottom half
+ if (lDetTLBR[3]) ConstrainR(i,iCh,fConstraintPR,2); // Right half
}
if (lVarXYT[3]) { // X-Z shearing
- if (lDetTLBR[0]) ConstrainT(i,iCh,fConstraintXZT,0,(lDetElemGloZ-lMeanZ)/lSigmaZ); // Top half
- if (lDetTLBR[1]) ConstrainL(i,iCh,fConstraintXZL,0,(lDetElemGloZ-lMeanZ)/lSigmaZ); // Left half
- if (lDetTLBR[2]) ConstrainB(i,iCh,fConstraintXZB,0,(lDetElemGloZ-lMeanZ)/lSigmaZ); // Bottom half
- if (lDetTLBR[3]) ConstrainR(i,iCh,fConstraintXZR,0,(lDetElemGloZ-lMeanZ)/lSigmaZ); // Right half
+ if (lDetTLBR[0]) ConstrainT(i,iCh,fConstraintXZT,0,(lDetElemGloZ-lMeanZ)/lSigmaZ); // Top half
+ if (lDetTLBR[1]) ConstrainL(i,iCh,fConstraintXZL,0,(lDetElemGloZ-lMeanZ)/lSigmaZ); // Left half
+ if (lDetTLBR[2]) ConstrainB(i,iCh,fConstraintXZB,0,(lDetElemGloZ-lMeanZ)/lSigmaZ); // Bottom half
+ if (lDetTLBR[3]) ConstrainR(i,iCh,fConstraintXZR,0,(lDetElemGloZ-lMeanZ)/lSigmaZ); // Right half
}
if (lVarXYT[4]) { // Y-Z shearing
- if (lDetTLBR[0]) ConstrainT(i,iCh,fConstraintYZT,1,(lDetElemGloZ-lMeanZ)/lSigmaZ); // Top half
- if (lDetTLBR[1]) ConstrainL(i,iCh,fConstraintYZL,1,(lDetElemGloZ-lMeanZ)/lSigmaZ); // Left half
- if (lDetTLBR[2]) ConstrainB(i,iCh,fConstraintYZB,1,(lDetElemGloZ-lMeanZ)/lSigmaZ); // Bottom half
- if (lDetTLBR[3]) ConstrainR(i,iCh,fConstraintYZR,1,(lDetElemGloZ-lMeanZ)/lSigmaZ); // Right half
+ if (lDetTLBR[0]) ConstrainT(i,iCh,fConstraintYZT,1,(lDetElemGloZ-lMeanZ)/lSigmaZ); // Top half
+ if (lDetTLBR[1]) ConstrainL(i,iCh,fConstraintYZL,1,(lDetElemGloZ-lMeanZ)/lSigmaZ); // Left half
+ if (lDetTLBR[2]) ConstrainB(i,iCh,fConstraintYZB,1,(lDetElemGloZ-lMeanZ)/lSigmaZ); // Bottom half
+ if (lDetTLBR[3]) ConstrainR(i,iCh,fConstraintYZR,1,(lDetElemGloZ-lMeanZ)/lSigmaZ); // Right half
}
if (lVarXYT[5]) { // P-Z rotation
- if (lDetTLBR[0]) ConstrainT(i,iCh,fConstraintPZT,2,(lDetElemGloZ-lMeanZ)/lSigmaZ); // Top half
- if (lDetTLBR[1]) ConstrainL(i,iCh,fConstraintPZL,2,(lDetElemGloZ-lMeanZ)/lSigmaZ); // Left half
- if (lDetTLBR[2]) ConstrainB(i,iCh,fConstraintPZB,2,(lDetElemGloZ-lMeanZ)/lSigmaZ); // Bottom half
- if (lDetTLBR[3]) ConstrainR(i,iCh,fConstraintPZR,2,(lDetElemGloZ-lMeanZ)/lSigmaZ); // Right half
+ if (lDetTLBR[0]) ConstrainT(i,iCh,fConstraintPZT,2,(lDetElemGloZ-lMeanZ)/lSigmaZ); // Top half
+ if (lDetTLBR[1]) ConstrainL(i,iCh,fConstraintPZL,2,(lDetElemGloZ-lMeanZ)/lSigmaZ); // Left half
+ if (lDetTLBR[2]) ConstrainB(i,iCh,fConstraintPZB,2,(lDetElemGloZ-lMeanZ)/lSigmaZ); // Bottom half
+ if (lDetTLBR[3]) ConstrainR(i,iCh,fConstraintPZR,2,(lDetElemGloZ-lMeanZ)/lSigmaZ); // Right half
}
if (lVarXYT[6]) { // X-Y shearing
- if (lDetTLBR[0]) ConstrainT(i,iCh,fConstraintXYT,0,(lDetElemGloY-lMeanY)/lSigmaY); // Top half
- if (lDetTLBR[1]) ConstrainL(i,iCh,fConstraintXYL,0,(lDetElemGloY-lMeanY)/lSigmaY); // Left half
- if (lDetTLBR[2]) ConstrainB(i,iCh,fConstraintXYB,0,(lDetElemGloY-lMeanY)/lSigmaY); // Bottom half
- if (lDetTLBR[3]) ConstrainR(i,iCh,fConstraintXYR,0,(lDetElemGloY-lMeanY)/lSigmaY); // Right half
+ if (lDetTLBR[0]) ConstrainT(i,iCh,fConstraintXYT,0,(lDetElemGloY-lMeanY)/lSigmaY); // Top half
+ if (lDetTLBR[1]) ConstrainL(i,iCh,fConstraintXYL,0,(lDetElemGloY-lMeanY)/lSigmaY); // Left half
+ if (lDetTLBR[2]) ConstrainB(i,iCh,fConstraintXYB,0,(lDetElemGloY-lMeanY)/lSigmaY); // Bottom half
+ if (lDetTLBR[3]) ConstrainR(i,iCh,fConstraintXYR,0,(lDetElemGloY-lMeanY)/lSigmaY); // Right half
}
if (lVarXYT[7]) { // Y-Y scaling
- if (lDetTLBR[0]) ConstrainT(i,iCh,fConstraintYYT,1,(lDetElemGloY-lMeanY)/lSigmaY); // Top half
- if (lDetTLBR[1]) ConstrainL(i,iCh,fConstraintYYL,1,(lDetElemGloY-lMeanY)/lSigmaY); // Left half
- if (lDetTLBR[2]) ConstrainB(i,iCh,fConstraintYYB,1,(lDetElemGloY-lMeanY)/lSigmaY); // Bottom half
- if (lDetTLBR[3]) ConstrainR(i,iCh,fConstraintYYR,1,(lDetElemGloY-lMeanY)/lSigmaY); // Right half
+ if (lDetTLBR[0]) ConstrainT(i,iCh,fConstraintYYT,1,(lDetElemGloY-lMeanY)/lSigmaY); // Top half
+ if (lDetTLBR[1]) ConstrainL(i,iCh,fConstraintYYL,1,(lDetElemGloY-lMeanY)/lSigmaY); // Left half
+ if (lDetTLBR[2]) ConstrainB(i,iCh,fConstraintYYB,1,(lDetElemGloY-lMeanY)/lSigmaY); // Bottom half
+ if (lDetTLBR[3]) ConstrainR(i,iCh,fConstraintYYR,1,(lDetElemGloY-lMeanY)/lSigmaY); // Right half
}
if (lVarXYT[8]) { // P-Y rotation
- if (lDetTLBR[0]) ConstrainT(i,iCh,fConstraintPYT,2,(lDetElemGloY-lMeanY)/lSigmaY); // Top half
- if (lDetTLBR[1]) ConstrainL(i,iCh,fConstraintPYL,2,(lDetElemGloY-lMeanY)/lSigmaY); // Left half
- if (lDetTLBR[2]) ConstrainB(i,iCh,fConstraintPYB,2,(lDetElemGloY-lMeanY)/lSigmaY); // Bottom half
- if (lDetTLBR[3]) ConstrainR(i,iCh,fConstraintPYR,2,(lDetElemGloY-lMeanY)/lSigmaY); // Right half
+ if (lDetTLBR[0]) ConstrainT(i,iCh,fConstraintPYT,2,(lDetElemGloY-lMeanY)/lSigmaY); // Top half
+ if (lDetTLBR[1]) ConstrainL(i,iCh,fConstraintPYL,2,(lDetElemGloY-lMeanY)/lSigmaY); // Left half
+ if (lDetTLBR[2]) ConstrainB(i,iCh,fConstraintPYB,2,(lDetElemGloY-lMeanY)/lSigmaY); // Bottom half
+ if (lDetTLBR[3]) ConstrainR(i,iCh,fConstraintPYR,2,(lDetElemGloY-lMeanY)/lSigmaY); // Right half
}
}
}
}
}
+//______________________________________________________________________
void AliMUONAlignment::ConstrainL(Int_t lDetElem, Int_t lCh, Double_t *lConstraintL, Int_t iVar, Double_t lWeight) const{
/// Set constrain equation for left half of spectrometer
Int_t lDetElemNumber = (lCh==1) ? lDetElem : lDetElem-fgSNDetElemCh[lCh-2];
}
}
+//______________________________________________________________________
void AliMUONAlignment::ConstrainB(Int_t lDetElem, Int_t lCh, Double_t *lConstraintB, Int_t iVar, Double_t /*lWeight*/) const{
/// Set constrain equation for bottom half of spectrometer
Int_t lDetElemNumber = (lCh==1) ? lDetElem : lDetElem-fgSNDetElemCh[lCh-2];
}
}
if (lCh>=5 && lCh<=6){
- if ((lDetElemNumber>=9&&lDetElemNumber<=17) ||
- (lDetElemNumber==0)){
+ if ((lDetElemNumber>=9&&lDetElemNumber<=17) ||
+ (lDetElemNumber==0)){
lConstraintB[lDetElem*fgNParCh+iVar]=1.0;
}
}
if (lCh>=7 && lCh<=10){
- if ((lDetElemNumber>=13&&lDetElemNumber<=25) ||
- (lDetElemNumber==0)){
+ if ((lDetElemNumber>=13&&lDetElemNumber<=25) ||
+ (lDetElemNumber==0)){
lConstraintB[lDetElem*fgNParCh+iVar]=1.0;
}
}
}
+//______________________________________________________________________
void AliMUONAlignment::ConstrainR(Int_t lDetElem, Int_t lCh, Double_t *lConstraintR, Int_t iVar, Double_t lWeight) const{
/// Set constrain equation for right half of spectrometer
Int_t lDetElemNumber = (lCh==1) ? lDetElem : lDetElem-fgSNDetElemCh[lCh-2];
}
}
if (lCh>=5 && lCh<=6){
- if ((lDetElemNumber>=0&&lDetElemNumber<=4) ||
- (lDetElemNumber>=14&&lDetElemNumber<=17)){
+ if ((lDetElemNumber>=0&&lDetElemNumber<=4) ||
+ (lDetElemNumber>=14&&lDetElemNumber<=17)){
lConstraintR[lDetElem*fgNParCh+iVar]=lWeight;
}
}
if (lCh>=7 && lCh<=10){
- if ((lDetElemNumber>=0&&lDetElemNumber<=6) ||
- (lDetElemNumber>=20&&lDetElemNumber<=25)){
+ if ((lDetElemNumber>=0&&lDetElemNumber<=6) ||
+ (lDetElemNumber>=20&&lDetElemNumber<=25)){
lConstraintR[lDetElem*fgNParCh+iVar]=lWeight;
}
}
}
+//______________________________________________________________________
void AliMUONAlignment::ResetConstraints(){
/// Reset all constraint equations
- for (Int_t i = 0; i < fgNDetElem; i++){
+ for (Int_t i = 0; i < fgNDetElem; i++){
fConstraintX[i*fgNParCh+0]=0.0;
fConstraintX[i*fgNParCh+1]=0.0;
fConstraintX[i*fgNParCh+2]=0.0;
}
}
+//______________________________________________________________________
void AliMUONAlignment::AddConstraint(Double_t *par, Double_t value) {
/// Constrain equation defined by par to value
fMillepede->SetGlobalConstraint(par, value);
AliInfo("Adding constraint");
}
+//______________________________________________________________________
void AliMUONAlignment::InitGlobalParameters(Double_t *par) {
/// Initialize global parameters with par array
fMillepede->SetGlobalParameters(par);
AliInfo("Init Global Parameters");
}
-
+
+//______________________________________________________________________
void AliMUONAlignment::FixParameter(Int_t iPar, Double_t value) {
- /// Parameter iPar is encourage to vary in [-value;value].
+ /// Parameter iPar is encourage to vary in [-value;value].
/// If value == 0, parameter is fixed
fMillepede->SetParSigma(iPar, value);
if (TMath::Abs(value)<1e-4) AliInfo(Form("Parameter %i Fixed", iPar));
}
+//______________________________________________________________________
void AliMUONAlignment::ResetLocalEquation()
{
/// Reset the derivative vectors
}
}
-void AliMUONAlignment::AllowVariations(const Bool_t *bChOnOff) {
+//______________________________________________________________________
+void AliMUONAlignment::AllowVariations(const Bool_t *bChOnOff)
+{
+
/// Set allowed variation for selected chambers based on fDoF and fAllowVar
- for (Int_t iCh=1; iCh<=10; iCh++) {
- if (bChOnOff[iCh-1]) {
- Int_t iDetElemFirst = (iCh>1) ? fgSNDetElemCh[iCh-2] : 0;
- Int_t iDetElemLast = fgSNDetElemCh[iCh-1];
- for (int i=0; i<fgNParCh; i++) {
- AliDebug(1,Form("fDoF[%d]= %d",i,fDoF[i]));
- if (fDoF[i]) {
- for (Int_t j=iDetElemFirst; j<iDetElemLast; j++){
- FixParameter(j*fgNParCh+i, fAllowVar[i]);
- }
- }
+ for (Int_t iCh=1; iCh<=10; iCh++)
+ {
+ if (bChOnOff[iCh-1])
+ {
+
+ Int_t iDetElemFirst = (iCh>1) ? fgSNDetElemCh[iCh-2] : 0;
+ Int_t iDetElemLast = fgSNDetElemCh[iCh-1];
+ for (int i=0; i<fgNParCh; i++)
+ {
+ AliDebug(1,Form("fDoF[%d]= %d",i,fDoF[i]));
+ if (fDoF[i])
+ {
+ for (Int_t j=iDetElemFirst; j<iDetElemLast; j++){
+ FixParameter(j*fgNParCh+i, fAllowVar[i]);
+ }
+ }
+
}
+
}
+
}
+
}
-void AliMUONAlignment::SetNonLinear(Int_t iPar /* set non linear flag */ ) {
+//______________________________________________________________________
+void AliMUONAlignment::SetNonLinear(Int_t iPar /* set non linear flag */ )
+{
/// Set nonlinear flag for parameter iPar
fMillepede->SetNonLinear(iPar);
AliInfo(Form("Parameter %i set to non linear", iPar));
}
+//______________________________________________________________________
+void AliMUONAlignment::SetSigmaXY(Double_t sigmaX, Double_t sigmaY)
+{
-void AliMUONAlignment::SetSigmaXY(Double_t sigmaX, Double_t sigmaY) {
/// Set expected measurement resolution
fSigma[0] = sigmaX; fSigma[1] = sigmaY;
AliInfo(Form("Using fSigma[0]=%f and fSigma[1]=%f",fSigma[0],fSigma[1]));
}
-void AliMUONAlignment::LocalEquationX() {
- /// Define local equation for current track and cluster in x coor. measurement
- // set local derivatives
- SetLocalDerivative(0, fCosPhi);
- SetLocalDerivative(1, fCosPhi * (fTrackPos[2] - fTrackPos0[2]));
- SetLocalDerivative(2, fSinPhi);
- SetLocalDerivative(3, fSinPhi * (fTrackPos[2] - fTrackPos0[2]));
+//______________________________________________________________________
+void AliMUONAlignment::LocalEquationX( Bool_t doAlignment )
+{
+
+
+ // local cluster record
+ AliMUONAlignmentClusterRecord clusterRecord;
+
+ // store detector and measurement
+ clusterRecord.SetDetElemId( fDetElemId );
+ clusterRecord.SetDetElemNumber( fDetElemNumber );
+ clusterRecord.SetMeas( fMeas[0] );
+ clusterRecord.SetSigma( fSigma[0] );
+
+ // store local derivatives
+ clusterRecord.SetLocalDerivative( 0, fCosPhi );
+ clusterRecord.SetLocalDerivative( 1, fCosPhi*(fTrackPos[2] - fTrackPos0[2]) );
+ clusterRecord.SetLocalDerivative( 2, fSinPhi );
+ clusterRecord.SetLocalDerivative( 3, fSinPhi*(fTrackPos[2] - fTrackPos0[2]) );
+
+ // store global derivatives
+ clusterRecord.SetGlobalDerivative( 0, -fCosPhi );
+ clusterRecord.SetGlobalDerivative( 1, -fSinPhi );
+
+ if (fBFieldOn)
+ {
+
+ clusterRecord.SetGlobalDerivative(
+ 2,
+ -fSinPhi*(fTrackPos[0]-fDetElemPos[0])
+ +fCosPhi*(fTrackPos[1]-fDetElemPos[1]) );
+
+ } else {
+
+ clusterRecord.SetGlobalDerivative(
+ 2,
+ -fSinPhi*(fTrackPos0[0]+fTrackSlope0[0]*(fTrackPos[2]-fTrackPos0[2])-fDetElemPos[0])
+ +fCosPhi*(fTrackPos0[1]+fTrackSlope0[1]*(fTrackPos[2]-fTrackPos0[2])-fDetElemPos[1]) );
- // set global derivatives
- SetGlobalDerivative(fDetElemNumber*fgNParCh+0, -fCosPhi);
- SetGlobalDerivative(fDetElemNumber*fgNParCh+1, -fSinPhi);
- if (fBFieldOn){
- SetGlobalDerivative(fDetElemNumber*fgNParCh+2,
- -fSinPhi*(fTrackPos[0]-fDetElemPos[0])
- +fCosPhi*(fTrackPos[1]-fDetElemPos[1]));
- }
- else {
- SetGlobalDerivative(fDetElemNumber*fgNParCh+2,
- -fSinPhi*(fTrackPos0[0]+fTrackSlope0[0]*
- (fTrackPos[2]-fTrackPos0[2])-fDetElemPos[0])
- +fCosPhi*(fTrackPos0[1]+fTrackSlope0[1]*
- (fTrackPos[2]-fTrackPos0[2])-fDetElemPos[1]));
}
- SetGlobalDerivative(fDetElemNumber*fgNParCh+3,
- fCosPhi*fTrackSlope0[0]+fSinPhi*fTrackSlope0[1]);
- fMillepede->SetLocalEquation(fGlobalDerivatives, fLocalDerivatives, fMeas[0], fSigma[0]);
+ clusterRecord.SetGlobalDerivative( 3, fCosPhi*fTrackSlope0[0]+fSinPhi*fTrackSlope0[1] );
+
+ // append to trackRecord
+ fTrackRecord.AddClusterRecord( clusterRecord );
+
+ // store local equation
+ if( doAlignment ) LocalEquation( clusterRecord );
+
}
-void AliMUONAlignment::LocalEquationY() {
- /// Define local equation for current track and cluster in y coor. measurement
- // set local derivatives
- SetLocalDerivative(0,-fSinPhi);
- SetLocalDerivative(1,-fSinPhi * (fTrackPos[2] - fTrackPos0[2]));
- SetLocalDerivative(2, fCosPhi);
- SetLocalDerivative(3, fCosPhi * (fTrackPos[2] - fTrackPos0[2]));
+//______________________________________________________________________
+void AliMUONAlignment::LocalEquationY(Bool_t doAlignment )
+{
+
+ // local cluster record
+ AliMUONAlignmentClusterRecord clusterRecord;
+
+ // store detector and measurement
+ clusterRecord.SetDetElemId( fDetElemId );
+ clusterRecord.SetDetElemNumber( fDetElemNumber );
+ clusterRecord.SetMeas( fMeas[1] );
+ clusterRecord.SetSigma( fSigma[1] );
+
+ // store local derivatives
+ clusterRecord.SetLocalDerivative( 0, -fSinPhi );
+ clusterRecord.SetLocalDerivative( 1, -fSinPhi*(fTrackPos[2] - fTrackPos0[2] ) );
+ clusterRecord.SetLocalDerivative( 2, fCosPhi );
+ clusterRecord.SetLocalDerivative( 3, fCosPhi*(fTrackPos[2] - fTrackPos0[2] ) );
// set global derivatives
- SetGlobalDerivative(fDetElemNumber*fgNParCh+0, fSinPhi);
- SetGlobalDerivative(fDetElemNumber*fgNParCh+1, -fCosPhi);
- if (fBFieldOn){
- SetGlobalDerivative(fDetElemNumber*fgNParCh+2,
- -fCosPhi*(fTrackPos[0]-fDetElemPos[0])
- -fSinPhi*(fTrackPos[1]-fDetElemPos[1]));
+ clusterRecord.SetGlobalDerivative( 0, fSinPhi);
+ clusterRecord.SetGlobalDerivative( 1, -fCosPhi);
+
+ if (fBFieldOn)
+ {
+
+ clusterRecord.SetGlobalDerivative(
+ 2,
+ -fCosPhi*(fTrackPos[0]-fDetElemPos[0])
+ -fSinPhi*(fTrackPos[1]-fDetElemPos[1]));
+
+ } else {
+
+ clusterRecord.SetGlobalDerivative(
+ 2,
+ -fCosPhi*(fTrackPos0[0]+fTrackSlope0[0]*(fTrackPos[2]-fTrackPos0[2])-fDetElemPos[0])
+ -fSinPhi*(fTrackPos0[1]+fTrackSlope0[1]*(fTrackPos[2]-fTrackPos0[2])-fDetElemPos[1]));
}
- else {
- SetGlobalDerivative(fDetElemNumber*fgNParCh+2,
- -fCosPhi*(fTrackPos0[0]+fTrackSlope0[0]*
- (fTrackPos[2]-fTrackPos0[2])-fDetElemPos[0])
- -fSinPhi*(fTrackPos0[1]+fTrackSlope0[1]*
- (fTrackPos[2]-fTrackPos0[2])-fDetElemPos[1]));
+
+ clusterRecord.SetGlobalDerivative( 3, -fSinPhi*fTrackSlope0[0]+fCosPhi*fTrackSlope0[1]);
+
+ // append to trackRecord
+ fTrackRecord.AddClusterRecord( clusterRecord );
+
+ // store local equation
+ if( doAlignment ) LocalEquation( clusterRecord );
+
+}
+
+//_____________________________________________________
+void AliMUONAlignment::LocalEquation( const AliMUONAlignmentClusterRecord& clusterRecord )
+{
+
+ // copy to local variables
+ for( Int_t index = 0; index < 4; ++index )
+ {
+ SetLocalDerivative( index, clusterRecord.GetLocalDerivative( index ) );
+ SetGlobalDerivative( clusterRecord.GetDetElemNumber()*fgNParCh + index, clusterRecord.GetGlobalDerivative( index ) );
}
- SetGlobalDerivative(fDetElemNumber*fgNParCh+3,
- -fSinPhi*fTrackSlope0[0]+fCosPhi*fTrackSlope0[1]);
- fMillepede->SetLocalEquation(fGlobalDerivatives, fLocalDerivatives, fMeas[1], fSigma[1]);
+ // pass equation parameters to millepede
+ fMillepede->SetLocalEquation( fGlobalDerivatives, fLocalDerivatives, clusterRecord.GetMeas(), clusterRecord.GetSigma() );
+
}
-void AliMUONAlignment::FillRecPointData() {
+//______________________________________________________________________
+void AliMUONAlignment::FillRecPointData()
+{
+
/// Get information of current cluster
fClustPos[0] = fCluster->GetX();
fClustPos[1] = fCluster->GetY();
fClustPos[2] = fCluster->GetZ();
- fTransform->Global2Local(fDetElemId,fClustPos[0],fClustPos[1],fClustPos[2],
- fClustPosLoc[0],fClustPosLoc[1],fClustPosLoc[2]);
+ fTransform->Global2Local(
+ fDetElemId,fClustPos[0],fClustPos[1],fClustPos[2],
+ fClustPosLoc[0],fClustPosLoc[1],fClustPosLoc[2]);
+
}
-void AliMUONAlignment::FillTrackParamData() {
+//______________________________________________________________________
+void AliMUONAlignment::FillTrackParamData()
+{
+
/// Get information of current track at current cluster
fTrackPos[0] = fTrackParam->GetNonBendingCoor();
fTrackPos[1] = fTrackParam->GetBendingCoor();
fTrackPos[2] = fTrackParam->GetZ();
fTrackSlope[0] = fTrackParam->GetNonBendingSlope();
fTrackSlope[1] = fTrackParam->GetBendingSlope();
- fTransform->Global2Local(fDetElemId,fTrackPos[0],fTrackPos[1],fTrackPos[2],
- fTrackPosLoc[0],fTrackPosLoc[1],fTrackPosLoc[2]);
+ fTransform->Global2Local(
+ fDetElemId,fTrackPos[0],fTrackPos[1],fTrackPos[2],
+ fTrackPosLoc[0],fTrackPosLoc[1],fTrackPosLoc[2]);
+
}
-void AliMUONAlignment::FillDetElemData() {
+//_____________________________________________________
+void AliMUONAlignment::FillDetElemData()
+{
+
/// Get information of current detection element
Double_t lDetElemLocX = 0.;
Double_t lDetElemLocY = 0.;
Double_t lDetElemLocZ = 0.;
fDetElemId = fCluster->GetDetElemId();
fDetElemNumber = fDetElemId%100;
- for (int iCh=0; iCh<fDetElemId/100-1; iCh++){
- fDetElemNumber += fgNDetElemCh[iCh];
- }
- fTransform->Local2Global(fDetElemId,lDetElemLocX,lDetElemLocY,lDetElemLocZ,
- fDetElemPos[0],fDetElemPos[1],fDetElemPos[2]);
+ for (int iCh=0; iCh<fDetElemId/100-1; iCh++)
+ { fDetElemNumber += fgNDetElemCh[iCh]; }
+
+ fTransform->Local2Global(
+ fDetElemId,lDetElemLocX,lDetElemLocY,lDetElemLocZ,
+ fDetElemPos[0],fDetElemPos[1],fDetElemPos[2]);
+
}
-void AliMUONAlignment::ProcessTrack(AliMUONTrack * track) {
- /// Process track; Loop over clusters and set local equations
+//_____________________________________________________
+AliMUONAlignmentTrackRecord* AliMUONAlignment::ProcessTrack( AliMUONTrack* track, Bool_t doAlignment )
+{
+
+ // store current track in running member.
fTrack = track;
- // get tclones arrays.
- fTrackParamAtCluster = fTrack->GetTrackParamAtCluster();
-
- // get size of arrays
- Int_t nTrackParam = fTrackParamAtCluster->GetEntries();
+
+ // clear track record
+ fTrackRecord.Clear();
+
+ // get number of tracks
+ Int_t nTrackParam = fTrack->GetTrackParamAtCluster()->GetEntries();
AliDebug(1,Form("Number of track param entries : %i ", nTrackParam));
- for(Int_t iCluster=0; iCluster<nTrackParam; iCluster++) {
- fTrackParam = (AliMUONTrackParam *) fTrack->GetTrackParamAtCluster()->At(iCluster);
- if (!fTrackParam) continue;
- fCluster = fTrackParam->GetClusterPtr();
- if (!fCluster) continue;
- // fill local variables for this position --> one measurement
- FillDetElemData();
- FillRecPointData();
- FillTrackParamData();
-// if (fDetElemId<500) continue;
- fTrackPos0[0] = fTrackPos[0];
- fTrackPos0[1] = fTrackPos[1];
- fTrackPos0[2] = fTrackPos[2];
- fTrackSlope0[0] = fTrackSlope[0];
- fTrackSlope0[1] = fTrackSlope[1];
- break;
- }
+ Bool_t first( kTRUE );
+ for(Int_t iCluster=0; iCluster<nTrackParam; iCluster++)
+ {
- for(Int_t iCluster=0; iCluster<nTrackParam; iCluster++) {
// and get new pointers
fTrackParam = (AliMUONTrackParam *) fTrack->GetTrackParamAtCluster()->At(iCluster);
- if (!fTrackParam) continue;
fCluster = fTrackParam->GetClusterPtr();
- if (!fCluster) continue;
+ if (!fCluster || !fTrackParam) continue;
+ // if (fDetElemId<500) continue;
+
// fill local variables for this position --> one measurement
- FillDetElemData();
+ FillDetElemData();
FillRecPointData();
FillTrackParamData();
-// if (fDetElemId<500) continue;
- AliDebug(1,Form("cluster: %i", iCluster));
- AliDebug(1,Form("x: %f\t y: %f\t z: %f\t DetElemID: %i\t ", fClustPos[0], fClustPos[1], fClustPos[2], fDetElemId));
- AliDebug(1,Form("fDetElemPos[0]: %f\t fDetElemPos[1]: %f\t fDetElemPos[2]: %f\t DetElemID: %i\t ", fDetElemPos[0],fDetElemPos[1],fDetElemPos[2], fDetElemId));
- AliDebug(1,Form("Track Parameter: %i", iCluster));
- AliDebug(1,Form("x: %f\t y: %f\t z: %f\t slopex: %f\t slopey: %f", fTrackPos[0], fTrackPos[1], fTrackPos[2], fTrackSlope[0], fTrackSlope[1]));
- AliDebug(1,Form("x0: %f\t y0: %f\t z0: %f\t slopex0: %f\t slopey0: %f", fTrackPos0[0], fTrackPos0[1], fTrackPos0[2], fTrackSlope0[0], fTrackSlope0[1]));
-
+ if( first )
+ {
+
+ // for first valid cluster, save track position as "starting" values
+ first = kFALSE;
+
+ fTrackPos0[0] = fTrackPos[0];
+ fTrackPos0[1] = fTrackPos[1];
+ fTrackPos0[2] = fTrackPos[2];
+ fTrackSlope0[0] = fTrackSlope[0];
+ fTrackSlope0[1] = fTrackSlope[1];
+
+ }
+
+ // calculate measurements
fCosPhi = TMath::Cos(fPhi);
fSinPhi = TMath::Sin(fPhi);
- if (fBFieldOn){
+ if( fBFieldOn )
+ {
+
fMeas[0] = fTrackPos[0] - fClustPos[0];
fMeas[1] = fTrackPos[1] - fClustPos[1];
- }
- else {
+
+ } else {
+
fMeas[0] = - fClustPos[0];
fMeas[1] = - fClustPos[1];
+
}
- AliDebug(1,Form("fMeas[0]: %f\t fMeas[1]: %f\t fSigma[0]: %f\t fSigma[1]: %f", fMeas[0], fMeas[1], fSigma[0], fSigma[1]));
+
+ // soùe debugging
+ AliDebug(1,Form("cluster: %i", iCluster));
+ AliDebug(1,Form("x: %f\t y: %f\t z: %f\t DetElemID: %i\t ", fClustPos[0], fClustPos[1], fClustPos[2], fDetElemId));
+ AliDebug(1,Form("fDetElemPos[0]: %f\t fDetElemPos[1]: %f\t fDetElemPos[2]: %f\t DetElemID: %i\t ", fDetElemPos[0],fDetElemPos[1],fDetElemPos[2], fDetElemId));
+
+ AliDebug(1,Form("Track Parameter: %i", iCluster));
+ AliDebug(1,Form("x: %f\t y: %f\t z: %f\t slopex: %f\t slopey: %f", fTrackPos[0], fTrackPos[1], fTrackPos[2], fTrackSlope[0], fTrackSlope[1]));
+ AliDebug(1,Form("x0: %f\t y0: %f\t z0: %f\t slopex0: %f\t slopey0: %f", fTrackPos0[0], fTrackPos0[1], fTrackPos0[2], fTrackSlope0[0], fTrackSlope0[1]));
+
+ AliDebug(1,Form("fMeas[0]: %f\t fMeas[1]: %f\t fSigma[0]: %f\t fSigma[1]: %f", fMeas[0], fMeas[1], fSigma[0], fSigma[1]));
+
// Set local equations
- LocalEquationX();
- LocalEquationY();
+ LocalEquationX( doAlignment );
+ LocalEquationY( doAlignment );
+
}
+
+ return &fTrackRecord;
+
}
-void AliMUONAlignment::LocalFit(Int_t iTrack, Double_t *lTrackParam, Int_t lSingleFit) {
+//______________________________________________________________________________
+void AliMUONAlignment::ProcessTrack( AliMUONAlignmentTrackRecord* track, Bool_t doAlignment )
+{
+
+ if( !( track && doAlignment ) ) return;
+
+ // loop over clusters
+ for( Int_t index = 0; index < track->GetNRecords(); ++index )
+ { if( AliMUONAlignmentClusterRecord* record = track->GetRecord( index ) ) LocalEquation( *record ); }
+
+ return;
+
+}
+
+//_____________________________________________________
+void AliMUONAlignment::LocalFit(Int_t iTrack, Double_t *lTrackParam, Int_t lSingleFit)
+{
+
/// Call local fit for this tracks
Int_t iRes = fMillepede->LocalFit(iTrack,lTrackParam,lSingleFit);
- if (iRes && !lSingleFit) {
- fMillepede->SetNLocalEquations(fMillepede->GetNLocalEquations()+1);
- }
+ if (iRes && !lSingleFit)
+ { fMillepede->SetNLocalEquations(fMillepede->GetNLocalEquations()+1); }
+
}
-void AliMUONAlignment::GlobalFit(Double_t *parameters,Double_t *errors,Double_t *pulls) {
+//_____________________________________________________
+void AliMUONAlignment::GlobalFit(Double_t *parameters,Double_t *errors,Double_t *pulls)
+{
+
/// Call global fit; Global parameters are stored in parameters
fMillepede->GlobalFit(parameters,errors,pulls);
AliInfo("Done fitting global parameters!");
- for (int iGlob=0; iGlob<fgNDetElem; iGlob++){
- printf("%d\t %f\t %f\t %f \n",iGlob,parameters[iGlob*fgNParCh+0],parameters[iGlob*fgNParCh+1],parameters[iGlob*fgNParCh+2]);
- }
+ for (int iGlob=0; iGlob<fgNDetElem; iGlob++)
+ { printf("%d\t %f\t %f\t %f\t %f \n",iGlob,parameters[iGlob*fgNParCh+0],parameters[iGlob*fgNParCh+1],parameters[iGlob*fgNParCh+3],parameters[iGlob*fgNParCh+2]); }
+
+
}
-Double_t AliMUONAlignment::GetParError(Int_t iPar) {
+//_____________________________________________________
+Double_t AliMUONAlignment::GetParError(Int_t iPar)
+{
/// Get error of parameter iPar
Double_t lErr = fMillepede->GetParError(iPar);
return lErr;
}
-void AliMUONAlignment::PrintGlobalParameters() {
+//_____________________________________________________
+void AliMUONAlignment::PrintGlobalParameters()
+{
/// Print global parameters
fMillepede->PrintGlobalParameters();
}
TGeoCombiTrans AliMUONAlignment::ReAlign(const TGeoCombiTrans & transform, const double *lMisAlignment) const
{
/// Realign given transformation by given misalignment and return the misaligned transformation
-
+
Double_t cartMisAlig[3] = {0,0,0};
Double_t angMisAlig[3] = {0,0,0};
// const Double_t *trans = transform.GetTranslation();
// TGeoRotation *rot;
// // check if the rotation we obtain is not NULL
-// if (transform.GetRotation()) {
+// if (transform.GetRotation()) {
// rot = transform.GetRotation();
// }
-// else {
+// else {
// rot = new TGeoRotation("rot");
// } // default constructor.
//______________________________________________________________________
AliMUONGeometryTransformer *
AliMUONAlignment::ReAlign(const AliMUONGeometryTransformer * transformer,
- const double *misAlignments, Bool_t verbose)
-
+ const double *misAlignments, Bool_t verbose)
+
{
/// Returns a new AliMUONGeometryTransformer with the found misalignments
- /// applied.
+ /// applied.
// Takes the internal geometry module transformers, copies them
// and gets the Detection Elements from them.
// Takes misalignment parameters and applies these
// to the local transform of the Detection Element
// Obtains the global transform by multiplying the module transformer
- // transformation with the local transformation
+ // transformation with the local transformation
// Applies the global transform to a new detection element
// Adds the new detection element to a new module transformer
// Adds the new module transformer to a new geometry transformer
AliMUONGeometryTransformer *newGeometryTransformer =
new AliMUONGeometryTransformer();
for (Int_t iMt = 0; iMt < transformer->GetNofModuleTransformers(); iMt++) {
- // module transformers
+ // module transformers
const AliMUONGeometryModuleTransformer *kModuleTransformer =
transformer->GetModuleTransformer(iMt, true);
-
+
AliMUONGeometryModuleTransformer *newModuleTransformer =
new AliMUONGeometryModuleTransformer(iMt);
newGeometryTransformer->AddModuleTransformer(newModuleTransformer);
-
+
TGeoCombiTrans moduleTransform =
TGeoCombiTrans(*kModuleTransformer->GetTransformation());
// New module transformation
TGeoCombiTrans newModuleTransform = ReAlign(moduleTransform,lModuleMisAlignment);
newModuleTransformer->SetTransformation(newModuleTransform);
-
- // Get delta transformation:
+
+ // Get delta transformation:
// Tdelta = Tnew * Told.inverse
- TGeoHMatrix deltaModuleTransform =
- AliMUONGeometryBuilder::Multiply(newModuleTransform,
- kModuleTransformer->GetTransformation()->Inverse());
+ TGeoHMatrix deltaModuleTransform =
+ AliMUONGeometryBuilder::Multiply(newModuleTransform,
+ kModuleTransformer->GetTransformation()->Inverse());
// Create module mis alignment matrix
newGeometryTransformer
->AddMisAlignModule(kModuleTransformer->GetModuleId(), deltaModuleTransform);
-
+
AliMpExMap *detElements = kModuleTransformer->GetDetElementStore();
-
+
if (verbose)
AliInfo(Form("%i DEs in old GeometryStore %i",detElements->GetSize(), iMt));
++iDe;
// make a new detection element
AliMUONGeometryDetElement *newDetElement =
- new AliMUONGeometryDetElement(detElement->GetId(),
- detElement->GetVolumePath());
+ new AliMUONGeometryDetElement(detElement->GetId(),
+ detElement->GetVolumePath());
TString lDetElemName(detElement->GetDEName());
lDetElemName.ReplaceAll("DE","");
iDetElemId = lDetElemName.Atoi();
iDetElemNumber = iDetElemId%100;
for (int iCh=0; iCh<iDetElemId/100-1; iCh++){
- iDetElemNumber += fgNDetElemCh[iCh];
+ iDetElemNumber += fgNDetElemCh[iCh];
}
for (int i=0; i<fgNParCh; i++) {
- lDetElemMisAlignment[i] = 0.0;
- if (iMt<fgNTrkMod) {
- AliInfo(Form("iMt %i, iCh %i, iDe %i, iDeId %i, iDeNb %i, iPar %i",iMt, iDetElemId/100, iDe, iDetElemId, iDetElemNumber, iDetElemNumber*fgNParCh+i));
- lDetElemMisAlignment[i] = misAlignments[iDetElemNumber*fgNParCh+i];
- }
+ lDetElemMisAlignment[i] = 0.0;
+ if (iMt<fgNTrkMod) {
+ AliInfo(Form("iMt %i, iCh %i, iDe %i, iDeId %i, iDeNb %i, iPar %i",iMt, iDetElemId/100, iDe, iDetElemId, iDetElemNumber, iDetElemNumber*fgNParCh+i));
+ lDetElemMisAlignment[i] = misAlignments[iDetElemNumber*fgNParCh+i];
+ }
}
// local transformation of this detection element.
TGeoCombiTrans localTransform
- = TGeoCombiTrans(*detElement->GetLocalTransformation());
+ = TGeoCombiTrans(*detElement->GetLocalTransformation());
TGeoCombiTrans newLocalTransform = ReAlign(localTransform,lDetElemMisAlignment);
- newDetElement->SetLocalTransformation(newLocalTransform);
+ newDetElement->SetLocalTransformation(newLocalTransform);
// global transformation
TGeoHMatrix newGlobalTransform =
- AliMUONGeometryBuilder::Multiply(newModuleTransform,
- newLocalTransform);
+ AliMUONGeometryBuilder::Multiply(newModuleTransform,
+ newLocalTransform);
newDetElement->SetGlobalTransformation(newGlobalTransform);
-
+
// add this det element to module
newModuleTransformer->GetDetElementStore()->Add(newDetElement->GetId(),
- newDetElement);
+ newDetElement);
// In the Alice Alignment Framework misalignment objects store
// global delta transformation
// Get detection "intermediate" global transformation
TGeoHMatrix newOldGlobalTransform = newModuleTransform * localTransform;
- // Get detection element global delta transformation:
+ // Get detection element global delta transformation:
// Tdelta = Tnew * Told.inverse
TGeoHMatrix deltaGlobalTransform
- = AliMUONGeometryBuilder::Multiply(newGlobalTransform,
- newOldGlobalTransform.Inverse());
-
+ = AliMUONGeometryBuilder::Multiply(newGlobalTransform,
+ newOldGlobalTransform.Inverse());
+
// Create mis alignment matrix
newGeometryTransformer
- ->AddMisAlignDetElement(detElement->GetId(), deltaGlobalTransform);
+ ->AddMisAlignDetElement(detElement->GetId(), deltaGlobalTransform);
}
-
+
if (verbose)
AliInfo(Form("Added module transformer %i to the transformer", iMt));
newGeometryTransformer->AddModuleTransformer(newModuleTransformer);
}
//______________________________________________________________________
-void AliMUONAlignment::SetAlignmentResolution(const TClonesArray* misAlignArray, Int_t rChId, Double_t rChResX, Double_t rChResY, Double_t rDeResX, Double_t rDeResY){
+void AliMUONAlignment::SetAlignmentResolution(const TClonesArray* misAlignArray, Int_t rChId, Double_t rChResX, Double_t rChResY, Double_t rDeResX, Double_t rDeResY)
+{
+
/// Set alignment resolution to misalign objects to be stored in CDB
Int_t chIdMin = (rChId<0)? 0 : rChId;
Int_t chIdMax = (rChId<0)? 9 : rChId;
chName1 = Form("GM%d",4+(chId-4)*2);
chName2 = Form("GM%d",4+(chId-4)*2+1);
}
-
+
for (int i=0; i<misAlignArray->GetEntries(); i++) {
alignMat = (AliAlignObjMatrix*)misAlignArray->At(i);
TString volName(alignMat->GetSymName());
if((volName.Contains(chName1)&&
- ((volName.Last('/')==volName.Index(chName1)+chName1.Length())||
- (volName.Length()==volName.Index(chName1)+chName1.Length())))||
- (volName.Contains(chName2)&&
- ((volName.Last('/')==volName.Index(chName2)+chName2.Length())||
- (volName.Length()==volName.Index(chName2)+chName2.Length())))){
- volName.Remove(0,volName.Last('/')+1);
- if (volName.Contains("GM")) {
- // alignMat->Print("NULL");
- alignMat->SetCorrMatrix(mChCorrMatrix);
- } else if (volName.Contains("DE")) {
- // alignMat->Print("NULL");
- alignMat->SetCorrMatrix(mDECorrMatrix);
- }
+ ((volName.Last('/')==volName.Index(chName1)+chName1.Length())||
+ (volName.Length()==volName.Index(chName1)+chName1.Length())))||
+ (volName.Contains(chName2)&&
+ ((volName.Last('/')==volName.Index(chName2)+chName2.Length())||
+ (volName.Length()==volName.Index(chName2)+chName2.Length())))){
+ volName.Remove(0,volName.Last('/')+1);
+ if (volName.Contains("GM")) {
+ // alignMat->Print("NULL");
+ alignMat->SetCorrMatrix(mChCorrMatrix);
+ } else if (volName.Contains("DE")) {
+ // alignMat->Print("NULL");
+ alignMat->SetCorrMatrix(mDECorrMatrix);
+ }
}
}
}
/// \author Javier Castillo, CEA/Saclay - Irfu/SPhN
//-----------------------------------------------------------------------------
-#include <fstream>
-
-#include <TString.h>
-#include <TError.h>
-#include <TGraphErrors.h>
-#include <TTree.h>
-#include <TChain.h>
-#include <TClonesArray.h>
-#include <TGeoGlobalMagField.h>
-#include <TGeoManager.h>
-#include <Riostream.h>
+// this class header must always come first
+#include "AliMUONAlignmentTask.h"
-#include "AliAnalysisTask.h"
+// local header must come before system headers
#include "AliAnalysisManager.h"
+#include "AliAODInputHandler.h"
+#include "AliAODHandler.h"
+#include "AliAODEvent.h"
+#include "AliAODHeader.h"
#include "AliESDInputHandler.h"
#include "AliESDEvent.h"
#include "AliESDMuonTrack.h"
#include "AliMagF.h"
+#include "AliCDBStorage.h"
#include "AliCDBManager.h"
#include "AliGRPManager.h"
#include "AliCDBMetaData.h"
#include "AliGeomManager.h"
#include "AliMpCDB.h"
+#include "AliMUONCDB.h"
#include "AliMUONAlignment.h"
#include "AliMUONTrack.h"
#include "AliMUONTrackExtrap.h"
#include "AliMUONGeometryTransformer.h"
#include "AliMUONESDInterface.h"
-#include "AliMUONAlignmentTask.h"
+// system headers
+#include <fstream>
+#include <TString.h>
+#include <TError.h>
+#include <TGraphErrors.h>
+#include <TTree.h>
+#include <TChain.h>
+#include <TClonesArray.h>
+#include <TGeoGlobalMagField.h>
+#include <TGeoManager.h>
+#include <Riostream.h>
-///\cond CLASSIMP
+///\cond CLASSIMP
ClassImp(AliMUONAlignmentTask)
///\endcond
-// //________________________________________________________________________
-// AliMUONAlignmentTask::AliMUONAlignmentTask(const char *name)
-// : AliAnalysisTask(name, ""),
-// fESD(0x0),
-// fAlign(0x0),
-// fGeoFilename(""),
-// fTransform(0x0),
-// fTrackTot(0),
-// fTrackOk(0),
-// fMSDEx(0x0),
-// fMSDEy(0x0),
-// fMSDEz(0x0),
-// fMSDEp(0x0),
-// fList(0x0)
-// {
-// /// Default Constructor
-// // Define input and output slots here
-// // Input slot #0 works with a TChain
-// DefineInput(0, TChain::Class());
-// // Output slot #0 writes NTuple/histos into a TList
-// DefineOutput(0, TList::Class());
-
-// // initialize parameters ...
-// for(Int_t k=0;k<4*156;k++) {
-// fParameters[k]=0.;
-// fErrors[k]=0.;
-// fPulls[k]=0.;
-// }
-
-// fAlign = new AliMUONAlignment();
-// fTransform = new AliMUONGeometryTransformer();
-// }
-
//________________________________________________________________________
-AliMUONAlignmentTask::AliMUONAlignmentTask(const char *name, const char *geofilename, const char *defaultocdb, const char *misalignocdb)
- : AliAnalysisTask(name, ""),
- fESD(0x0),
- fAlign(0x0),
- fGeoFilename(geofilename),
- fMisAlignOCDB(misalignocdb),
- fDefaultOCDB(defaultocdb),
- fTransform(0x0),
- fTrackTot(0),
- fTrackOk(0),
- fLastRunNumber(-1),
- fMSDEx(0x0),
- fMSDEy(0x0),
- fMSDEz(0x0),
- fMSDEp(0x0),
- fList(0x0)
+AliMUONAlignmentTask::AliMUONAlignmentTask(const char *name, const char *newalignocdb, const char *oldalignocdb, const char *defaultocdb, const char *geofilename):
+ AliAnalysisTaskSE(name),
+ fReadRecords( kFALSE ),
+ fWriteRecords( kFALSE ),
+ fDoAlignment( kTRUE ),
+ fAlign(0x0),
+ fGeoFilename(geofilename),
+ fDefaultStorage(defaultocdb),
+ fOldAlignStorage(oldalignocdb),
+ fNewAlignStorage(newalignocdb),
+ fOldGeoTransformer(NULL),
+ fNewGeoTransformer(NULL),
+ fLoadOCDBOnce(kFALSE),
+ fOCDBLoaded(kFALSE),
+ fTrackTot(0),
+ fTrackOk(0),
+ fLastRunNumber(-1),
+ fMSDEx(0x0),
+ fMSDEy(0x0),
+ fMSDEz(0x0),
+ fMSDEp(0x0),
+ fList(0x0),
+ fRecords(0x0),
+ fRecordCount(0)
{
/// Default Constructor
// Define input and output slots here
// Input slot #0 works with a TChain
DefineInput(0, TChain::Class());
+
// Output slot #0 writes NTuple/histos into a TList
- DefineOutput(0, TList::Class());
+ DefineOutput(1, TList::Class());
// initialize parameters ...
- for(Int_t k=0;k<4*156;k++) {
+ for(Int_t k=0;k<4*156;k++)
+ {
fParameters[k]=0.;
fErrors[k]=0.;
fPulls[k]=0.;
}
- fAlign = new AliMUONAlignment();
- fTransform = new AliMUONGeometryTransformer();
-}
+ fOldGeoTransformer = new AliMUONGeometryTransformer();
+}
//________________________________________________________________________
-AliMUONAlignmentTask::AliMUONAlignmentTask(const AliMUONAlignmentTask& obj)
- : AliAnalysisTask(obj),
- fESD(0x0),
- fAlign(0x0),
- fGeoFilename(""),
- fMisAlignOCDB(""),
- fDefaultOCDB(""),
- fTransform(0x0),
- fTrackTot(0),
- fTrackOk(0),
- fLastRunNumber(-1),
- fMSDEx(0x0),
- fMSDEy(0x0),
- fMSDEz(0x0),
- fMSDEp(0x0),
- fList(0x0)
+AliMUONAlignmentTask::AliMUONAlignmentTask(const AliMUONAlignmentTask& other):
+ AliAnalysisTaskSE(other),
+ fReadRecords( other.fReadRecords ),
+ fWriteRecords( other.fWriteRecords ),
+ fDoAlignment( other.fDoAlignment ),
+ fAlign( fAlign ),
+ fGeoFilename( other.fGeoFilename ),
+ fDefaultStorage( other.fDefaultStorage ),
+ fOldAlignStorage( other.fOldAlignStorage ),
+ fNewAlignStorage( other.fNewAlignStorage ),
+ fOldGeoTransformer( other.fOldGeoTransformer ),
+ fNewGeoTransformer( other.fNewGeoTransformer ),
+ fLoadOCDBOnce( other.fLoadOCDBOnce ),
+ fOCDBLoaded( other.fOCDBLoaded ),
+ fTrackTot( other.fTrackTot ),
+ fTrackOk( other.fTrackOk ),
+ fLastRunNumber( other.fLastRunNumber ),
+ fMSDEx( other.fMSDEx ),
+ fMSDEy( other.fMSDEy ),
+ fMSDEz( other.fMSDEz ),
+ fMSDEp( other.fMSDEp ),
+ fList( other.fList ),
+ fRecords( other.fRecords ),
+ fRecordCount( other.fRecordCount )
{
- /// Copy constructor
- fESD = obj.fESD;
- fAlign = obj.fAlign;
- fGeoFilename = obj.fGeoFilename;
- fTransform = obj.fTransform;
- fTrackTot = obj.fTrackTot;
- fTrackOk = obj.fTrackOk;
- fLastRunNumber = obj.fLastRunNumber;
- fMSDEx = obj.fMSDEx;
- fMSDEy = obj.fMSDEy;
- fMSDEz = obj.fMSDEz;
- fMSDEp = obj.fMSDEp;
- fList = obj.fList;
- // initialize parameters ...
- for(Int_t k=0;k<4*156;k++) {
- fParameters[k]=obj.fParameters[k];
- fErrors[k]=obj.fErrors[k];
- fPulls[k]=obj.fPulls[k];
+ // initialize parameters
+ for(Int_t k=0;k<4*156;k++)
+ {
+
+ fParameters[k]=other.fParameters[k];
+ fErrors[k]=other.fErrors[k];
+ fPulls[k]=other.fPulls[k];
+
}
}
//________________________________________________________________________
-AliMUONAlignmentTask& AliMUONAlignmentTask::operator=(const AliMUONAlignmentTask& other)
+AliMUONAlignmentTask& AliMUONAlignmentTask::operator=(const AliMUONAlignmentTask& other)
{
/// Assignment
- AliAnalysisTask::operator=(other);
- fESD = other.fESD;
+ AliAnalysisTaskSE::operator=(other);
+
+ fReadRecords = other.fReadRecords;
+ fWriteRecords = other.fWriteRecords;
+ fDoAlignment = other.fDoAlignment;
+
+ // this breaks in destructor
fAlign = other.fAlign;
+
fGeoFilename = other.fGeoFilename;
- fMisAlignOCDB = other.fMisAlignOCDB;
- fDefaultOCDB = other.fDefaultOCDB;
- fTransform = other.fTransform;
- fTrackTot = other.fTrackTot;
- fTrackOk = other.fTrackOk;
+ fDefaultStorage = other.fDefaultStorage;
+ fOldAlignStorage = other.fOldAlignStorage;
+ fNewAlignStorage = other.fNewAlignStorage;
+
+ // this breaks in destructor
+ fOldGeoTransformer = other.fOldGeoTransformer;
+ fNewGeoTransformer = other.fNewGeoTransformer;
+
+ fLoadOCDBOnce = other.fLoadOCDBOnce;
+ fOCDBLoaded = other.fOCDBLoaded;
+ fTrackTot = other.fTrackTot;
+ fTrackOk = other.fTrackOk;
fLastRunNumber = other.fLastRunNumber;
- fMSDEx = other.fMSDEx;
- fMSDEy = other.fMSDEy;
+ fMSDEx = other.fMSDEx;
+ fMSDEy = other.fMSDEy;
fMSDEz = other.fMSDEz;
fMSDEp = other.fMSDEp;
- fList = other.fList;
+ fList = other.fList;
+ fRecords = other.fRecords;
+ fRecordCount = other.fRecordCount;
// initialize parameters ...
- for(Int_t k=0;k<4*156;k++) {
+ for( Int_t k=0; k<4*156; ++k)
+ {
+
fParameters[k]=other.fParameters[k];
fErrors[k]=other.fErrors[k];
fPulls[k]=other.fPulls[k];
+
}
-
+
return *this;
}
//________________________________________________________________________
-AliMUONAlignmentTask::~AliMUONAlignmentTask()
-{
- /// Destructor
- if (fAlign) delete fAlign;
- if (fTransform) delete fTransform;
+AliMUONAlignmentTask::~AliMUONAlignmentTask()
+{
+ /*
+ it is safe to delete NULL pointers, so
+ there is no need to test their validity here.
+ However, it crashes here, because of incorrect assignment operator
+ and copy constructor, resulting in double deletion of objects.
+ Would require deep copy instead.
+ */
+ delete fAlign;
+ delete fOldGeoTransformer;
+ delete fNewGeoTransformer;
}
//________________________________________________________________________
-void AliMUONAlignmentTask::LocalInit()
+void AliMUONAlignmentTask::LocalInit()
{
- /// Local initialization, called once per task on the client machine
+ /// Local initialization, called once per task on the client machine
/// where the analysis train is assembled
- fLastRunNumber = 0;
- // Prepare(fGeoFilename.Data(),fDefaultOCDB.Data(),fMisAlignOCDB.Data());
- Prepare(fGeoFilename.Data(),"local://$ALICE_ROOT/OCDB",fMisAlignOCDB.Data());
- fLastRunNumber = -1;
+
+ // print configuration
+ AliInfo( Form( "fReadRecords: %s", (fReadRecords ? "kTRUE":"kFALSE" ) ) );
+ AliInfo( Form( "fWriteRecords: %s", (fWriteRecords ? "kTRUE":"kFALSE" ) ) );
+ AliInfo( Form( "fDoAlignment: %s", (fDoAlignment ? "kTRUE":"kFALSE" ) ) );
+
+ // consistency checks between flags
+ /* need at least one of the flags set to true */
+ if( !( fReadRecords || fWriteRecords || fDoAlignment ) )
+ { AliFatal( "Need at least one of the three flags fReadRecords, fWriteRecords and fDoAlignment set to kTRUE" ); }
+
+ /* cannot read and write records at the same time */
+ if( fReadRecords && fWriteRecords )
+ { AliFatal( "Cannot set both fReadRecords and fWriteRecords to kTRUE at the same time" ); }
+
+ /* must run alignment when reading records */
+ if( fReadRecords && !fDoAlignment )
+ {
+
+ AliInfo( "Automatically setting fDoAlignment to kTRUE because fReadRecords is kTRUE" );
+ SetDoAlignment( kTRUE );
+
+ }
// Set initial values here, good guess may help convergence
- // St 1
+ // St 1
// Int_t iPar = 0;
- // fParameters[iPar++] = 0.010300 ; fParameters[iPar++] = 0.010600 ; fParameters[iPar++] = 0.000396 ;
-
+ // fParameters[iPar++] = 0.010300 ; fParameters[iPar++] = 0.010600 ; fParameters[iPar++] = 0.000396 ;
+ fAlign = new AliMUONAlignment();
fAlign->InitGlobalParameters(fParameters);
-
- fTransform->LoadGeometryData();
- fAlign->SetGeometryTransformer(fTransform);
+// AliCDBManager::Instance()->Print();
+//
+// fAlign->SetGeometryTransformer(fOldGeoTransformer);
// Do alignment with magnetic field off
fAlign->SetBFieldOn(kFALSE);
-
+
// Set tracking station to use
// Bool_t bStOnOff[5] = {kTRUE,kTRUE,kTRUE,kTRUE,kTRUE};
- Bool_t bChOnOff[10] = {kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kFALSE,kTRUE,kTRUE,kTRUE,kTRUE};
+ Bool_t bChOnOff[10] = {kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE};
// Set degrees of freedom to align (see AliMUONAlignment)
fAlign->AllowVariations(bChOnOff);
+ // Set expected resolution (see AliMUONAlignment)
+ fAlign->SetSigmaXY(0.15,0.10);
+
// Fix parameters or add constraints here
// for (Int_t iSt=0; iSt<5; iSt++)
// if (!bStOnOff[iSt]) fAlign->FixStation(iSt+1);
for (Int_t iCh=0; iCh<10; iCh++)
if (!bChOnOff[iCh]) fAlign->FixChamber(iCh+1);
- // Left and right sides of the detector are independent, one can choose to align
+ // Left and right sides of the detector are independent, one can choose to align
// only one side
Bool_t bSpecLROnOff[2] = {kTRUE,kTRUE};
fAlign->FixHalfSpectrometer(bChOnOff,bSpecLROnOff);
fAlign->SetSpecLROnOff(bChOnOff);
// Here we can fix some detection elements
- fAlign->FixDetElem(908);
- fAlign->FixDetElem(1020);
+ // fAlign->FixDetElem(908);
+ // fAlign->FixDetElem(1012);
+ fAlign->FixDetElem(608);
// Set predifined global constrains: X, Y, P, XvsZ, YvsZ, PvsZ, XvsY, YvsY, PvsY
// Bool_t bVarXYT[9] = {kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE};
}
//________________________________________________________________________
-void AliMUONAlignmentTask::ConnectInputData(Option_t *)
+void AliMUONAlignmentTask::UserCreateOutputObjects()
{
- /// Connect ESD here. Called on each input data change.
-
- // Connect ESD here
- TTree* esdTree = dynamic_cast<TTree*> (GetInputData(0));
- if (!esdTree) {
- Printf("ERROR: Could not read chain from input slot 0");
- }
- else {
- AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
- if (!esdH) {
- Printf("ERROR: Could not get ESDInputHandler");
- }
- else {
- fESD = esdH->GetEvent();
- }
+
+ if( fDoAlignment )
+ {
+
+ // Creating graphs
+ fMSDEx = new TGraphErrors(156);
+ fMSDEy = new TGraphErrors(156);
+ fMSDEz = new TGraphErrors(156);
+ fMSDEp = new TGraphErrors(156);
+
+ // Add Ntuples to the list
+ fList = new TList();
+ fList->Add(fMSDEx);
+ fList->Add(fMSDEy);
+ fList->Add(fMSDEz);
+ fList->Add(fMSDEp);
+
+ fList->SetOwner(kTRUE);
+ PostData(1, fList);
+ }
+
+ // connect AOD output
+ if( fWriteRecords )
+ {
+ AliAODHandler* handler = dynamic_cast<AliAODHandler*>( AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler() );
+ if( handler )
+ {
+
+ // get AOD output handler and add Branch
+ fRecords = new TClonesArray( "AliMUONAlignmentTrackRecord", 0 );
+ fRecords->SetName( "records" );
+ handler->AddBranch("TClonesArray", &fRecords);
+
+ fRecordCount = 0;
+
+ } else AliInfo( "Error: invalid output event handler" );
+
}
-}
-//________________________________________________________________________
-void AliMUONAlignmentTask::CreateOutputObjects()
-{
- /// Executed once on each worker (machine actually running the analysis code)
- //
- // This method has to be called INSIDE the user redefined CreateOutputObjects
- // method, before creating each object corresponding to the output containers
- // that are to be written to a file. This need to be done in general for the big output
- // objects that may not fit memory during processing.
- // OpenFile(0);
-
- // Creating graphs
- fMSDEx = new TGraphErrors(156);
- fMSDEy = new TGraphErrors(156);
- fMSDEz = new TGraphErrors(156);
- fMSDEp = new TGraphErrors(156);
-
- // Add Ntuples to the list
- fList = new TList();
- fList->Add(fMSDEx);
- fList->Add(fMSDEy);
- fList->Add(fMSDEz);
- fList->Add(fMSDEp);
}
//________________________________________________________________________
-void AliMUONAlignmentTask::Exec(Option_t *)
+void AliMUONAlignmentTask::UserExec(Option_t *)
{
- /// Main loop, called for each event
- if (!fESD) {
- Printf("ERROR: fESD not available");
- return;
+ // print configuration
+// AliInfo( Form( "fReadRecords: %s", (fReadRecords ? "kTRUE":"kFALSE" ) ) );
+// AliInfo( Form( "fWriteRecords: %s", (fWriteRecords ? "kTRUE":"kFALSE" ) ) );
+// AliInfo( Form( "fDoAlignment: %s", (fDoAlignment ? "kTRUE":"kFALSE" ) ) );
+ // clear array
+ if( fWriteRecords && fRecords )
+ {
+ fRecords->Clear();
+ fRecordCount = 0;
}
+ // local track parameters
Double_t trackParams[8] = {0.,0.,0.,0.,0.,0.,0.,0.};
- if (fESD->GetRunNumber()!=fLastRunNumber){
- fLastRunNumber = fESD->GetRunNumber();
- Prepare(fGeoFilename.Data(),fDefaultOCDB.Data(),fMisAlignOCDB.Data());
- }
-
- Int_t nTracks = Int_t(fESD->GetNumberOfMuonTracks());
- // if (!event%100) cout << " there are " << nTracks << " tracks in event " << event << endl;
- for (Int_t iTrack = 0; iTrack < nTracks; iTrack++) {
- AliESDMuonTrack* esdTrack = fESD->GetMuonTrack(iTrack);
- if (!esdTrack->ClustersStored()) continue;
- Double_t invBenMom = esdTrack->GetInverseBendingMomentum();
-// fInvBenMom->Fill(invBenMom);
-// fBenMom->Fill(1./invBenMom);
- if (TMath::Abs(invBenMom)<=1.04) {
- AliMUONTrack track;
- AliMUONESDInterface::ESDToMUON(*esdTrack, track);
- fAlign->ProcessTrack(&track);
- fAlign->LocalFit(fTrackOk++,trackParams,0);
+
+ if( fReadRecords ) {
+
+ AliAODEvent* lAOD( dynamic_cast<AliAODEvent*>(InputEvent() ) );
+
+ // check AOD
+ if( !lAOD )
+ {
+ AliInfo("Error: AOD event not available");
+ return;
+ }
+
+ // read records
+ TClonesArray* records = static_cast<TClonesArray*>( lAOD->FindListObject( "records" ) );
+ if( !records )
+ {
+ AliInfo( "Unable to read object name \"records\"" );
+ return;
+ }
+
+ // get number of records and print
+ const Int_t lRecordsCount( records->GetEntriesFast() );
+ fTrackTot += lRecordsCount;
+
+ // loop over records
+ for( Int_t index = 0; index < lRecordsCount; ++index )
+ {
+
+ if( AliMUONAlignmentTrackRecord* record = dynamic_cast<AliMUONAlignmentTrackRecord*>( records->UncheckedAt(index) ) )
+ {
+
+ fAlign->ProcessTrack( record, fDoAlignment );
+ if( fDoAlignment )
+ { fAlign->LocalFit( fTrackOk++, trackParams, 0 ); }
+
+ } else AliInfo( Form( "Invalid record at %i", index ) );
+
+ }
+
+ } else {
+
+ /// Main loop, called for each event
+ AliESDEvent* lESD( dynamic_cast<AliESDEvent*>(InputEvent()) );
+
+ // check ESD
+ if( !lESD )
+ {
+ AliInfo("Error: ESD event not available");
+ return;
}
- fTrackTot++;
- cout << "Processed " << fTrackTot << " Tracks." << endl;
+
+// HUGO: Comment out check on run number, to be able to run on MC
+// if (!lESD->GetRunNumber())
+// {
+// AliInfo( Form( "Current Run Number: %i", lESD->GetRunNumber() ) );
+// return;
+// }
+
+ // if (lESD->GetRunNumber()!=fLastRunNumber){
+ // fLastRunNumber = lESD->GetRunNumber();
+ // Prepare(fGeoFilename.Data(),fDefaultOCDB.Data(),fMisAlignOCDB.Data());
+ // }
+
+ Int_t nTracks = Int_t(lESD->GetNumberOfMuonTracks());
+// cout << " there are " << nTracks << " tracks" << endl;
+ for( Int_t iTrack = 0; iTrack < nTracks; iTrack++ )
+ {
+
+ AliESDMuonTrack* esdTrack = lESD->GetMuonTrack(iTrack);
+ if (!esdTrack->ClustersStored()) continue;
+ if (!esdTrack->ContainTriggerData()) continue;
+
+ Double_t invBenMom = esdTrack->GetInverseBendingMomentum();
+ // fInvBenMom->Fill(invBenMom);
+ // fBenMom->Fill(1./invBenMom);
+ if (TMath::Abs(invBenMom)<=1.04)
+ {
+
+ AliMUONTrack track;
+ AliMUONESDInterface::ESDToMUON(*esdTrack, track);
+
+ // process track and retrieve corresponding records, for storage
+ const AliMUONAlignmentTrackRecord* lRecords( fAlign->ProcessTrack( &track, fDoAlignment ) );
+
+ // do the fit, if required
+ if( fDoAlignment ) fAlign->LocalFit(fTrackOk++,trackParams,0);
+ else fTrackOk++;
+
+ // store in array
+ if( fWriteRecords && fRecords )
+ {
+ new((*fRecords)[fRecordCount]) AliMUONAlignmentTrackRecord( *lRecords );
+ ++fRecordCount;
+ }
+
+ }
+
+ ++fTrackTot;
+
+ if(!(fTrackTot%1000))
+ { AliInfo( Form( "Processed %i Tracks and %i were fitted.", fTrackTot, fTrackOk ) ); }
+
+// cout << "Processed " << fTrackTot << " Tracks." << endl;
+ // Post final data. Write histo list to a file with option "RECREATE"
+ PostData(1,fList);
+
+ }
+
+ // save AOD
+ if( fWriteRecords && fRecordCount > 0 ) {
+ AliAODHandler* handler = dynamic_cast<AliAODHandler*>( AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler() );
+// printf("handler: %p\n",handler);
+ AliAODEvent* aod = handler->GetAOD();
+// printf("aod: %p\n",aod);
+ AliAODHeader* header = aod->GetHeader();
+// printf("header: %p\n",header);
+ header->SetRunNumber(lESD->GetRunNumber());
+// printf("RunNumber: %d\n",lESD->GetRunNumber());
+ AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()->SetFillAOD(kTRUE); }
+
}
-
- // Post final data. Write histo list to a file with option "RECREATE"
- PostData(0,fList);
-}
+
+}
+
+//________________________________________________________________________
+void AliMUONAlignmentTask::Terminate(const Option_t*)
+{ return; }
//________________________________________________________________________
void AliMUONAlignmentTask::FinishTaskOutput()
{
+
/// Called once per task on the client machine at the end of the analysis.
+ AliInfo( Form( "Processed %i tracks.", fTrackTot ) );
+ AliInfo( Form( "Accepted %i tracks.", fTrackOk ) );
+
+ // stop here if no alignment is to be performed
+ if( !fDoAlignment ) return;
+
+ AliLog::SetGlobalLogLevel(AliLog::kInfo);
- cout << "Processed " << fTrackTot << " Tracks." << endl;
// Perform global fit
fAlign->GlobalFit(fParameters,fErrors,fPulls);
- cout << "Done with GlobalFit " << endl;
-
// // Update pointers reading them from the output slot
// fList = (TList*)GetOutputData(0);
// fMSDEx = (TGraphErrors*)fList->At(0);
// fMSDEp = (TGraphErrors*)fList->At(3);
// Store results
- Double_t deId[156] = {0};
- Double_t msdEx[156] = {0};
- Double_t msdEy[156] = {0};
- Double_t msdEz[156] = {0};
- Double_t msdEp[156] = {0};
- Double_t deIdErr[156] = {0};
- Double_t msdExErr[156] = {0};
- Double_t msdEyErr[156] = {0};
- Double_t msdEzErr[156] = {0};
- Double_t msdEpErr[156] = {0};
+ Double_t DEid[156] = {0};
+ Double_t MSDEx[156] = {0};
+ Double_t MSDEy[156] = {0};
+ Double_t MSDEz[156] = {0};
+ Double_t MSDEp[156] = {0};
+ Double_t DEidErr[156] = {0};
+ Double_t MSDExErr[156] = {0};
+ Double_t MSDEyErr[156] = {0};
+ Double_t MSDEzErr[156] = {0};
+ Double_t MSDEpErr[156] = {0};
Int_t lNDetElem = 4*2+4*2+18*2+26*2+26*2;
Int_t lNDetElemCh[10] = {4,4,4,4,18,18,26,26,26,26};
// Int_t lSNDetElemCh[10] = {4,8,12,16,34,52,78,104,130,156};
Int_t idOffset = 0; // 400
Int_t lSDetElemCh = 0;
- for(Int_t iDE=0; iDE<lNDetElem; iDE++){
- deIdErr[iDE] = 0.;
- deId[iDE] = idOffset+100;
- deId[iDE] += iDE;
+ for(Int_t iDE=0; iDE<lNDetElem; iDE++)
+ {
+
+ DEidErr[iDE] = 0.;
+ DEid[iDE] = idOffset+100;
+ DEid[iDE] += iDE;
lSDetElemCh = 0;
- for(Int_t iCh=0; iCh<9; iCh++){
+ for(Int_t iCh=0; iCh<9; iCh++)
+ {
lSDetElemCh += lNDetElemCh[iCh];
- if (iDE>=lSDetElemCh) {
- deId[iDE] += 100;
- deId[iDE] -= lNDetElemCh[iCh];
+
+ if (iDE>=lSDetElemCh)
+ {
+ DEid[iDE] += 100;
+ DEid[iDE] -= lNDetElemCh[iCh];
}
+
}
- msdEx[iDE]=fParameters[3*iDE+0];
- msdEy[iDE]=fParameters[3*iDE+1];
- msdEz[iDE]=fParameters[3*iDE+3];
- msdEp[iDE]=fParameters[3*iDE+2];
- msdExErr[iDE]=(Double_t)fAlign->GetParError(3*iDE+0);
- msdEyErr[iDE]=(Double_t)fAlign->GetParError(3*iDE+1);
- msdEzErr[iDE]=(Double_t)fAlign->GetParError(3*iDE+3);
- msdEpErr[iDE]=(Double_t)fAlign->GetParError(3*iDE+2);
- fMSDEx->SetPoint(iDE,deId[iDE],fParameters[3*iDE+0]);
- fMSDEx->SetPoint(iDE,deIdErr[iDE],(Double_t)fAlign->GetParError(3*iDE+0));
- fMSDEy->SetPoint(iDE,deId[iDE],fParameters[3*iDE+1]);
- fMSDEy->SetPoint(iDE,deIdErr[iDE],(Double_t)fAlign->GetParError(3*iDE+1));
- fMSDEp->SetPoint(iDE,deId[iDE],fParameters[3*iDE+2]);
- fMSDEp->SetPoint(iDE,deIdErr[iDE],(Double_t)fAlign->GetParError(3*iDE+2));
- fMSDEz->SetPoint(iDE,deId[iDE],fParameters[3*iDE+3]);
- fMSDEz->SetPoint(iDE,deIdErr[iDE],(Double_t)fAlign->GetParError(3*iDE+3));
+
+ MSDEx[iDE]=fParameters[4*iDE+0];
+ MSDEy[iDE]=fParameters[4*iDE+1];
+ MSDEz[iDE]=fParameters[4*iDE+3];
+ MSDEp[iDE]=fParameters[4*iDE+2];
+ MSDExErr[iDE]=(Double_t)fAlign->GetParError(4*iDE+0);
+ MSDEyErr[iDE]=(Double_t)fAlign->GetParError(4*iDE+1);
+ MSDEzErr[iDE]=(Double_t)fAlign->GetParError(4*iDE+3);
+ MSDEpErr[iDE]=(Double_t)fAlign->GetParError(4*iDE+2);
+ fMSDEx->SetPoint(iDE,DEid[iDE],fParameters[4*iDE+0]);
+ fMSDEx->SetPointError(iDE,DEidErr[iDE],(Double_t)fAlign->GetParError(4*iDE+0));
+ fMSDEy->SetPoint(iDE,DEid[iDE],fParameters[4*iDE+1]);
+ fMSDEy->SetPointError(iDE,DEidErr[iDE],(Double_t)fAlign->GetParError(4*iDE+1));
+ fMSDEz->SetPoint(iDE,DEid[iDE],fParameters[4*iDE+3]);
+ fMSDEz->SetPointError(iDE,DEidErr[iDE],(Double_t)fAlign->GetParError(4*iDE+3));
+ fMSDEp->SetPoint(iDE,DEid[iDE],fParameters[4*iDE+2]);
+ fMSDEp->SetPointError(iDE,DEidErr[iDE],(Double_t)fAlign->GetParError(4*iDE+2));
+
}
// Post final data. Write histo list to a file with option "RECREATE"
- PostData(0,fList);
+ PostData(1,fList);
+
+ // HUGO: stop here to test reproducibility
+ // return;
// Re Align
- AliMUONGeometryTransformer *newTransform = fAlign->ReAlign(fTransform,fParameters,true);
- newTransform->WriteTransformations("transform2ReAlign.dat");
-
+ fNewGeoTransformer = fAlign->ReAlign(fOldGeoTransformer,fParameters,true);
+ // newTransform->WriteTransformations("transform2ReAlign.dat");
+
// Generate realigned data in local cdb
- const TClonesArray* array = newTransform->GetMisAlignmentData();
+ const TClonesArray* array = fNewGeoTransformer->GetMisAlignmentData();
// 100 mum residual resolution for chamber misalignments?
fAlign->SetAlignmentResolution(array,-1,0.01,0.01,0.004,0.003);
-
+
// CDB manager
- AliCDBManager* cdbManager = AliCDBManager::Instance();
- cdbManager->SetDefaultStorage(fDefaultOCDB.Data());
- cdbManager->SetSpecificStorage("MUON/Align/Data",fMisAlignOCDB.Data());
-
+ AliLog::SetGlobalDebugLevel(2);
+ AliCDBManager* cdbm = AliCDBManager::Instance();
+ // cdbManager->SetDefaultStorage(fDefaultOCDB.Data());
+
+ // recover default storage full name (raw:// cannot be used to set specific storage)
+ TString defaultStorage(cdbm->GetDefaultStorage()->GetType());
+ if (defaultStorage == "alien") defaultStorage += Form("://folder=%s", cdbm->GetDefaultStorage()->GetBaseFolder().Data());
+ else defaultStorage += Form("://%s", cdbm->GetDefaultStorage()->GetBaseFolder().Data());
+
+ if (fOldAlignStorage != "none") cdbm->UnloadFromCache("MUON/Align/Data");
+ if (!fNewAlignStorage.IsNull()) cdbm->SetSpecificStorage("MUON/Align/Data",fNewAlignStorage.Data());
+ // else cdbm->SetSpecificStorage("MUON/Align/Data",fDefaultStorage.Data());
+ else cdbm->SetSpecificStorage("MUON/Align/Data",defaultStorage.Data());
+
+
AliCDBMetaData* cdbData = new AliCDBMetaData();
cdbData->SetResponsible("Dimuon Offline project");
cdbData->SetComment("MUON alignment objects with residual misalignment");
- AliCDBId id("MUON/Align/Data", 0, AliCDBRunRange::Infinity());
- cdbManager->Put(const_cast<TClonesArray*>(array), id, cdbData);
+ AliCDBId id("MUON/Align/Data", 0, AliCDBRunRange::Infinity());
+ cdbm->Put(const_cast<TClonesArray*>(array), id, cdbData);
+
+ gSystem->Exec("cp -a ReAlignOCDB/MUON/Align/Data/Run0_999999999_v0_s0.root Run0_999999999_v0_s0.root");
+ gSystem->Exec("ls -l");
}
//________________________________________________________________________
-void AliMUONAlignmentTask::Terminate(const Option_t*)
+void AliMUONAlignmentTask::NotifyRun()
{
- /// Called once per task on the client machine at the end of the analysis.
+
+ if (fOCDBLoaded && fLoadOCDBOnce) { AliError(Form("OCDB already loaded %d %d",fOCDBLoaded,fLoadOCDBOnce));
+ return;
+ }
-}
+ AliCDBManager* cdbm = AliCDBManager::Instance();
+ cdbm->SetDefaultStorage(fDefaultStorage.Data());
+// cdbm->SetDefaultStorage("local://$ALICE_ROOT/OCDB");
-//-----------------------------------------------------------------------
-void AliMUONAlignmentTask::Prepare(const char* geoFilename, const char* defaultOCDB, const char* misAlignOCDB)
-{
- /// Set the geometry, the magnetic field, the mapping and the reconstruction parameters
-
- // Load mapping
- AliCDBManager* man = AliCDBManager::Instance();
- man->SetDefaultStorage(defaultOCDB);
- man->SetSpecificStorage("MUON/Align/Data",misAlignOCDB);
- man->Print();
- man->SetRun(fLastRunNumber);
- if ( ! AliMpCDB::LoadDDLStore() ) {
- Error("MUONRefit","Could not access mapping from OCDB !");
- exit(-1);
- }
+ // HUGO: add to comment out the specific settings below in order to be able to run.
+ //cdbm->SetSpecificStorage("GRP/GRP/Data",fDefaultStorage.Data());
+ //cdbm->SetSpecificStorage("GRP/Geometry/data",fDefaultStorage.Data());
+ //cdbm->SetSpecificStorage("MUON/Calib/MappingData",fDefaultStorage.Data());
+ cdbm->SetRun(fCurrentRunNumber);
- // Import TGeo geometry (needed by AliMUONTrackExtrap::ExtrapToVertex)
- if (!gGeoManager) {
- AliGeomManager::LoadGeometry(geoFilename);
- if (!gGeoManager) {
- Error("AliMUONReAlignTask", "getting geometry from file %s failed", "generated/galice.root");
- return;
- }
- }
+ if (!AliMUONCDB::LoadField()) { AliError("Problem field"); return;}
- // set mag field
- if (!TGeoGlobalMagField::Instance()->GetField()) {
- printf("Loading field map...\n");
- AliGRPManager *grpMan = new AliGRPManager();
- grpMan->ReadGRPEntry();
- grpMan->SetMagField();
- delete grpMan;
- }
// set the magnetic field for track extrapolations
AliMUONTrackExtrap::SetField();
-
+
+ if (!AliMUONCDB::LoadMapping(kTRUE)) { AliError("Problem mapping"); return;}
+
+ // recover default storage full name (raw:// cannot be used to set specific storage)
+ TString defaultStorage(cdbm->GetDefaultStorage()->GetType());
+ if (defaultStorage == "alien") defaultStorage += Form("://folder=%s", cdbm->GetDefaultStorage()->GetBaseFolder().Data());
+ else defaultStorage += Form("://%s", cdbm->GetDefaultStorage()->GetBaseFolder().Data());
+
+ // reset existing geometry/alignment if any
+ if (cdbm->GetEntryCache()->Contains("GRP/Geometry/Data")) cdbm->UnloadFromCache("GRP/Geometry/Data");
+ if (cdbm->GetEntryCache()->Contains("MUON/Align/Data")) cdbm->UnloadFromCache("MUON/Align/Data");
+ if (AliGeomManager::GetGeometry()) AliGeomManager::GetGeometry()->UnlockGeometry();
+
+ // get original geometry transformer
+ AliGeomManager::LoadGeometry();
+ if (!AliGeomManager::GetGeometry()) { AliError("Problem geometry"); return;}
+ if (fOldAlignStorage != "none")
+ {
+
+ if (!fOldAlignStorage.IsNull()) cdbm->SetSpecificStorage("MUON/Align/Data",fOldAlignStorage.Data());
+ else cdbm->SetSpecificStorage("MUON/Align/Data",defaultStorage.Data());
+ // else cdbm->SetSpecificStorage("MUON/Align/Data",fDefaultStorage.Data());
+
+ AliGeomManager::ApplyAlignObjsFromCDB("MUON");
+
+ }
+
+ // fOldGeoTransformer = new AliMUONGeometryTransformer();
+ fOldGeoTransformer->LoadGeometryData();
+ fAlign->SetGeometryTransformer(fOldGeoTransformer);
+
+ fOCDBLoaded = kTRUE;
+
}