+/**************************************************************************
+ * 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. *
+ **************************************************************************/
+
+//
+// This implementation AliTPCExB is using an aproximate calculation of the ExB
+// effect. Therefore the drift ODE is Taylor expanded and only the first
+// order part is taken.
+//
+//
+// The ExB correction map is stored in the calib DB
+// To test current version:
+/*
+ char *storage = "local://OCDBres"
+ Int_t RunNumber=0;
+ AliCDBManager::Instance()->SetDefaultStorage(storage);
+ AliCDBManager::Instance()->SetRun(RunNumber)
+ AliTPCExBFirst * exb = AliTPCcalibDB::Instance()->GetExB();
+
+ // exb->TestExB("exb.root");
+ // TFile f("exb.root");
+ //positions->Draw("drphi");
+*/
+
+
+
#include "TMath.h"
+//#include "AliFieldMap.h"
+#include "AliMagF.h"
#include "TTreeStream.h"
#include "AliTPCExBFirst.h"
ClassImp(AliTPCExBFirst)
const Double_t AliTPCExBFirst::fgkEM=1.602176487e-19/9.10938215e-31;
-const Double_t AliTPCExBFirst::fgkDriftField=40.e3;
+const Double_t AliTPCExBFirst::fgkDriftField=-40.e3;
AliTPCExBFirst::AliTPCExBFirst()
: fDriftVelocity(0),
//
// purely for I/O
//
+ SetInstance(this);
}
AliTPCExBFirst::AliTPCExBFirst(const AliMagF *bField,
// drift velocity. Since some kind of lookuptable is created the
// number of its meshpoints can be supplied.
//
- ConstructCommon(0,bField);
+ // ConstructCommon(0,bField);
+ ConstructCommon(bField);
+ SetInstance(this);
}
+/*
AliTPCExBFirst::AliTPCExBFirst(const AliFieldMap *bFieldMap,
Double_t driftVelocity)
: fDriftVelocity(driftVelocity),
// The constructor. One has to supply a field map and an (initial)
// drift velocity.
//
-
+ SetInstance(this);
fkXMin=bFieldMap->Xmin()
-TMath::Ceil( (bFieldMap->Xmin()+250.0)/bFieldMap->DelX())
*bFieldMap->DelX();
ConstructCommon(bFieldMap,0);
}
+*/
AliTPCExBFirst::~AliTPCExBFirst() {
//
//
// correct for the distortion
//
- Double_t Bx,By;
- GetMeanFields(position[0],position[1],position[2],&Bx,&By);
+ Double_t bx,by;
+ GetMeanFields(position[0],position[1],position[2],&bx,&by);
if (position[2]>0.) {
- Double_t Bxe,Bye;
- GetMeanFields(position[0],position[1],250.,&Bxe,&Bye);
+ Double_t bxe,bye;
+ GetMeanFields(position[0],position[1],250.,&bxe,&bye);
if (position[2]!=250.) {
- Bx=(500.*Bxe-(position[2]+250.)*Bx)/(250.-position[2]);
- By=(500.*Bye-(position[2]+250.)*By)/(250.-position[2]);
+ bx=(500.*bxe-(position[2]+250.)*bx)/(250.-position[2]);
+ by=(500.*bye-(position[2]+250.)*by)/(250.-position[2]);
}
else {
- Bx=Bxe;
- By=Bye;
+ bx=bxe;
+ by=bye;
}
}
Double_t mu=fDriftVelocity/fgkDriftField;
Double_t wt=mu*fkMeanBz;
- corrected[0]=mu*(wt*Bx-By)/(1.+wt*wt);
- corrected[1]=mu*(wt*By+Bx)/(1.+wt*wt);
+ corrected[0]=mu*(wt*bx-by)/(1.+wt*wt);
+ corrected[1]=mu*(wt*by+bx)/(1.+wt*wt);
if (position[2]>0.) {
corrected[0]*=(250.-position[2]);
}
}
-void AliTPCExBFirst::ConstructCommon(const AliFieldMap *bFieldMap,
+
+void AliTPCExBFirst::ConstructCommon(//const AliFieldMap *bFieldMap,
const AliMagF *bField) {
//
// THIS IS PRIVATE! (a helper for the constructor)
//
fkNMean=fkNX*fkNY*fkNZ;
- fkMeanBx=new Float_t[fkNMean];
- fkMeanBy=new Float_t[fkNMean];
+ fkMeanBx=new Double_t[fkNMean];
+ fkMeanBy=new Double_t[fkNMean];
Double_t x[3];
Double_t nBz=0;
x[0]=fkXMin+i*(fkXMax-fkXMin)/(fkNX-1);
for (int j=0;j<fkNY;++j) {
x[1]=fkYMin+j*(fkYMax-fkYMin)/(fkNY-1);
- Double_t R=TMath::Sqrt(x[0]*x[0]+x[1]*x[1]);
- Double_t Bx=0.,By=0.;
+ Double_t r=TMath::Sqrt(x[0]*x[0]+x[1]*x[1]);
+ Double_t bx=0.,by=0.;
for (int k=0;k<fkNZ;++k) {
x[2]=fkZMin+k*(fkZMax-fkZMin)/(fkNZ-1);
- Float_t B[3];
+ Double_t b[3];
// the x is not const in the Field function...
- Float_t xt[3];
+ Double_t xt[3];
for (int l=0;l<3;++l) xt[l]=x[l];
// that happens due to the lack of a sophisticated class design:
- if (bFieldMap!=0)
- bFieldMap->Field(xt,B);
- else
- bField->Field(xt,B);
- Bx+=B[0]/10.;
- By+=B[1]/10.;
- /*old
- fkMeanBx[(k*fkNY+j)*fkNX+i]=Bx/(k+1);
- fkMeanBy[(k*fkNY+j)*fkNX+i]=By/(k+1);
- */
- fkMeanBx[(k*fkNY+j)*fkNX+i]=Bx;
- fkMeanBy[(k*fkNY+j)*fkNX+i]=By;
- if (90.<=R&&R<=250.) {
- fkMeanBz+=B[2]/10.;
+ // if (bFieldMap!=0)
+ // bFieldMap->Field(xt,b);
+ // else
+ ((AliMagF*)bField)->Field(xt,b);
+ bx+=b[0]/10.;
+ by+=b[1]/10.;
+ fkMeanBx[(k*fkNY+j)*fkNX+i]=bx;
+ fkMeanBy[(k*fkNY+j)*fkNX+i]=by;
+ if (90.<=r&&r<=250.) {
+ fkMeanBz+=b[2]/10.;
++nBz;
}
}
fkMeanBz/=nBz;
}
+
void AliTPCExBFirst::GetMeanFields(Double_t rx,Double_t ry,Double_t rz,
Double_t *Bx,Double_t *By) const {
//
double s0x=fkMeanBx[yi1*fkNX+xi1]*dx1*dy1
+fkMeanBx[yi2*fkNX+xi1]*dx1*dy
- +fkMeanBx[yi1*fkNX+xi2]*dx *dy
- +fkMeanBx[yi2*fkNX+xi2]*dx *dy1;
+ +fkMeanBx[yi1*fkNX+xi2]*dx *dy1
+ +fkMeanBx[yi2*fkNX+xi2]*dx *dy;
double s0y=fkMeanBy[yi1*fkNX+xi1]*dx1*dy1
+fkMeanBy[yi2*fkNX+xi1]*dx1*dy
- +fkMeanBy[yi1*fkNX+xi2]*dx *dy
- +fkMeanBy[yi2*fkNX+xi2]*dx *dy1;
+ +fkMeanBy[yi1*fkNX+xi2]*dx *dy1
+ +fkMeanBy[yi2*fkNX+xi2]*dx *dy;
Int_t zi0=zi1-1;
double snmx,snmy;
if (zi0>=0) {
snmx=fkMeanBx[(zi0*fkNY+yi1)*fkNX+xi1]*dx1*dy1
+fkMeanBx[(zi0*fkNY+yi2)*fkNX+xi1]*dx1*dy
- +fkMeanBx[(zi0*fkNY+yi1)*fkNX+xi2]*dx *dy
- +fkMeanBx[(zi0*fkNY+yi2)*fkNX+xi2]*dx *dy1;
+ +fkMeanBx[(zi0*fkNY+yi1)*fkNX+xi2]*dx *dy1
+ +fkMeanBx[(zi0*fkNY+yi2)*fkNX+xi2]*dx *dy;
snmy=fkMeanBy[(zi0*fkNY+yi1)*fkNX+xi1]*dx1*dy1
+fkMeanBy[(zi0*fkNY+yi2)*fkNX+xi1]*dx1*dy
- +fkMeanBy[(zi0*fkNY+yi1)*fkNX+xi2]*dx *dy
- +fkMeanBy[(zi0*fkNY+yi2)*fkNX+xi2]*dx *dy1;
+ +fkMeanBy[(zi0*fkNY+yi1)*fkNX+xi2]*dx *dy1
+ +fkMeanBy[(zi0*fkNY+yi2)*fkNX+xi2]*dx *dy;
}
else
snmx=snmy=0.;
double snx=fkMeanBx[(zi1*fkNY+yi1)*fkNX+xi1]*dx1*dy1
+fkMeanBx[(zi1*fkNY+yi2)*fkNX+xi1]*dx1*dy
- +fkMeanBx[(zi1*fkNY+yi1)*fkNX+xi2]*dx *dy
- +fkMeanBx[(zi1*fkNY+yi2)*fkNX+xi2]*dx *dy1;
+ +fkMeanBx[(zi1*fkNY+yi1)*fkNX+xi2]*dx *dy1
+ +fkMeanBx[(zi1*fkNY+yi2)*fkNX+xi2]*dx *dy;
double sny=fkMeanBy[(zi1*fkNY+yi1)*fkNX+xi1]*dx1*dy1
+fkMeanBy[(zi1*fkNY+yi2)*fkNX+xi1]*dx1*dy
- +fkMeanBy[(zi1*fkNY+yi1)*fkNX+xi2]*dx *dy
- +fkMeanBy[(zi1*fkNY+yi2)*fkNX+xi2]*dx *dy1;
+ +fkMeanBy[(zi1*fkNY+yi1)*fkNX+xi2]*dx *dy1
+ +fkMeanBy[(zi1*fkNY+yi2)*fkNX+xi2]*dx *dy;
double snpx=fkMeanBx[(zi2*fkNY+yi1)*fkNX+xi1]*dx1*dy1
+fkMeanBx[(zi2*fkNY+yi2)*fkNX+xi1]*dx1*dy
- +fkMeanBx[(zi2*fkNY+yi1)*fkNX+xi2]*dx *dy
- +fkMeanBx[(zi2*fkNY+yi2)*fkNX+xi2]*dx *dy1;
+ +fkMeanBx[(zi2*fkNY+yi1)*fkNX+xi2]*dx *dy1
+ +fkMeanBx[(zi2*fkNY+yi2)*fkNX+xi2]*dx *dy;
double snpy=fkMeanBy[(zi2*fkNY+yi1)*fkNX+xi1]*dx1*dy1
+fkMeanBy[(zi2*fkNY+yi2)*fkNX+xi1]*dx1*dy
- +fkMeanBy[(zi2*fkNY+yi1)*fkNX+xi2]*dx *dy
- +fkMeanBy[(zi2*fkNY+yi2)*fkNX+xi2]*dx *dy1;
+ +fkMeanBy[(zi2*fkNY+yi1)*fkNX+xi2]*dx *dy1
+ +fkMeanBy[(zi2*fkNY+yi2)*fkNX+xi2]*dx *dy;
+
+
+
*Bx=0.5*(((snpx-2.*snx+snmx)*dz+2.*(snx-snmx))*dz+snx-s0x+snmx);
*By=0.5*(((snpy-2.*sny+snmy)*dz+2.*(sny-snmy))*dz+sny-s0y+snmy);
//TODO: make this nice
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff