LCOV - code coverage report
Current view: top level - MUON/MUONgeometry - AliMUONSurveyObj.cxx (source / functions) Hit Total Coverage
Test: coverage.info Lines: 1 615 0.2 %
Date: 2016-06-14 17:26:59 Functions: 1 37 2.7 %

          Line data    Source code
       1             : /**************************************************************************
       2             :  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
       3             :  *                                                                        *
       4             :  * Author: The ALICE Off-line Project.                                    *
       5             :  * Contributors are mentioned in the code where appropriate.              *
       6             :  *                                                                        *
       7             :  * Permission to use, copy, modify and distribute this software and its   *
       8             :  * documentation strictly for non-commercial purposes is hereby granted   *
       9             :  * without fee, provided that the above copyright notice appears in all   *
      10             :  * copies and that both the copyright notice and this permission notice   *
      11             :  * appear in the supporting documentation. The authors make no claims     *
      12             :  * about the suitability of this software for any purpose. It is          *
      13             :  * provided "as is" without express or implied warranty.                  *
      14             :  **************************************************************************/
      15             : 
      16             : //-----------------------------------------------------------------------------
      17             : /// \class AliMUONSurveyObj
      18             : /// Base class for the survey processing of the ALICE DiMuon spectrometer 
      19             : ///
      20             : /// This base object provides methods to process the survey+photogrammetry
      21             : /// data of the Chambers (frames) and Detection Elements of the DiMuon
      22             : /// Spectrometer and calculate their misalignments. 
      23             : ///
      24             : /// \author Javier Castillo
      25             : //-----------------------------------------------------------------------------
      26             : 
      27             : #include <fstream>
      28             : 
      29             : #include "TMath.h"
      30             : #include "TVector3.h"
      31             : #include "TGeoMatrix.h"
      32             : #include "TFitter.h"
      33             : #include "TMinuit.h"
      34             : #include "TString.h"
      35             : #include "TH2.h"
      36             : #include "TF2.h"
      37             : #include "TGraph2DErrors.h"
      38             : #include "TArrayD.h"
      39             : 
      40             : #include "AliLog.h"
      41             : #include "AliSurveyPoint.h"
      42             : 
      43             : #include "AliMUONSurveyObj.h"
      44             : #include "AliMUONSurveyUtil.h"
      45             : 
      46             : void SurveyFcn(int &npar, double *g, double &f, double *par, int iflag);
      47             : 
      48             : /// \cond CLASSIMP
      49          18 : ClassImp(AliMUONSurveyObj)
      50             : /// \endcond
      51             : 
      52             : AliMUONSurveyObj::AliMUONSurveyObj() 
      53           0 :   : TObject() 
      54           0 :   , fSTargets(0x0)
      55           0 :   , fGBTargets(0x0)
      56           0 :   , fLBTargets(0x0)
      57           0 :   , fLocalTrf(0x0)
      58           0 :   , fAlignTrf(0x0)
      59           0 :   , fBaseTrf(0x0)
      60           0 :   , fOwnerLocalTrf(kFALSE)
      61           0 :   , fOwnerAlignTrf(kTRUE)
      62           0 :   , fOwnerBaseTrf(kFALSE)
      63           0 :   , fUseCM(kTRUE)
      64           0 :   , fPlane(0x0)
      65           0 :   , fFitter(0x0)
      66           0 :   , fXMin(-400.)
      67           0 :   , fXMax(400.)
      68           0 :   , fYMin(-400.)
      69           0 :   , fYMax(400.)
      70           0 :   , fZMin(-2000.)
      71           0 :   , fZMax(2000.)
      72           0 : {
      73             :   /// Default constructor
      74             : 
      75           0 :   fSTargets = new TObjArray();  
      76           0 :   fSTargets->SetOwner(kFALSE);
      77           0 :   fGBTargets = new TObjArray();  
      78           0 :   fGBTargets->SetOwner(kFALSE);
      79           0 :   fLBTargets = new TObjArray();  
      80           0 :   fLBTargets->SetOwner(kFALSE);  
      81             : 
      82           0 :   fAlignTrf = new TGeoCombiTrans();
      83             : 
      84           0 :   fFitter = new TFitter(100);
      85           0 : }
      86             : 
      87           0 : AliMUONSurveyObj::~AliMUONSurveyObj() {
      88             :   /// Destructor
      89           0 :   if(fSTargets) {
      90           0 :     fSTargets->Delete();
      91           0 :     fSTargets = 0x0;
      92           0 :   }
      93           0 :   if(fGBTargets) {
      94           0 :     fGBTargets->Delete();
      95           0 :     fGBTargets = 0x0;
      96           0 :   }
      97           0 :   if(fLBTargets) {
      98           0 :     fLBTargets->Delete();
      99           0 :     fLBTargets = 0x0;
     100           0 :   }
     101           0 :   if (fPlane) {
     102           0 :     fPlane->Delete();
     103           0 :     fPlane = 0x0;
     104           0 :   }
     105           0 :   if(fOwnerLocalTrf && fLocalTrf) {
     106           0 :     fLocalTrf->Delete();
     107           0 :     fLocalTrf = 0x0;
     108           0 :   }
     109           0 :   if(fOwnerAlignTrf && fAlignTrf) {
     110           0 :     fAlignTrf->Delete();
     111           0 :     fAlignTrf = 0x0;
     112           0 :   }
     113           0 :   if(fOwnerBaseTrf && fBaseTrf) {
     114           0 :     fBaseTrf->Delete();
     115           0 :     fBaseTrf = 0x0;
     116           0 :   }
     117           0 :   if (fFitter){
     118           0 :     fFitter->Delete();
     119           0 :     fFitter = 0x0;
     120           0 :   }
     121           0 : }
     122             : 
     123             : void AliMUONSurveyObj::AddStickerTarget(AliSurveyPoint *stPoint){
     124             :   /// Add sticker target
     125           0 :   if (fUseCM) {
     126           0 :     fSTargets->Add(ConvertPointUnits(stPoint,0.1));
     127           0 :   } else {
     128           0 :     fSTargets->Add(stPoint);
     129             :   }
     130           0 : }
     131             : 
     132             : void AliMUONSurveyObj::AddGButtonTarget(AliSurveyPoint *btPoint){
     133             :   /// Add global button target
     134           0 :   if (fUseCM) {
     135           0 :     fGBTargets->Add(ConvertPointUnits(btPoint,0.1));
     136           0 :   } else {
     137           0 :     fGBTargets->Add(btPoint);
     138             :   }  
     139           0 : }
     140             : 
     141             : void AliMUONSurveyObj::AddLButtonTarget(AliSurveyPoint *btPoint){
     142             :   /// Add local button target target; AliSurveyPoint
     143           0 :   if (fUseCM) {
     144           0 :     fLBTargets->Add(ConvertPointUnits(btPoint,0.1));
     145           0 :   } else {
     146           0 :     fLBTargets->Add(btPoint);
     147             :   }  
     148           0 : }
     149             : 
     150             : void AliMUONSurveyObj::AddLButtonTarget(TVector3 *btVector){
     151             :   /// Add local button target target; TVector3
     152           0 :   fLBTargets->Add(btVector);
     153           0 : }
     154             : 
     155             : Int_t AliMUONSurveyObj::AddStickerTargets(TObjArray *pArray, TString stBaseName, Int_t lTargetMax){
     156             :   /// Add a maximum of lTargetMax sticker targets with stBaseName from pArray 
     157           0 :   if (!pArray) {
     158           0 :     AliError(Form("Survey points array is empty %p!",pArray));
     159           0 :     return 0;
     160             :   }
     161           0 :   if (stBaseName.IsNull()) {
     162           0 :     AliError(Form("Need base name for sticker targets %s!",stBaseName.Data()));
     163           0 :     return 0;
     164             :   }
     165             :   
     166             :   Int_t stIndex = 0;
     167             :   AliSurveyPoint *pointSST = 0x0;
     168             : 
     169           0 :   TString stNumber;
     170             :   
     171           0 :   for (int iPoint=0; iPoint<lTargetMax; iPoint++) {
     172           0 :     TString stFullName(stBaseName);
     173           0 :     stNumber = Form("%d",iPoint+1);
     174           0 :     if(lTargetMax>9&&iPoint+1<10) {
     175           0 :       stFullName+="0";
     176             :     }
     177           0 :     stFullName+=stNumber;
     178             : 
     179           0 :     pointSST = (AliSurveyPoint *)pArray->FindObject(stFullName.Data());
     180             : 
     181           0 :     if(pointSST) {
     182           0 :       AddStickerTarget(pointSST);
     183           0 :       AliInfo(Form("Added survey sticker target %s at index %d",pointSST->GetName(),stIndex));
     184           0 :       stIndex++;
     185           0 :     }
     186           0 :   }
     187             : 
     188           0 :   AliInfo(Form("Found %d sticker targets with base name %s",fSTargets->GetEntries(),stBaseName.Data()));
     189             :   return stIndex;
     190           0 : }
     191             : 
     192             : Int_t AliMUONSurveyObj::AddGButtonTargets(TObjArray *pArray, TString btBaseName, Int_t lTargetMax){
     193             :   /// Add a maximum of lTargetMax global button targets with stBaseName from pArray 
     194           0 :   printf("%s \n",btBaseName.Data());
     195           0 :   if (!pArray) {
     196           0 :     AliError(Form("Survey points array is empty %p!",pArray));
     197           0 :     return 0;
     198             :   }
     199           0 :   if (btBaseName.IsNull()) {
     200           0 :     AliError(Form("Need base name for button targets %s!",btBaseName.Data()));
     201           0 :     return 0;
     202             :   }
     203             :   
     204             :   Int_t btIndex = 0;
     205             :   AliSurveyPoint *pointSBT = 0x0;
     206             : 
     207           0 :   TString btNumber;
     208             :   
     209           0 :   for (int iPoint=0; iPoint<lTargetMax; iPoint++) {
     210           0 :     TString btFullName(btBaseName);
     211           0 :     btNumber = Form("%d",iPoint+1);
     212           0 :     if(lTargetMax>9&&iPoint+1<10) {
     213           0 :       btFullName+="0";
     214             :     }
     215           0 :     btFullName+=btNumber;
     216           0 :     printf("%s \n",btFullName.Data());
     217           0 :     pointSBT = (AliSurveyPoint *)pArray->FindObject(btFullName.Data());
     218             : 
     219           0 :     if(pointSBT) {
     220           0 :       AddGButtonTarget(pointSBT);
     221           0 :       AliInfo(Form("Added survey button target %s at index %d",pointSBT->GetName(),btIndex));
     222           0 :       btIndex++;
     223           0 :     }
     224           0 :   }
     225             : 
     226           0 :   AliInfo(Form("Found %d button targets with base name %s",fGBTargets->GetEntries(),btBaseName.Data()));
     227             :   return btIndex;
     228           0 : }
     229             : 
     230             : Int_t AliMUONSurveyObj::AddLButtonTargets(TObjArray *pArray, TString btBaseName, Int_t lTargetMax){
     231             :   /// Add a maximum of lTargetMax local button targets with stBaseName from pArray 
     232           0 :   printf("%s \n",btBaseName.Data());
     233           0 :   if (!pArray) {
     234           0 :     AliError(Form("Local points array is empty %p!",pArray));
     235           0 :     return 0;
     236             :   }
     237           0 :   if (btBaseName.IsNull()) {
     238           0 :     AliError(Form("Need base name for button targets %s!",btBaseName.Data()));
     239           0 :     return 0;
     240             :   }
     241             :   
     242             :   Int_t btIndex = 0;
     243             :   AliSurveyPoint *pointSBT = 0x0;
     244             : 
     245           0 :   TString btNumber;
     246             :   
     247           0 :   for (int iPoint=0; iPoint<lTargetMax; iPoint++) {
     248           0 :     TString btFullName(btBaseName);
     249           0 :     btNumber = Form("%d",iPoint+1);
     250           0 :     if(lTargetMax>9&&iPoint+1<10) {
     251           0 :       btFullName+="0";
     252             :     }
     253           0 :     btFullName+=btNumber;
     254           0 :     printf("%s \n",btFullName.Data());
     255           0 :     pointSBT = (AliSurveyPoint *)pArray->FindObject(btFullName.Data());
     256             : 
     257           0 :     if(pointSBT) {
     258           0 :       AddLButtonTarget(pointSBT);
     259           0 :       AliInfo(Form("Added local button target %s at index %d",pointSBT->GetName(),btIndex));
     260           0 :       btIndex++;
     261           0 :     }
     262           0 :   }
     263             : 
     264           0 :   AliInfo(Form("Found %d local button targets with base name %s",fLBTargets->GetEntries(),btBaseName.Data()));
     265             :   return btIndex;
     266           0 : }
     267             : 
     268             : Int_t AliMUONSurveyObj::GetNStickerTargets() {
     269             :   /// return number of sticker targets
     270           0 :   return fSTargets->GetEntriesFast();
     271             : }
     272             : 
     273             : AliSurveyPoint* AliMUONSurveyObj::GetStickerTarget(Int_t stIndex){
     274             :   /// return sticker target at stIndex
     275           0 :   if (stIndex<0||stIndex>=fSTargets->GetEntriesFast()) {
     276           0 :     AliError(Form("No sticker target at index %d",stIndex));
     277           0 :     return 0x0;
     278             :   }
     279             :   else {
     280           0 :     return (AliSurveyPoint*)fSTargets->At(stIndex);
     281             :   }
     282           0 : }
     283             : 
     284             : Int_t AliMUONSurveyObj::GetNGButtonTargets() {
     285             :   /// return number of global button targets
     286           0 :   return fGBTargets->GetEntriesFast();
     287             : }
     288             : 
     289             : AliSurveyPoint* AliMUONSurveyObj::GetGButtonTarget(Int_t btIndex){
     290             :   /// return global button target at btIndex
     291           0 :   if (btIndex<0||btIndex>=fGBTargets->GetEntriesFast()) {
     292           0 :     AliError(Form("No surveyed button target at index %d",btIndex));
     293           0 :     return 0x0;
     294             :   }
     295             :   else {
     296           0 :     return (AliSurveyPoint*)fGBTargets->At(btIndex);
     297             :   }
     298           0 : }
     299             : 
     300             : Int_t AliMUONSurveyObj::GetNLButtonTargets() {
     301             :   /// return number of local button targets
     302           0 :   return fGBTargets->GetEntriesFast();
     303             : }
     304             : 
     305             : AliSurveyPoint* AliMUONSurveyObj::GetLButtonTarget(Int_t btIndex){
     306             :   /// return local button target at btIndex
     307           0 :   if (btIndex<0||btIndex>=fLBTargets->GetEntriesFast()) {
     308           0 :     AliError(Form("No surveyed button target at index %d",btIndex));
     309           0 :     return 0x0;
     310             :   }
     311             :   else {
     312           0 :     if(fLBTargets->At(btIndex)->IsA()==TVector3::Class()){
     313           0 :       TVector3 *lBT = (TVector3*)fLBTargets->At(btIndex);
     314           0 :       TString str("B");
     315           0 :       return (new AliSurveyPoint(TString("local"),(float)lBT->X(),(float)lBT->Y(),(float)lBT->Z(),(float)0.,(float)0.,(float)0.,'B',kTRUE));
     316           0 :     } else if(fLBTargets->At(btIndex)->IsA()==AliSurveyPoint::Class()) {
     317           0 :       AliSurveyPoint *lBT = (AliSurveyPoint*)fLBTargets->At(btIndex);
     318             :       return lBT;
     319             :     } else {
     320           0 :       AliError(Form("Unexpected class %s ! Valid classes are TVector3 or AliSurveyPoint",fLBTargets->At(btIndex)->ClassName()));
     321           0 :       return 0;
     322             :     }
     323             :   }
     324           0 : }
     325             : 
     326             : void AliMUONSurveyObj::SetPlane(TString pName, Double_t xMin, Double_t xMax, Double_t yMin, Double_t yMax){
     327             :   /// Set the plane function for the plane fitting
     328           0 :   if(fPlane) {
     329           0 :     fPlane->Delete();
     330           0 :     fPlane = 0x0;
     331           0 :   }
     332           0 :   fPlane = new TF2(pName,this,&AliMUONSurveyObj::EqPlane,xMin,xMax,yMin,yMax,3,"AliMUONSurveyObj","EqPlane");
     333           0 : }
     334             : 
     335             : void AliMUONSurveyObj::SetPlaneParameters(Double_t p0, Double_t p1, Double_t p2) {
     336             :   /// Set the parameters of plane function for the plane fitting
     337           0 :   if (!fPlane) {
     338           0 :     AliError("Must use SetPlane before SetPlaneParameters!!!");
     339           0 :   }
     340             :   else {
     341           0 :     fPlane->SetParameter(0,p0);
     342           0 :     fPlane->SetParameter(1,p1);
     343           0 :     fPlane->SetParameter(2,p2);
     344             :   }
     345           0 : }
     346             : 
     347             : void AliMUONSurveyObj::DrawSTargets() {
     348             :   /// Draw a graph of the sticker targets
     349           0 :   TGraph2DErrors *gST = new TGraph2DErrors(3);
     350             :   AliSurveyPoint *pST = 0x0;
     351           0 :   for (Int_t iPoint=0; iPoint<GetNStickerTargets(); iPoint++) {
     352           0 :     pST = GetStickerTarget(iPoint);
     353             :     //    pST->PrintPoint();
     354           0 :     gST->SetPoint(iPoint,pST->GetX(),pST->GetY(),pST->GetZ());
     355           0 :     gST->SetPointError(iPoint,pST->GetPrecisionX(),pST->GetPrecisionY(),pST->GetPrecisionZ());
     356             :   }
     357           0 :   gST->DrawClone("P0");
     358             : 
     359           0 :   delete gST;
     360           0 : }
     361             : 
     362             : Double_t AliMUONSurveyObj::FitPlane() {
     363             :   /// Fit plane to sticker targets
     364           0 :   if (!fPlane) {
     365           0 :     AliError("Must use SetPlane before FitPlane!!!");
     366           0 :     return 0.;
     367             :   }
     368           0 :   if (fSTargets->GetEntriesFast()<3) {
     369           0 :     AliError("Not enough sticker targets (%d) for plane fitting!!!");
     370           0 :     return 0.;
     371             :   }
     372             : 
     373           0 :   Double_t pl[3] = {0};
     374           0 :   Double_t pg[3] = {0};
     375             :     
     376           0 :   TGraph2DErrors *gST = new TGraph2DErrors(3);
     377             :   AliSurveyPoint *pST = 0x0;
     378           0 :   for (Int_t iPoint=0; iPoint<GetNStickerTargets(); iPoint++) {
     379           0 :     pST = GetStickerTarget(iPoint);
     380             :     //    pST->PrintPoint();
     381           0 :     pg[0] =  pST->GetX(); 
     382           0 :     pg[1] =  pST->GetY(); 
     383           0 :     pg[2] =  pST->GetZ(); 
     384           0 :     fBaseTrf->MasterToLocal(pg,pl);
     385           0 :     gST->SetPoint(iPoint,pl[0],pl[1],pl[2]);
     386           0 :     printf("%d %f %f %f\n",iPoint,pl[0],pl[1],pl[2]);
     387           0 :     gST->SetPointError(iPoint,pST->GetPrecisionX(),pST->GetPrecisionY(),pST->GetPrecisionZ());
     388             :   }
     389           0 :   gST->Fit(fPlane);
     390             : 
     391           0 :   delete gST;
     392             : 
     393           0 :   return fPlane->GetChisquare();
     394           0 : }
     395             : 
     396             : Double_t AliMUONSurveyObj::SurveyChi2(Double_t *par){
     397             :   /// Returns the chisquare between local2global transform of local button targets and their surveyed position
     398           0 :   TGeoTranslation transTemp;
     399           0 :   TGeoRotation rotTemp;
     400           0 :   TGeoCombiTrans trfTemp;
     401             : 
     402             :   Double_t lChi2=0.;
     403             : 
     404           0 :   trfTemp.SetTranslation(transTemp);
     405           0 :   trfTemp.SetRotation(rotTemp);
     406           0 :   trfTemp.Clear();
     407           0 :   trfTemp.RotateZ(TMath::RadToDeg()*par[5]);
     408           0 :   trfTemp.RotateY(TMath::RadToDeg()*par[4]);
     409           0 :   trfTemp.RotateX(TMath::RadToDeg()*par[3]);
     410           0 :   trfTemp.SetTranslation(par[0],par[1],par[2]);
     411             : 
     412           0 :   TGeoHMatrix matGlo = (*fBaseTrf)*trfTemp;
     413           0 :   TGeoCombiTrans trfGlo(matGlo);
     414             :                 
     415           0 :   Double_t pl[3] = {0};
     416           0 :   Double_t pg[3] = {0};
     417             :   
     418           0 :   for(Int_t iPoint=0; iPoint<fGBTargets->GetEntries(); iPoint++){
     419           0 :     AliSurveyPoint *gBT = (AliSurveyPoint*)fGBTargets->At(iPoint);
     420           0 :     if(fLBTargets->At(iPoint)->IsA()==TVector3::Class()){
     421           0 :       TVector3 *lBT = (TVector3*)fLBTargets->At(iPoint);
     422           0 :       pl[0]=lBT->X();
     423           0 :       pl[1]=lBT->Y();
     424           0 :       pl[2]=lBT->Z();
     425           0 :     } else if(fLBTargets->At(iPoint)->IsA()==AliSurveyPoint::Class()) {
     426           0 :       AliSurveyPoint *lBT = (AliSurveyPoint*)fLBTargets->At(iPoint);
     427           0 :       pl[0]=lBT->GetX();
     428           0 :       pl[1]=lBT->GetY();
     429           0 :       pl[2]=lBT->GetZ();
     430             :     } else {
     431           0 :       AliError(Form("Unexpected class %s ! Valid classes are TVector3 or AliSurveyPoint",fLBTargets->At(iPoint)->ClassName()));
     432           0 :       return 0;
     433             :     }
     434             : 
     435           0 :     trfGlo.LocalToMaster(pl, pg);
     436             :     //    printf("%d %f %f %f\n",iPoint,pg[0],pg[1],pg[2]);
     437           0 :     if(fLBTargets->At(iPoint)->IsA()==TVector3::Class()){
     438           0 :       lChi2 += (pg[0]-gBT->GetX())*(pg[0]-gBT->GetX())/(gBT->GetPrecisionX()*gBT->GetPrecisionX());
     439           0 :       lChi2 += (pg[1]-gBT->GetY())*(pg[1]-gBT->GetY())/(gBT->GetPrecisionY()*gBT->GetPrecisionY());
     440           0 :       lChi2 += (pg[2]-gBT->GetZ())*(pg[2]-gBT->GetZ())/(gBT->GetPrecisionZ()*gBT->GetPrecisionZ());
     441           0 :     } else if(fLBTargets->At(iPoint)->IsA()==AliSurveyPoint::Class()) {
     442           0 :       AliSurveyPoint *lBT = (AliSurveyPoint*)fLBTargets->At(iPoint);
     443           0 :       lChi2 += (pg[0]-gBT->GetX())*(pg[0]-gBT->GetX())/(gBT->GetPrecisionX()*gBT->GetPrecisionX()+lBT->GetPrecisionX()*lBT->GetPrecisionX());
     444           0 :       lChi2 += (pg[1]-gBT->GetY())*(pg[1]-gBT->GetY())/(gBT->GetPrecisionY()*gBT->GetPrecisionY()+lBT->GetPrecisionY()*lBT->GetPrecisionY());
     445           0 :       lChi2 += (pg[2]-gBT->GetZ())*(pg[2]-gBT->GetZ())/(gBT->GetPrecisionZ()*gBT->GetPrecisionZ()+lBT->GetPrecisionZ()*lBT->GetPrecisionZ());
     446             :     } else {
     447           0 :       AliError(Form("Unexpected class %s ! Valid classes are TVector3 or AliSurveyPoint",fLBTargets->At(iPoint)->ClassName()));
     448           0 :       return 0;
     449             :     }
     450           0 :   }
     451             : 
     452           0 :   return lChi2;
     453           0 : }
     454             : 
     455             : //_____________________________________________________________________________
     456             : void SurveyFcn(int &npar, double *g, double &f, double *par, int iflag) {
     457             :   /// 
     458             :   /// Standard function as needed by Minuit-like minimization procedures. 
     459             :   /// For the set of parameters par calculates and returns chi-squared.
     460             :   ///
     461             : 
     462             :   // smuggle a C++ object into a C function
     463           0 :   AliMUONSurveyObj *aSurveyObj = (AliMUONSurveyObj*) gMinuit->GetObjectFit(); 
     464             : 
     465           0 :   f = aSurveyObj->SurveyChi2(par);
     466             :   if (iflag==3) {}
     467           0 :   if (npar) {} 
     468             :   if (g) {} // no warnings about unused stuff...
     469             : 
     470           0 : }
     471             : 
     472             : //_____________________________________________________________________________
     473             : Int_t AliMUONSurveyObj::SurveyToAlign(TGeoCombiTrans &quadTransf, Double_t *parErr, Double_t psi, Double_t tht, Double_t epsi, Double_t etht) {
     474             :   /// Main function to obtain the misalignments from the surveyed position of the button targets; 
     475           0 :   if (fGBTargets->GetEntries()!=fLBTargets->GetEntries()){
     476           0 :     AliError(Form("Different number of button targets: %d survey points and %d local coord!",
     477             :                   fGBTargets->GetEntries(),fLBTargets->GetEntries()));
     478           0 :     return 0;
     479             :   }
     480             : 
     481           0 :   TFitter fitter(100);
     482           0 :   gMinuit->SetObjectFit(this);
     483           0 :   fitter.SetFCN(SurveyFcn);
     484           0 :   fitter.SetParameter(0,"dx",0,0.1,-20,20);
     485           0 :   fitter.SetParameter(1,"dy",0,0.1,-20,20);
     486           0 :   fitter.SetParameter(2,"dz",0,0.1,0,0);
     487           0 :   fitter.SetParameter(3,"rx",psi,0.0001,-0.05,0.05);
     488           0 :   fitter.SetParameter(4,"ry",tht,0.0001,-0.05,0.05);
     489             : //   fitter.SetParameter(3,"rx",psi,0.0001,psi-5*epsi,psi+5*epsi);
     490             : //   fitter.SetParameter(4,"ry",tht,0.0001,tht-5*etht,tht+5*etht);
     491           0 :   fitter.SetParameter(5,"rz",0,0.0001,0,0);
     492             : 
     493           0 :   if(psi) fitter.FixParameter(3);
     494           0 :   if(tht) fitter.FixParameter(4);
     495             : 
     496           0 :   double arglist[100];
     497           0 :   arglist[0] = 2;
     498           0 :   fitter.ExecuteCommand("SET PRINT", arglist, 1);
     499           0 :   fitter.ExecuteCommand("SET ERR", arglist, 1);
     500           0 :   arglist[0]=0;
     501             :   //fitter.ExecuteCommand("SIMPLEX", arglist, 1);
     502             :   //  fitter.ExecuteCommand("MINIMIZE", arglist, 1);
     503           0 :   fitter.ExecuteCommand("MIGRAD", arglist, 1);
     504           0 :   fitter.ExecuteCommand("IMPROVE", arglist, 1);
     505             :   //  fitter.ExecuteCommand("MINOS", arglist, 1);
     506             :   //  fitter.ExecuteCommand("CALL 3", arglist,0);
     507             : 
     508           0 :   for (int j=0; j<3; j++) printf("%10.3f ",fitter.GetParameter(j));   
     509           0 :   for (int j=3; j<6; j++) printf("%10.3f ",1000*fitter.GetParameter(j));   
     510           0 :   printf("\n");
     511           0 :   for (int j=0; j<3; j++) printf("%10.3f ",fitter.GetParError(j));
     512           0 :   for (int j=3; j<6; j++) printf("%10.3f ",1000*fitter.GetParError(j));
     513           0 :   printf("\n");
     514             : 
     515           0 :   quadTransf.Clear();
     516           0 :   quadTransf.RotateZ(TMath::RadToDeg()*fitter.GetParameter(5));
     517           0 :   quadTransf.RotateY(TMath::RadToDeg()*fitter.GetParameter(4));
     518           0 :   quadTransf.RotateX(TMath::RadToDeg()*fitter.GetParameter(3));
     519           0 :   quadTransf.SetTranslation(fitter.GetParameter(0),fitter.GetParameter(1),fitter.GetParameter(2));
     520             : 
     521           0 :   for(Int_t iPar=0; iPar<6; iPar++){
     522           0 :     parErr[iPar] = fitter.GetParError(iPar);
     523             :   }
     524           0 :   if(epsi) parErr[3] = epsi;
     525           0 :   if(etht) parErr[4] = etht;
     526             : 
     527             :   return 1;
     528             : 
     529           0 : }
     530             : 
     531             : //_____________________________________________________________________________
     532             : Int_t AliMUONSurveyObj::SurveyToAlign(Double_t psi, Double_t tht, Double_t epsi, Double_t etht) {
     533             :   /// Main function to obtain the misalignments from the surveyed position of the button targets; 
     534           0 :   if (fGBTargets->GetEntries()!=fLBTargets->GetEntries()){
     535           0 :     AliError(Form("Different number of button targets: %d survey points and %d local coord!",
     536             :                   fGBTargets->GetEntries(),fLBTargets->GetEntries()));
     537           0 :     return 0;
     538             :   }
     539             : 
     540             :   //  TFitter fitter(100);
     541           0 :   gMinuit->SetObjectFit(this);
     542           0 :   fFitter->SetFCN(SurveyFcn);
     543           0 :   fFitter->SetParameter(0,"dx",0,0.1,-20,20);
     544           0 :   fFitter->SetParameter(1,"dy",0,0.1,-20,20);
     545           0 :   fFitter->SetParameter(2,"dz",0,0.1,0,0);
     546           0 :   if(psi)
     547           0 :       fFitter->SetParameter(3,"rx",psi,epsi,psi-5*epsi,psi+5*epsi);
     548             :   else
     549           0 :     fFitter->SetParameter(3,"rx",psi,0.0001,-0.05,0.05);
     550           0 :   if(tht)
     551           0 :     fFitter->SetParameter(4,"ry",tht,etht,tht-5*etht,tht+5*etht);
     552             :   else
     553           0 :     fFitter->SetParameter(4,"ry",tht,0.0001,-0.05,0.05);
     554             : //   fFitter->SetParameter(3,"rx",psi,0.0001,psi-5*epsi,psi+5*epsi);
     555             : //   fFitter->SetParameter(4,"ry",tht,0.0001,tht-5*etht,tht+5*etht);
     556           0 :   fFitter->SetParameter(5,"rz",0,0.0001,0,0);
     557             : 
     558           0 :   if(psi) fFitter->FixParameter(3);
     559           0 :   if(tht) fFitter->FixParameter(4);
     560             : 
     561           0 :   double arglist[100];
     562           0 :   arglist[0] = 2;
     563           0 :   fFitter->ExecuteCommand("SET PRINT", arglist, 1);
     564           0 :   fFitter->ExecuteCommand("SET ERR", arglist, 1);
     565           0 :   arglist[0]=0;
     566             :   //fFitter->ExecuteCommand("SIMPLEX", arglist, 1);
     567             :   //  fFitter->ExecuteCommand("MINIMIZE", arglist, 1);
     568           0 :   fFitter->ExecuteCommand("MIGRAD", arglist, 1);
     569           0 :   fFitter->ExecuteCommand("IMPROVE", arglist, 1);
     570             : //   fFitter->ExecuteCommand("MINOS", arglist, 1);
     571             : //   fFitter->ExecuteCommand("CALL 3", arglist,0);
     572             : 
     573           0 :   for (int j=0; j<3; j++) printf("%10.3f ",fFitter->GetParameter(j));   
     574           0 :   for (int j=3; j<6; j++) printf("%10.3f ",1000*fFitter->GetParameter(j));   
     575           0 :   printf("\n");
     576           0 :   for (int j=0; j<3; j++) printf("%10.3f ",fFitter->GetParError(j));
     577           0 :   for (int j=3; j<6; j++) printf("%10.3f ",1000*fFitter->GetParError(j));
     578           0 :   printf("\n");
     579             : 
     580           0 :   fAlignTrf->Clear();
     581           0 :   fAlignTrf->RotateZ(TMath::RadToDeg()*fFitter->GetParameter(5));
     582           0 :   fAlignTrf->RotateY(TMath::RadToDeg()*fFitter->GetParameter(4));
     583           0 :   fAlignTrf->RotateX(TMath::RadToDeg()*fFitter->GetParameter(3));
     584           0 :   fAlignTrf->SetTranslation(fFitter->GetParameter(0),fFitter->GetParameter(1),fFitter->GetParameter(2));
     585             : 
     586           0 :   if(epsi) fFitter->ReleaseParameter(3); // To get error
     587           0 :   if(etht) fFitter->ReleaseParameter(4); // To get error
     588             : 
     589           0 :   TGeoCombiTrans lGlobalTrf = TGeoCombiTrans((*fBaseTrf)*(*fAlignTrf));
     590             :   AliSurveyPoint *pointGBT;
     591             :   AliSurveyPoint *pointLBT;
     592           0 :   Double_t pl[3] = {0};
     593           0 :   Double_t pg[3] = {0};
     594           0 :   for (int iPoint=0; iPoint<GetNGButtonTargets(); iPoint++){
     595           0 :     pointGBT=GetGButtonTarget(iPoint);
     596           0 :     pointLBT=GetLButtonTarget(iPoint);
     597           0 :     pl[0] = pointLBT->GetX();
     598           0 :     pl[1] = pointLBT->GetY();
     599           0 :     pl[2] = pointLBT->GetZ();
     600           0 :     lGlobalTrf.LocalToMaster(pl,pg);
     601           0 :     printf("Point %d  local: %.3f %.3f %.3f\n",iPoint,pl[0],pl[1],pl[2]);
     602           0 :     printf("Point %d global: %.3f %.3f %.3f\n",iPoint,pg[0],pg[1],pg[2]);
     603           0 :     printf("Point %d survey: %.3f %.3f %.3f\n",iPoint,pointGBT->GetX(),pointGBT->GetY(),pointGBT->GetZ());
     604             :   }
     605             : 
     606             :   return 1;
     607             : 
     608           0 : }
     609             : 
     610             : Double_t AliMUONSurveyObj::EvalFunction(const TF2 *lFunction, Int_t iP1, Int_t iP2, const Char_t *lCoord) {
     611             :   /// Evaluate the given function at the given points for the given coordinate
     612           0 :   if (!lFunction) {
     613           0 :     AliError("No function given!!!");
     614           0 :     return 0;
     615             :   }
     616           0 :   AliSurveyPoint *gP1 = GetGButtonTarget(iP1);
     617           0 :   AliSurveyPoint *gP2 = GetGButtonTarget(iP2);
     618             :   
     619           0 :   if(!gP1||!gP2){
     620           0 :     AliError("Missing global button target!!!");
     621           0 :     return 0;
     622             :   }
     623             : 
     624             :   //  AliInfo(Form("Function %s parameters %f %f %f %f %f %f",lFunction->GetName(),lFunction->GetParameter(0),lFunction->GetParameter(1),lFunction->GetParameter(2),lFunction->GetParameter(3),lFunction->GetParameter(4),lFunction->GetParameter(5)));
     625           0 :   Double_t pl1[3] = {0};
     626           0 :   Double_t pl2[3] = {0};
     627           0 :   Double_t pg1[3] = {0};
     628           0 :   Double_t pg2[3] = {0};
     629             :     
     630           0 :   pg1[0] =  gP1->GetX(); 
     631           0 :   pg1[1] =  gP1->GetY(); 
     632           0 :   pg1[2] =  gP1->GetZ(); 
     633           0 :   pg2[0] =  gP2->GetX(); 
     634           0 :   pg2[1] =  gP2->GetY(); 
     635           0 :   pg2[2] =  gP2->GetZ(); 
     636             : 
     637           0 :   fBaseTrf->MasterToLocal(pg1,pl1);
     638           0 :   fBaseTrf->MasterToLocal(pg2,pl2);
     639             : 
     640             :   Double_t lVal = 0.;
     641           0 :   switch (lCoord[0]) {
     642             :   case 'X':
     643             :     {
     644           0 :       lVal = lFunction->Eval(pl1[0],pl2[0]);
     645             :       //      lVal = lFunction->Eval(gP1->GetX(),gP2->GetX());
     646             :       //      AliInfo(Form("case X, lVal = %f",lVal));
     647           0 :       return lVal;
     648             :     }
     649             :   case 'Y':
     650             :     {
     651           0 :       lVal = lFunction->Eval(pl1[1],pl2[1]);
     652             :       //      lVal = lFunction->Eval(gP1->GetY(),gP2->GetY());
     653             :       //      AliInfo(Form("case Y, lVal = %f",lVal));
     654           0 :       return lVal;
     655             :     }
     656             :   case 'Z':
     657             :     {
     658           0 :       lVal = lFunction->Eval(pl1[2],pl2[2]);
     659             :       //      lVal = lFunction->Eval(gP1->GetZ(),gP2->GetZ());
     660             :       //      AliInfo(Form("case Z, lVal = %f",lVal));
     661           0 :       return lVal;
     662             :     }
     663             :   default:
     664             :     {
     665           0 :       AliError(Form("Coordinate %s is not valid, options are X Y Z",lCoord));
     666           0 :       return 0;
     667             :     }
     668             :   }
     669           0 : }
     670             : 
     671             : void AliMUONSurveyObj::CalculateTranslation(TF2 *xFunc, TF2 *yFunc, TF2 *zFunc, Int_t iP1, Int_t iP2, Double_t *lCenTemp) {
     672             :   /// Calculate the center translation using analytic functions
     673           0 :   lCenTemp[0] = EvalFunction(xFunc,iP1,iP2,"X");
     674           0 :   lCenTemp[1] = EvalFunction(yFunc,iP1,iP2,"Y");
     675           0 :   lCenTemp[2] = EvalFunction(zFunc,iP1,iP2,"Z");
     676             : 
     677           0 : }
     678             : 
     679             : Double_t AliMUONSurveyObj::CalculateGlobalDiff(TGeoCombiTrans &lTransf, Int_t nPoints, TArrayD &lDiff){
     680             :   /// By hand computation of distance between local2global transform of target position and its surveyed position
     681           0 :   if (nPoints > GetNGButtonTargets()) {
     682           0 :     nPoints = GetNGButtonTargets();
     683           0 :   }
     684             : 
     685           0 :   for(Int_t iVal=0; iVal<nPoints*(3+1)+1; iVal++){
     686           0 :     lDiff[iVal] = 0.;
     687             :   }
     688             : 
     689           0 :   Double_t pl[3] = {0};
     690           0 :   Double_t pg[3] = {0};
     691           0 :   Double_t pml[3] = {0};
     692           0 :   Double_t pmg[3] = {0};
     693             :   AliSurveyPoint *gBT = 0x0;
     694           0 :   for(Int_t iPoint=0; iPoint<nPoints; iPoint++){
     695           0 :     gBT = GetGButtonTarget(iPoint);
     696           0 :     if(!gBT||!fLBTargets->At(iPoint)){
     697           0 :       AliError(Form("The local or global target %d is missing!",iPoint));
     698           0 :       lDiff[nPoints*(3+1)] = 1.e7;
     699           0 :       return lDiff[nPoints*(3+1)];
     700             :     }
     701           0 :     if(fLBTargets->At(iPoint)->IsA()==TVector3::Class()){
     702           0 :       TVector3 *lBT = (TVector3*)fLBTargets->At(iPoint);
     703           0 :       pl[0]=lBT->X();
     704           0 :       pl[1]=lBT->Y();
     705           0 :       pl[2]=lBT->Z();
     706           0 :     } else if(fLBTargets->At(iPoint)->IsA()==AliSurveyPoint::Class()) {
     707           0 :       AliSurveyPoint *lBT = (AliSurveyPoint*)fLBTargets->At(iPoint);
     708           0 :       pl[0]=lBT->GetX();
     709           0 :       pl[1]=lBT->GetY();
     710           0 :       pl[2]=lBT->GetZ();
     711             :     } else {
     712           0 :       AliError(Form("Unexpected class %s ! Valid classes are TVector3 or AliSurveyPoint",fLBTargets->At(iPoint)->ClassName()));
     713           0 :       return 0;
     714             :     }
     715             : 
     716           0 :     lTransf.LocalToMaster(pl,pg);
     717           0 :     pmg[0] = gBT->GetX();
     718           0 :     pmg[1] = gBT->GetY();
     719           0 :     pmg[2] = gBT->GetZ();
     720           0 :     fBaseTrf->MasterToLocal(pmg,pml);
     721             : //     printf("l %d %f %f %f\n",iPoint,pl[0],pl[1],pl[2]);
     722             : //     printf("g %d %f %f %f\n",iPoint,pg[0],pg[1],pg[2]);
     723             : //     printf("ml %d %f %f %f\n",iPoint,pml[0],pml[1],pml[2]);
     724             : //     printf("mg %d %f %f %f\n",iPoint,gBT->GetX(),gBT->GetY(),gBT->GetZ());
     725           0 :     lDiff[iPoint*(3+1)+0] = (pml[0]-pg[0]);
     726           0 :     lDiff[iPoint*(3+1)+1] = (pml[1]-pg[1]);
     727           0 :     lDiff[iPoint*(3+1)+2] = (pml[2]-pg[2]);
     728             :     
     729           0 :     lDiff[iPoint*(3+1)+3] = TMath::Sqrt(lDiff[iPoint*(3+1)+0]*lDiff[iPoint*(3+1)+0]+
     730           0 :                                         lDiff[iPoint*(3+1)+1]*lDiff[iPoint*(3+1)+1]+
     731           0 :                                         lDiff[iPoint*(3+1)+2]*lDiff[iPoint*(3+1)+2]);
     732             : 
     733           0 :     lDiff[nPoints*(3+1)] += lDiff[iPoint*(3+1)+3]*lDiff[iPoint*(3+1)+3];
     734             :   }
     735             :                 
     736           0 :   lDiff[nPoints*(3+1)] = TMath::Sqrt(lDiff[nPoints*(3+1)]);
     737           0 :   return lDiff[nPoints*(3+1)];
     738           0 : }
     739             : 
     740             : Int_t AliMUONSurveyObj::CalculateBestTransf(Int_t iP1, Int_t iP2, Double_t *lXYZ, Double_t *lPTP) {
     741             :   /// By hand calculation of the best local to global transform using 2 button targets
     742           0 :   Double_t lPsi = lPTP[0];
     743           0 :   Double_t lTht = lPTP[1];
     744             : 
     745             :   Double_t pl1[3] = {0};
     746             :   Double_t pl2[3] = {0};
     747             : 
     748           0 :   if(!fLBTargets->At(iP1)||!fLBTargets->At(iP2)){
     749           0 :     AliError(Form("Local target %d or %d is missing!",iP1,iP2));
     750           0 :     return 0;
     751             :   }
     752             : 
     753           0 :   if(fLBTargets->At(iP1)->IsA()==TVector3::Class()){
     754           0 :     TVector3 *lBT1 = (TVector3*)fLBTargets->At(iP1);
     755           0 :     pl1[0]=lBT1->X();
     756           0 :     pl1[1]=lBT1->Y();
     757           0 :     pl1[2]=lBT1->Z();
     758           0 :   } else if(fLBTargets->At(iP1)->IsA()==AliSurveyPoint::Class()) {
     759           0 :     AliSurveyPoint *lBT1 = (AliSurveyPoint*)fLBTargets->At(iP1);
     760           0 :     pl1[0]=lBT1->GetX();
     761           0 :     pl1[1]=lBT1->GetY();
     762           0 :     pl1[2]=lBT1->GetZ();
     763             :   } else {
     764           0 :     AliError(Form("Unexpected class %s ! Valid classes are TVector3 or AliSurveyPoint",fLBTargets->At(iP1)->ClassName()));
     765           0 :     return 0;
     766             :   }
     767           0 :   if(fLBTargets->At(iP2)->IsA()==TVector3::Class()){
     768           0 :     TVector3 *lBT2 = (TVector3*)fLBTargets->At(iP2);
     769           0 :     pl2[0]=lBT2->X();
     770           0 :     pl2[1]=lBT2->Y();
     771           0 :     pl2[2]=lBT2->Z();
     772           0 :   } else if(fLBTargets->At(iP2)->IsA()==AliSurveyPoint::Class()) {
     773           0 :     AliSurveyPoint *lBT2 = (AliSurveyPoint*)fLBTargets->At(iP2);
     774           0 :     pl2[0]=lBT2->GetX();
     775           0 :     pl2[1]=lBT2->GetY();
     776           0 :     pl2[2]=lBT2->GetZ();
     777             :   } else {
     778           0 :     AliError(Form("Unexpected class %s ! Valid classes are TVector3 or AliSurveyPoint",fLBTargets->At(iP2)->ClassName()));
     779           0 :     return 0;
     780             :   }
     781             : 
     782             :   
     783           0 :   AliMUONSurveyUtil *surveyUtil = AliMUONSurveyUtil::Instance();
     784             : 
     785             :   // Xcenter functions
     786             :   const char *fxcName = "fXcn00"; 
     787           0 :   TF2 **fXc = new TF2*[2];
     788             :   fxcName = "fXcn";
     789           0 :   fXc[0] = new TF2(fxcName,surveyUtil,&AliMUONSurveyUtil::XnCenter,fXMin,fXMax,fYMin,fYMax,7,"AliMUONSurveyUtil","XnCenter");
     790             :   fxcName = "fXcp";
     791           0 :   fXc[1] = new TF2(fxcName,surveyUtil,&AliMUONSurveyUtil::XpCenter,fXMin,fXMax,fYMin,fYMax,7,"AliMUONSurveyUtil","XpCenter");
     792             : 
     793             :   // Ycenter functions
     794             :   const char *fycName = "fYcn00"; 
     795           0 :   TF2 **fYc = new TF2*[2];
     796             :   fycName = "fYcn";
     797           0 :   fYc[0] = new TF2(fycName,surveyUtil,&AliMUONSurveyUtil::YnCenter,fYMin,fYMax,fYMin,fYMax,8,"AliMUONSurveyUtil","YnCenter");
     798             :   fycName = "fYcp";
     799           0 :   fYc[1] = new TF2(fycName,surveyUtil,&AliMUONSurveyUtil::YpCenter,fYMin,fYMax,fYMin,fYMax,8,"AliMUONSurveyUtil","YpCenter");   
     800             : 
     801             :   // Zcenter functions
     802             :   const char *fzcName = "fZcn00"; 
     803           0 :   TF2 **fZc = new TF2*[2];
     804             :   fzcName = "fZcn";
     805           0 :   fZc[0] = new TF2(fzcName,surveyUtil,&AliMUONSurveyUtil::ZnCenter,fZMin,fZMax,fZMin,fZMax,8,"AliMUONSurveyUtil","ZnCenter");
     806             :   fzcName = "fZcp";
     807           0 :   fZc[1] = new TF2(fzcName,surveyUtil,&AliMUONSurveyUtil::ZpCenter,fZMin,fZMax,fZMin,fZMax,8,"AliMUONSurveyUtil","ZpCenter");   
     808             : 
     809             :   // Phi rotation using xglobal coords functions
     810             :   const char *fphixName = "fPhiXnn00"; 
     811           0 :   TF2 ***fPhiX = new TF2**[2];
     812           0 :   for (Int_t iX =0; iX<2; iX++) {
     813           0 :     fPhiX[iX] = new TF2*[2];
     814             :   }
     815             :   fphixName = "fPhiXnn";
     816           0 :   fPhiX[0][0] = new TF2(fphixName,surveyUtil,&AliMUONSurveyUtil::PhiXnn,fXMin,fXMax,fXMin,fXMax,7,"AliMUONSurveyUtil","PhiXnn");
     817             :   fphixName = "fPhiXnp";
     818           0 :   fPhiX[0][1] = new TF2(fphixName,surveyUtil,&AliMUONSurveyUtil::PhiXnp,fXMin,fXMax,fXMin,fXMax,7,"AliMUONSurveyUtil","PhiXnp");   
     819             :   fphixName = "fPhiXpn";
     820           0 :   fPhiX[1][0] = new TF2(fphixName,surveyUtil,&AliMUONSurveyUtil::PhiXpn,fXMin,fXMax,fXMin,fXMax,7,"AliMUONSurveyUtil","PhiXpn");
     821             :   fphixName = "fPhiXpp";
     822           0 :   fPhiX[1][1] = new TF2(fphixName,surveyUtil,&AliMUONSurveyUtil::PhiXpp,fXMin,fXMax,fXMin,fXMax,7,"AliMUONSurveyUtil","PhiXpp");   
     823             : 
     824             :   // Phi rotation using yglobal coords functions
     825             :   const char *fphiyName = "fPhiYnn00"; 
     826           0 :   TF2 ***fPhiY = new TF2**[2];
     827           0 :   for (Int_t iY =0; iY<2; iY++) {
     828           0 :     fPhiY[iY] = new TF2*[2];
     829             :   }
     830             :   fphiyName = "fPhiYnn";
     831           0 :   fPhiY[0][0] = new TF2(fphiyName,surveyUtil,&AliMUONSurveyUtil::PhiYnn,fYMin,fYMax,fYMin,fYMax,8,"AliMUONSurveyUtil","PhiYnn");
     832             :   fphiyName = "fPhiYnp";
     833           0 :   fPhiY[0][1] = new TF2(fphiyName,surveyUtil,&AliMUONSurveyUtil::PhiYnp,fYMin,fYMax,fYMin,fYMax,8,"AliMUONSurveyUtil","PhiYnp");   
     834             :   fphiyName = "fPhiYpn";
     835           0 :   fPhiY[1][0] = new TF2(fphiyName,surveyUtil,&AliMUONSurveyUtil::PhiYpn,fYMin,fYMax,fYMin,fYMax,8,"AliMUONSurveyUtil","PhiYpn");
     836             :   fphiyName = "fPhiYpp";
     837           0 :   fPhiY[1][1] = new TF2(fphiyName,surveyUtil,&AliMUONSurveyUtil::PhiYpp,fYMin,fYMax,fYMin,fYMax,8,"AliMUONSurveyUtil","PhiYpp");   
     838             :   
     839             : 
     840             :   // Set Parameters of functions
     841           0 :   for(Int_t iS=0; iS<2; iS++){
     842           0 :     fXc[iS]->SetParameters(pl1[0],pl1[1],pl1[2],pl2[0],pl2[1],pl2[2],lTht);
     843           0 :     fYc[iS]->SetParameters(pl1[0],pl1[1],pl1[2],pl2[0],pl2[1],pl2[2],lPsi,lTht);
     844           0 :     fZc[iS]->SetParameters(pl1[0],pl1[1],pl1[2],pl2[0],pl2[1],pl2[2],lPsi,lTht);
     845             : //     fXc[iS]->SetParameters(lBT1->X(),lBT1->Y(),lBT1->Z(),lBT2->X(),lBT2->Y(),lBT2->Z(),lTht);
     846             : //     fYc[iS]->SetParameters(lBT1->X(),lBT1->Y(),lBT1->Z(),lBT2->X(),lBT2->Y(),lBT2->Z(),lPsi,lTht);
     847             : //     fZc[iS]->SetParameters(lBT1->X(),lBT1->Y(),lBT1->Z(),lBT2->X(),lBT2->Y(),lBT2->Z(),lPsi,lTht);
     848           0 :     for(Int_t jS=0; jS<2; jS++){
     849           0 :       fPhiX[iS][jS]->SetParameters(pl1[0],pl1[1],pl1[2],pl2[0],pl2[1],pl2[2],lTht);
     850           0 :       fPhiY[iS][jS]->SetParameters(pl1[0],pl1[1],pl1[2],pl2[0],pl2[1],pl2[2],lPsi,lTht);
     851             : //     fPhiX[iS][jS]->SetParameters(lBT1->X(),lBT1->Y(),lBT1->Z(),lBT2->X(),lBT2->Y(),lBT2->Z(),lTht);
     852             : //     fPhiY[iS][jS]->SetParameters(lBT1->X(),lBT1->Y(),lBT1->Z(),lBT2->X(),lBT2->Y(),lBT2->Z(),lPsi,lTht);
     853             :     }
     854             :   }
     855             : 
     856           0 :   Double_t lCenTemp[3];
     857             :   Double_t lRotTemp[3];
     858             : 
     859           0 :   TGeoCombiTrans trfTemp;
     860             : 
     861           0 :   Int_t nPoints = GetNGButtonTargets();
     862             : 
     863           0 :   TArrayD lDiffTemp(nPoints*(3+1)+1);
     864           0 :   TArrayD lDiffMin(nPoints*(3+1)+1);
     865             : 
     866           0 :   for(Int_t i=0; i<nPoints*(3+1)+1; i++){
     867           0 :     lDiffMin[i]=1000000.;
     868           0 :     lDiffTemp[i]=0.;
     869             :   }
     870             : 
     871             :   //
     872             :   // Calculate Detection Element Center from button targets
     873             :   //    
     874             :   
     875             :   // Trying 2x*2y*2z*(2phi+2phi) possibilities
     876           0 :   for(Int_t iX=0; iX<2; iX++){
     877           0 :     for(Int_t iY=0; iY<2; iY++){
     878           0 :       for(Int_t iZ=0; iZ<2; iZ++){
     879           0 :         CalculateTranslation(fXc[iX],fYc[iY],fZc[iZ],iP1,iP2,lCenTemp);
     880             :         
     881             :         lRotTemp[0] = lPsi;
     882             :         lRotTemp[1] = lTht;
     883           0 :         for(Int_t iP=0; iP<2; iP++){
     884           0 :           lRotTemp[2] = EvalFunction(fPhiX[iX][iP],iP1,iP2,"X");
     885             :           
     886           0 :           trfTemp.Clear();
     887           0 :           trfTemp.RotateZ(TMath::RadToDeg()*lRotTemp[2]);
     888           0 :           trfTemp.RotateY(TMath::RadToDeg()*lRotTemp[1]);
     889           0 :           trfTemp.RotateX(TMath::RadToDeg()*lRotTemp[0]);
     890           0 :           trfTemp.SetTranslation(lCenTemp);
     891             :           
     892           0 :           if(CalculateGlobalDiff(trfTemp,nPoints,lDiffTemp)<lDiffMin[nPoints*(3+1)]){
     893           0 :             printf("Diffs");
     894           0 :             for(Int_t i=0; i<nPoints*(3+1)+1; i++){
     895           0 :               printf(" %f",lDiffTemp[i]); 
     896             :             }
     897           0 :             printf("\n");
     898           0 :             printf(" : mycenX%dY%dZ%d(%f,%f,%f); rotx%d(%f,%f,%f)\n",iX,iY,iZ,lCenTemp[0],lCenTemp[1],lCenTemp[2],iP,lRotTemp[0],lRotTemp[1],lRotTemp[2]);
     899           0 :             printf("Transformation improved ...\n");
     900           0 :             for (int i=0; i<3; i++) {
     901           0 :               lXYZ[i] = lCenTemp[i];
     902             :             } 
     903           0 :             lPTP[2] = lRotTemp[2];
     904           0 :             for(Int_t i=0; i<nPoints*(3+1)+1; i++){
     905           0 :               lDiffMin[i]=lDiffTemp[i];
     906             :             }
     907           0 :           }
     908             :         }
     909           0 :         for(Int_t iP=0; iP<2; iP++){
     910           0 :           lRotTemp[2] = EvalFunction(fPhiY[iY][iP],iP1,iP2,"Y");
     911             :           
     912           0 :           trfTemp.Clear();
     913           0 :           trfTemp.RotateZ(TMath::RadToDeg()*lRotTemp[2]);
     914           0 :           trfTemp.RotateY(TMath::RadToDeg()*lRotTemp[1]);
     915           0 :           trfTemp.RotateX(TMath::RadToDeg()*lRotTemp[0]);
     916           0 :           trfTemp.SetTranslation(lCenTemp);
     917             :           
     918           0 :           if(CalculateGlobalDiff(trfTemp,nPoints,lDiffTemp)<lDiffMin[nPoints*(3+1)]){
     919           0 :             printf("Diffs");
     920           0 :             for(Int_t i=0; i<nPoints*(3+1)+1; i++){
     921           0 :               printf(" %f",lDiffTemp[i]); 
     922             :             }
     923           0 :             printf("\n");
     924           0 :             printf(" : mycenX%dY%dZ%d(%f,%f,%f); roty%d(%f,%f,%f)\n",iX,iY,iZ,lCenTemp[0],lCenTemp[1],lCenTemp[2],iP,lRotTemp[0],lRotTemp[1],lRotTemp[2]);
     925           0 :             printf("Transformation improved ...\n");
     926           0 :             for (int i=0; i<3; i++) {
     927           0 :               lXYZ[i] = lCenTemp[i];
     928             :             } 
     929           0 :             lPTP[2] = lRotTemp[2];
     930           0 :             for(Int_t i=0; i<nPoints*(3+1)+1; i++){
     931           0 :               lDiffMin[i]=lDiffTemp[i];
     932             :             }
     933           0 :           }
     934             :         }
     935             :       }
     936             :     }
     937             :   }
     938             : 
     939           0 :   for (Int_t i=0; i<2; i++) {
     940           0 :     delete fXc[i];
     941           0 :     delete fYc[i];
     942           0 :     delete fZc[i];
     943           0 :     for (Int_t j=0; j<2; j++) {
     944           0 :       delete fPhiX[i][j];
     945           0 :       delete fPhiY[i][j];
     946             :     }
     947           0 :     delete[] fPhiX[i];
     948           0 :     delete[] fPhiY[i];
     949             :   }
     950           0 :   delete[] fXc;
     951           0 :   delete[] fYc;
     952           0 :   delete[] fZc;
     953           0 :   delete[] fPhiX;
     954           0 :   delete[] fPhiY;
     955             : 
     956           0 :   if (lDiffMin[nPoints*(3+1)]>20) return 0;
     957             : 
     958           0 :   return 1;
     959           0 : }
     960             : 
     961             : void AliMUONSurveyObj::CalculateMeanTransf(Double_t *lXYZ, Double_t *lPTP) {
     962             :   /// By hand calculation of the mean (for nPairs of targets) of the best local to global transform using 2 button targets
     963             :     Double_t xce=0.;
     964             :     Double_t yce=0.;
     965             :     Double_t zce=0.;
     966             :     Double_t phi=0.;
     967             :     
     968             :     Int_t nPairs = 0;
     969           0 :     Int_t nPoints = GetNGButtonTargets();
     970             :     // Loop over all possible pairs of button tragets
     971           0 :     for(Int_t iP1=0; iP1<nPoints; iP1++){
     972           0 :       for(Int_t iP2=iP1+1; iP2<nPoints; iP2++){
     973           0 :         printf("%d and %d\n",iP1,iP2);
     974             : 
     975           0 :         if(CalculateBestTransf(iP1,iP2,lXYZ,lPTP)) {
     976           0 :           nPairs++;
     977             :         
     978           0 :           xce+=lXYZ[0];
     979           0 :           yce+=lXYZ[1];
     980           0 :           zce+=lXYZ[2];
     981           0 :           phi+=lPTP[2];
     982           0 :         }
     983             :       }
     984             :     }
     985             : 
     986           0 :     if (!nPairs) return;
     987             :     
     988           0 :     lXYZ[0]=xce/nPairs;
     989           0 :     lXYZ[1]=yce/nPairs;
     990           0 :     lXYZ[2]=zce/nPairs;
     991           0 :     lPTP[2]=phi/nPairs;
     992           0 : }
     993             : 
     994             : void AliMUONSurveyObj::PrintLocalTrf() {
     995             :   /// Print the local transformation
     996           0 :   Double_t lRotTemp[3];
     997           0 :   AliMUONSurveyUtil::MatrixToAngles(fLocalTrf->GetRotationMatrix(),lRotTemp);
     998           0 :   printf("(%.3f %.3f %.3f), (%.6f %.6f %.6f)\n",fLocalTrf->GetTranslation()[0],fLocalTrf->GetTranslation()[1],fLocalTrf->GetTranslation()[2],lRotTemp[0],lRotTemp[1],lRotTemp[2]);
     999           0 : }
    1000             : 
    1001             : void AliMUONSurveyObj::PrintAlignTrf() {
    1002             :   /// Print the alignment transformation
    1003           0 :   Double_t lRotTemp[3];
    1004           0 :   AliMUONSurveyUtil::MatrixToAngles(fAlignTrf->GetRotationMatrix(),lRotTemp);
    1005           0 :   printf("(%.3f %.3f %.3f), (%.6f %.6f %.6f)\n",fAlignTrf->GetTranslation()[0],fAlignTrf->GetTranslation()[1],fAlignTrf->GetTranslation()[2],lRotTemp[0],lRotTemp[1],lRotTemp[2]);
    1006           0 : }
    1007             : 
    1008             : void AliMUONSurveyObj::FillSTHistograms(TString baseNameC, TH2 *hSTc, TString baseNameA, TH2 *hSTa) {
    1009             :   /// Fill sticker target histograms for monitoring
    1010           0 :   if(baseNameC.IsNull()||!hSTc){
    1011           0 :     AliError("Need base name for points on side C and/or a histogram for them!");
    1012           0 :     return;
    1013             :   }
    1014             :   AliSurveyPoint *pointST = 0x0;
    1015           0 :   for (Int_t iPoint=0; iPoint<GetNStickerTargets(); iPoint++) {
    1016           0 :     pointST = GetStickerTarget(iPoint);
    1017           0 :     if (!pointST) continue;
    1018           0 :     if (pointST->GetPointName().Contains(baseNameC)){
    1019           0 :       hSTc->Fill(pointST->GetX(),pointST->GetY(),-pointST->GetZ());
    1020           0 :     } else if ((!baseNameA.IsNull()) && 
    1021           0 :               (pointST->GetPointName().Contains(baseNameA))) {
    1022           0 :       if (!hSTa){
    1023           0 :         AliError("Base name for points on side A provided but no histogram for them!");
    1024           0 :         continue;
    1025             :       }
    1026           0 :       hSTa->Fill(pointST->GetX(),pointST->GetY(),-pointST->GetZ());
    1027           0 :     }
    1028             :   }
    1029           0 : }
    1030             : 
    1031             : Double_t AliMUONSurveyObj::GetAlignResX() {
    1032             :   /// Returns the uncertainty of the x translation parameter 
    1033           0 :   if(!fFitter) {
    1034           0 :     AliError("There is no fitter for this object! X resolution will be 0.");
    1035           0 :     return 0.;
    1036             :   }
    1037           0 :   return fFitter->GetParError(0);
    1038           0 : }
    1039             : 
    1040             : Double_t AliMUONSurveyObj::GetAlignResY() {
    1041             :   /// Returns the uncertainty of the y translation parameter 
    1042           0 :   if(!fFitter) {
    1043           0 :     AliError("There is no fitter for this object! Y resolution will be 0.");
    1044           0 :     return 0.;
    1045             :   }
    1046           0 :   return fFitter->GetParError(1);
    1047           0 : }
    1048             : 
    1049             : AliSurveyPoint* AliMUONSurveyObj::ConvertPointUnits(AliSurveyPoint *stPoint, Float_t lFactor) {
    1050             :   /// Return the AliSurveyPoint with new units. Default is from mm -> cm 
    1051           0 :   return new AliSurveyPoint(stPoint->GetPointName(),
    1052           0 :                             lFactor*stPoint->GetX(),lFactor*stPoint->GetY(),lFactor*stPoint->GetZ(),
    1053           0 :                             lFactor*stPoint->GetPrecisionX(),lFactor*stPoint->GetPrecisionY(),lFactor*stPoint->GetPrecisionZ(),
    1054           0 :                             stPoint->GetType(), stPoint->GetTarget());
    1055           0 : }

Generated by: LCOV version 1.11