--- /dev/null
+//Macros to calculate PHOS alignment matrixes from survay
+//To run select module number (mod=1,2,3,4) and uncomment correcponding
+//pairs of matrixes: survay points in crystal plane and in ALICE global system
+//These data one takes from syrvay results.
+// Note that for Mod.4 there was now survay repers in crystal surface frame.
+// So we choose module with closest ditances between repers (mod2==mod4)
+// and assume that mod4 has similr reper coordinates wrt crytal surface.
+//
+
+TGeoHMatrix *mPHOS[6] ;
+Int_t mod=1 ; //Which module minimize
+
+
+//Module 1, points 500x survey
+ const Float_t srv5000[3]={-217.61,-386.44,-77.98} ;
+ const Float_t srv5001[3]={ -79.32,-439.05,-77.54} ;
+ const Float_t srv5002[3]={ -79.40,-439.29, 80.80} ;
+ const Float_t srv5003[3]={-217.64,-386.49, 80.04} ;
+
+//Module 1, points 500x wrt crystals (copied from Mod3)
+ const Float_t p5000[3]={ -3.08,28.97,-16.43} ;
+ const Float_t p5001[3]={145.94,28.90,-16.03} ;
+ const Float_t p5002[3]={145.58,29.01,142.17} ;
+ const Float_t p5003[3]={ -3.58,28.88,141.88} ;
+
+
+
+
+
+/*
+//Module 2, points 500x survey
+ const Float_t srv5000[3]={-74.67,-429.97,-90.95} ;
+ const Float_t srv5001[3]={ 72.39,-430.96,-91.25} ;
+ const Float_t srv5002[3]={ 64.43,-431.06, 92.01} ;
+ const Float_t srv5003[3]={-70.56,-430.27, 92.43} ;
+*/
+/*
+//Module 3, points 500x survey
+ const Float_t srv5000[3]={ 77.95,-430.77,-78.85} ;
+ const Float_t srv5001[3]={217.78,-380.47,-78.76} ;
+ const Float_t srv5002[3]={217.86,-380.25, 79.28} ;
+ const Float_t srv5003[3]={ 77.57,-430.93, 79.58} ;
+*/
+/*
+//Module 4, points 500x survey
+ const Float_t srv5000[3]={222.58,-377.81,-78.96} ;
+ const Float_t srv5001[3]={333.99,-284.63,-78.21} ;
+ const Float_t srv5002[3]={335.38,-283.62, 80.39} ;
+ const Float_t srv5003[3]={220.94,-379.38, 79.89} ;
+*/
+
+
+/*
+ //Module 2 reference points on crystal surface (in cm)
+ const Float_t p5000[3]={0.0000,0.0000, 0.00} ;
+ const Float_t p5001[3]={0.0000,0.0000,122.62} ;
+ const Float_t p5002[3]={142.42,0.0000, 1.95} ;
+ const Float_t p5003[3]={142.49,0.0001,124.69} ;
+
+ //Module 2 survey: position of these points in global ALICE system
+ const Float_t srv5000[3]={-71.06,-459.00,-61.82} ;
+ const Float_t srv5001[3]={-71.16,-459.15, 61.06} ;
+ const Float_t srv5002[3]={ 71.65,-459.99,-59.76} ;
+ const Float_t srv5003[3]={ 71.63,-460.13, 63.24} ;
+*/
+/*
+ //Module 3 reference points on crystal surface
+ const Float_t p5000[3]={0.0000,0.0000,0.0000 } ;
+ const Float_t p5001[3]={0.0000,0.0000,122.63 } ;
+ const Float_t p5002[3]={142.53,0.0000, 2.37 } ;
+ const Float_t p5003[3]={142.28, -0.11,125.06 } ;
+
+ //Module 3 survey
+ const Float_t srv5000[3]={ 90.51,-457.01,-62.35 } ;
+ const Float_t srv5001[3]={ 90.77,-456.86, 60.19 } ;
+ const Float_t srv5002[3]={224.51,-408.73,-60.33 } ;
+ const Float_t srv5003[3]={224.57,-408.77, 62.27 } ;
+
+*/
+/*
+ //Module 4 reference points on crystal surface
+ const Float_t p5000[3]={0.0000,0.0000,0.0000 } ;
+ const Float_t p5001[3]={0.0000,0.0000,122.98 } ;
+ const Float_t p5002[3]={142.39,0.0000, 2.05 } ;
+ const Float_t p5003[3]={142.41,0.0019,124.82 } ;
+
+ //Module 4 survey
+ const Float_t srv5000[3]={244.57,-399.74,-62.67 } ;
+ const Float_t srv5001[3]={244.10,-400.79, 60.38 } ;
+ const Float_t srv5002[3]={353.95,-308.47,-59.42 } ;
+ const Float_t srv5003[3]={353.38,-309.36, 63.42 } ;
+*/
+
+
+
+
+
+Double_t DistV(TVector3 &a, Float_t * b){
+
+ return sqrt((a[0]-b[0])*(a[0]-b[0])+ (a[1]-b[1])*(a[1]-b[1])+ (a[2]-b[2])*(a[2]-b[2]) ) ;
+
+}
+
+Double_t Dist(Float_t *a, Float_t *b){
+
+ return sqrt((a[0]-b[0])*(a[0]-b[0])+ (a[1]-b[1])*(a[1]-b[1])+ (a[2]-b[2])*(a[2]-b[2]) ) ;
+
+}
+void fcn(Int_t &npar, Double_t *gin, Double_t &f, Double_t *par, Int_t iflag)
+{
+ //Function used in Minuit minimization:
+ //Minimize distances between measred reper coordinates
+ //and calculate after PHOS rotaion/shifts.
+ //Uses dedicated function AliPHOSGeoUtils::TestSurvey()
+
+ Double_t chisq = 0;
+
+ //Create current rotation matrix
+ TGeoRotation r1;
+ r1.SetAngles(par[0],par[1],par[2]) ;
+ mPHOS[mod]->SetRotation(r1.GetRotationMatrix()) ;
+ mPHOS[mod]->SetDx(par[3]) ;
+ mPHOS[mod]->SetDy(par[4]) ;
+ mPHOS[mod]->SetDz(par[5]) ;
+
+ AliPHOSGeoUtils * geom = new AliPHOSGeoUtils("PHOS") ;
+ geom->SetMisalMatrix(mPHOS[mod],mod-1) ;
+
+ TVector3 globaPos ;
+ //Compare ditances for all 4 points
+ geom->TestSurvey(mod,p5000,globaPos) ;
+ chisq+=DistV(globaPos,srv5000) ;
+
+ geom->TestSurvey(mod,p5001,globaPos) ;
+ chisq+=DistV(globaPos,srv5001) ;
+
+ geom->TestSurvey(mod,p5002,globaPos) ;
+ chisq+=DistV(globaPos,srv5002) ;
+
+ geom->TestSurvey(mod,p5003,globaPos) ;
+ chisq+=DistV(globaPos,srv5003) ;
+
+ delete geom ;
+ f = chisq;
+}
+
+
+void CheckPHOSAlignment(){
+ //Checks if PHOS alignment agree with results of photogrammetry
+
+
+ mPHOS[mod] = new TGeoHMatrix() ;
+
+ //Initial conditions for minimization, close to ideal PHOS positions
+ TGeoRotation r1;
+ r1.SetAngles(0,90+(mod-2)*20,0) ;
+ r1.Print() ;
+
+ mPHOS[mod]->SetRotation(r1.GetRotationMatrix()) ;
+ mPHOS[mod]->SetDx(-480.*TMath::Cos((210+mod*20)*TMath::DegToRad())) ;
+ mPHOS[mod]->SetDy(-480.*TMath::Sin((210+mod*20)*TMath::DegToRad())) ;
+ mPHOS[mod]->SetDz(0.) ;
+
+
+
+ //Each module has 4 reference points
+ //Define them w.r.t. center of top right crystal (see more details in EDH comenents
+ // 992344, 1002387, 1012391.
+
+ //Print distances between repers in current module
+//printf("%6.2f,\t %6.2f,\t %6.2f,\t %6.2f\n",Dist(srv5000,srv5001),Dist(srv5000,srv5002),Dist(srv5000,srv5003),Dist(srv5001,srv5003)) ;
+//return ;
+
+
+
+ //Module 3
+ const Float_t m3p5000[3]={-30.76, 191.15, -164.32} ;
+ const Float_t m3p5001[3]={1459.44, 190.47, -160.25} ;
+ const Float_t m3p5002[3]={1455.82, 191.57, 1421.70} ;
+ const Float_t m3p5003[3]={-35.83, 190.32, 1418.78} ;
+
+ //Module 3 survay (in meters!)
+ const Float_t m3srv5000[3]={0.7795,-4.3077,-0.7885} ;
+ const Float_t m3srv5001[3]={2.1778,-3.8047,-0.7876} ;
+ const Float_t m3srv5002[3]={2.1786,-3.8025,0.7928} ;
+ const Float_t m3srv5003[3]={0.7757,-4.3093,0.7958} ;
+
+
+ //Module 4
+ const Float_t m4p5000[3]={-34.59, 210.66, -165.23} ;
+ const Float_t m4p5001[3]={1456.73, 210.71, -164.91} ;
+ const Float_t m4p5002[3]={1454.76, 213.62, 1418.01} ;
+ const Float_t m4p5003[3]={-36.60, 210.14, 1419.84} ;
+
+
+ const Float_t m4srv5000[3]={2.2258,-3.7781,-0.7896} ;
+ const Float_t m4srv5001[3]={3.3399,-2.8463,-0.7821} ;
+ const Float_t m4srv5002[3]={3.3538,-2.8362, 0.8039} ;
+ const Float_t m4srv5003[3]={2.2094,-3.7938, 0.7989} ;
+
+
+
+
+ // prepare minuit
+ const Double_t nPar = 6; // number of fit parameters
+ TMinuit* gMinuit = new TMinuit(nPar);
+ gMinuit->SetFCN(fcn);
+ gMinuit->SetPrintLevel(1); // -1 quiet, 0 normal, 1 verbose
+
+ Double_t arglist[10];
+ Int_t ierflg = 0;
+
+ // start minimization
+ Double_t chi2 = 0;
+
+ ierflg = 0;
+ arglist[0] = 1;
+ gMinuit->mnexcm("SET ERR", arglist, 1, ierflg);
+ gMinuit->SetMaxIterations(5000);
+
+ // parameters: parameter no., name, start value, step size, parameter range min., parameter range max., status
+ gMinuit->mnparm(0,"A1",0., 1., 0, 0, ierflg);
+ gMinuit->mnparm(1,"A2",90., 1., 0, 0, ierflg);
+ gMinuit->mnparm(2,"A3",0., 1., 0, 0, ierflg);
+ gMinuit->mnparm(3,"dx",0., 1., 0, 0, ierflg);
+ gMinuit->mnparm(4,"dy",-460., 1., 0, 0, ierflg);
+ gMinuit->mnparm(5,"dz",0., 1., 0, 0, ierflg);
+
+ // Now ready for minimization step
+ arglist[0] = 5000;
+ arglist[1] = 1.;
+ gMinuit->mnexcm("MIGRAD", arglist, 2, ierflg);
+
+ AliPHOSGeometry * geom = AliPHOSGeometry::GetInstance("PHOS_mod1234") ;
+ geom->SetMisalMatrix(mPHOS[mod],mod-1) ;
+
+
+
+ TVector3 globaPos ;
+
+
+ geom->TestSurvey(mod,p5000,globaPos) ;
+ printf("P5000: (%f, %f, %f) \n",globaPos.X(),globaPos.Y(),globaPos.Z()) ;
+ printf(" srv: (%f, %f, %f) \n",srv5000[0],srv5000[1],srv5000[2]) ;
+ geom->TestSurvey(mod,p5001,globaPos) ;
+ printf("P5001: (%f, %f, %f) \n",globaPos.X(),globaPos.Y(),globaPos.Z()) ;
+ printf(" srv: (%f, %f, %f) \n",srv5001[0],srv5001[1],srv5001[2]) ;
+ geom->TestSurvey(mod,p5002,globaPos) ;
+ printf("P5002: (%f, %f, %f) \n",globaPos.X(),globaPos.Y(),globaPos.Z()) ;
+ printf(" srv: (%f, %f, %f) \n",srv5002[0],srv5002[1],srv5002[2]) ;
+ geom->TestSurvey(mod,p5003,globaPos) ;
+ printf("P5003: (%f, %f, %f) \n",globaPos.X(),globaPos.Y(),globaPos.Z()) ;
+ printf(" srv: (%f, %f, %f) \n",srv5003[0],srv5003[1],srv5003[2]) ;
+
+ mPHOS[mod]->Print() ;
+}
+
+
+
--- /dev/null
+/*
+**********No misalignment*******
+i=0
+matrix global_1 - tr=1 rot=1 refl=0 scl=0
+ 0.766044 0.000000 0.642788 Tx = 310.215729
+ 0.642788 0.000000 -0.766044 Ty = -369.700684
+ 0.000000 1.000000 0.000000 Tz = 0.000000
+i=1
+matrix global_1 - tr=1 rot=1 refl=0 scl=0
+ 0.939693 0.000000 0.342020 Tx = 165.062332
+ 0.342020 0.000000 -0.939693 Ty = -453.505035
+ 0.000000 1.000000 0.000000 Tz = 0.000000
+i=2
+matrix global_1 - tr=1 rot=1 refl=0 scl=0
+ 1.000000 0.000000 0.000000 Tx = 0.000000
+ 0.000000 0.000000 -1.000000 Ty = -482.609985
+ 0.000000 1.000000 0.000000 Tz = 0.000000
+i=3
+matrix global_1 - tr=1 rot=1 refl=0 scl=0
+ 0.939693 0.000000 -0.342020 Tx = -165.062332
+ -0.342020 0.000000 -0.939693 Ty = -453.505035
+ 0.000000 1.000000 0.000000 Tz = 0.000000
+
+*************
+=========Should be========
+
+mod 4 (i=0)
+matrix - tr=1 rot=1 refl=0 scl=0
+ 0.767337 -0.004312 0.641229 Tx = 313.381746
+ 0.641189 -0.007942 -0.767342 Ty = -372.024764
+ 0.008402 0.999959 -0.003330 Tz = -2.354855
+
+mod 3 (i=1)
+matrix - tr=1 rot=1 refl=0 scl=0
+ 0.940798 0.002136 0.338961 Tx = 165.078518
+ 0.338959 0.001227 -0.940800 Ty = -454.148853
+ -0.002426 0.999997 0.000431 Tz = -2.801675
+
+mod 3 (i=2)
+matrix - tr=1 rot=1 refl=0 scl=0
+ 0.999976 -0.000756 -0.006876 Tx = -0.083626
+ -0.006877 -0.001172 -0.999976 Ty = -482.172358
+ 0.000748 0.999999 -0.001178 Tz = -1.906440
+
+mod 3 (i=3)
+matrix - tr=1 rot=1 refl=0 scl=0
+ 0.934365 -0.002092 -0.356311 Tx = -167.233337
+ -0.356313 -0.000801 -0.934366 Ty = -460.901723
+ 0.001670 0.999997 -0.001494 Tz = -1.821245
+
+------------Correct angles, but zero offsets------
+ 0.767338 -0.004312 0.641229 Tx = 309.463421
+ 0.641188 -0.007942 -0.767343 Ty = -370.327157
+ 0.008401 0.999959 -0.003329 Tz = -1.606798
+i=1
+matrix global_1 - tr=1 rot=1 refl=0 scl=0
+ 0.940798 0.002136 0.338961 Tx = 163.586079
+ 0.338959 0.001227 -0.940800 Ty = -454.039582
+ -0.002425 0.999997 0.000430 Tz = 0.207687
+i=2
+matrix global_1 - tr=1 rot=1 refl=0 scl=0
+ 0.999976 -0.000756 -0.006876 Tx = -3.318427
+ -0.006877 -0.001172 -0.999976 Ty = -482.598242
+ 0.000748 0.999999 -0.001177 Tz = -0.568114
+i=3
+matrix global_1 - tr=1 rot=1 refl=0 scl=0
+ 0.934365 -0.002092 -0.356311 Tx = -171.959374
+ -0.356314 -0.000801 -0.934366 Ty = -450.934413
+ 0.001669 0.999997 -0.001494 Tz = -0.720939
+
+
+
+
+*/
+
+void MakeFinalAlignment(){
+ // Create ideal (no misalignment) object for PHOS
+
+
+ const AliPHOSGeometry *phosGeom = AliPHOSGeometry::GetInstance("Run2", "Run2");
+ if (!phosGeom) {
+ Error(macroName, "Cannot obtain AliPHOSGeometry singleton.\n");
+ return;
+ }
+
+ //Activate CDB storage and load geometry from CDB
+ //[Part of code, taken from ITS version of MakexxxFullMisalignment
+ AliCDBManager * cdb = AliCDBManager::Instance();
+ cdb->SetDefaultStorage("local://./OCDB");
+ cdb->SetRun(0);
+
+ AliPHOSEMCAGeometry * emca = phosGeom->GetEMCAGeometry();
+ TClonesArray alobj("AliAlignObjParams", 16);
+
+ const Double_t dpsi = 0., dtheta = 0., dphi = 0.;
+ Int_t iIndex = 0; //let all modules have index=0 in a layer with no LUT
+ const AliGeomManager::ELayerID iLayer = AliGeomManager::kInvalidLayer;
+ UShort_t volid = AliGeomManager::LayerToVolUID(iLayer, iIndex);
+ Int_t i = 0;
+
+ // Alignment for 5 PHOS modules
+ Double_t ideal1[9] = {0.766044, 0.000000, 0.642788,
+ 0.642788, 0.000000, -0.766044,
+ 0.000000, 1.000000, 0.000000} ;
+ Double_t final1[9] = {0.767337, -0.004312, 0.641229,
+ 0.641189, -0.007942, -0.767342,
+ 0.008402, 0.999959, -0.003330} ;
+ TGeoRotation rot1;
+ rot1.SetMatrix(ideal1);
+ TGeoRotation inv =rot1.Inverse() ;
+ TGeoRotation rotF1;
+ rotF1.SetMatrix(final1);
+ rotF1.MultiplyBy(&inv) ;
+ Double_t dX= 313.381746-309.463421;
+ Double_t dY=-372.024764+370.327157;
+ Double_t dZ=-2.354855 +1.606798 ;
+
+ AliAlignObjParams * mod1 = new(alobj[i++]) AliAlignObjParams("PHOS/Module1",volid, dX, dY, dZ, 0., 0., 0., kTRUE);
+ mod1->SetRotation(rotF1);
+
+
+ Double_t ideal2[9] = {0.939693, 0.000000, 0.342020,
+ 0.342020, 0.000000, -0.939693,
+ 0.000000, 1.000000, 0.000000} ;
+ Double_t final2[9] = {0.940798, 0.002136, 0.338961,
+ 0.338959, 0.001227, -0.940800,
+ -0.002426, 0.999997, 0.000431} ;
+ TGeoRotation rot2;
+ rot2.SetMatrix(ideal2);
+ TGeoRotation inv2 =rot2.Inverse() ;
+ TGeoRotation rotF2;
+ rotF2.SetMatrix(final2);
+ rotF2.MultiplyBy(&inv2) ;
+ dX= 165.078518-163.586079;
+ dY=-454.148853+454.039582;
+ dZ=-2.801675-0.207687 ;
+
+ AliAlignObjParams * mod2 = new(alobj[i++]) AliAlignObjParams("PHOS/Module2",volid, dX, dY, dZ, 0., 0., 0., kTRUE);
+ mod2->SetRotation(rotF2);
+
+ Double_t ideal3[9] = {1.000000, 0.000000, 0.000000,
+ 0.000000, 0.000000, -1.000000,
+ 0.000000, 1.000000, 0.000000} ;
+ Double_t final3[9] = {0.999976, -0.000756, -0.006876,
+ -0.006877, -0.001172, -0.999976,
+ 0.000748, 0.999999, -0.001178} ;
+ TGeoRotation rot3;
+ rot3.SetMatrix(ideal3);
+ TGeoRotation inv3 =rot3.Inverse() ;
+ TGeoRotation rotF3;
+ rotF3.SetMatrix(final3);
+ rotF3.MultiplyBy(&inv3) ;
+ dX= -0.083626+3.318427;
+ dY=-482.172358+482.598242;
+ dZ=-1.906440+0.568114 ;
+
+ AliAlignObjParams * mod3 = new(alobj[i++]) AliAlignObjParams("PHOS/Module3",volid, dX, dY, dZ, 0., 0., 0., kTRUE);
+ mod3->SetRotation(rotF3);
+
+ Double_t ideal4[9] = {0.939693, 0.000000, -0.342020,
+ -0.342020, 0.000000, -0.939693,
+ 0.000000, 1.000000, 0.000000} ;
+ Double_t final4[9] = {0.934365, -0.002092, -0.356311,
+ -0.356313, -0.000801, -0.934366,
+ 0.001670, 0.999997, -0.001494} ;
+ TGeoRotation rot4;
+ rot4.SetMatrix(ideal4);
+ TGeoRotation inv4 =rot4.Inverse() ;
+ TGeoRotation rotF4;
+ rotF4.SetMatrix(final4);
+ rotF4.MultiplyBy(&inv4) ;
+ dX=-167.233337+171.959374;
+ dY=-460.901723+450.934413;
+ dZ=-1.821245+0.720939 ;
+
+ AliAlignObjParams * mod4 = new(alobj[i++]) AliAlignObjParams("PHOS/Module4",volid, dX, dY, dZ, 0., 0., 0., kTRUE);
+ mod4->SetRotation(rotF4);
+
+ new(alobj[i++]) AliAlignObjParams("PHOS/Module5",
+ volid, 0., 0., 0., 0., 0., 0., kTRUE);
+
+ const Double_t dx = 0., dy = 0., dz = 0. ;
+ // Alignment of CPV modules
+ new(alobj[i++]) AliAlignObjParams("PHOS/Module1/CPV",
+ volid, dx, dy, dz, dpsi, dtheta, dphi, kTRUE);
+ new(alobj[i++]) AliAlignObjParams("PHOS/Module2/CPV",
+ volid, dx, dy, dz, dpsi, dtheta, dphi, kTRUE);
+ new(alobj[i++]) AliAlignObjParams("PHOS/Module3/CPV",
+ volid, dx, dy, dz, dpsi, dtheta, dphi, kTRUE);
+ new(alobj[i++]) AliAlignObjParams("PHOS/Module4/CPV",
+ volid, dx, dy, dz, dpsi, dtheta, dphi, kTRUE);
+ new(alobj[i++]) AliAlignObjParams("PHOS/Module5/CPV",
+ volid, dx, dy, dz, dpsi, dtheta, dphi, kTRUE);
+
+ Double_t displacement=0 ;
+ // Alignment for PHOS cradle
+ new(alobj[i++]) AliAlignObjParams("PHOS/Cradle0",
+ volid, 0., 0., -displacement, dpsi, dtheta, dphi, kTRUE);
+ new(alobj[i++]) AliAlignObjParams("PHOS/Cradle1",
+ volid, 0., 0., +displacement, dpsi, dtheta, dphi, kTRUE);
+
+ // Alignment for cradle wheels
+ new(alobj[i++]) AliAlignObjParams("PHOS/Wheel0",
+ volid, 0., 0., -displacement, dpsi, dtheta, dphi, kTRUE);
+ new(alobj[i++]) AliAlignObjParams("PHOS/Wheel1",
+ volid, 0., 0., -displacement, dpsi, dtheta, dphi, kTRUE);
+ new(alobj[i++]) AliAlignObjParams("PHOS/Wheel2",
+ volid, 0., 0., +displacement, dpsi, dtheta, dphi, kTRUE);
+ new(alobj[i++]) AliAlignObjParams("PHOS/Wheel3",
+ volid, 0., 0., +displacement, dpsi, dtheta, dphi, kTRUE);
+
+ AliCDBMetaData md;
+ md.SetResponsible("Dmitri Peressounko");
+ md.SetComment("Ideal alignment objects for PHOS");
+ md.SetAliRootVersion(gSystem->Getenv("ARVERSION"));
+ AliCDBId id("PHOS/Align/Data",0,AliCDBRunRange::Infinity());
+ cdb->Put(&alobj, id, &md);
+
+ alobj.Delete();
+}
+
+
\ No newline at end of file