]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - FMD/AliFMDSurveyToAlignObjs.cxx
Configurable when to add internally
[u/mrichter/AliRoot.git] / FMD / AliFMDSurveyToAlignObjs.cxx
index f236cfd3113da6d4ed12f5ca0a91dfd767690d23..cbfb403e71613b51b30541a204c1a856cfcdb916 100644 (file)
@@ -1,3 +1,9 @@
+//
+// Class to take survey data and 
+// transform that to alignment objects. 
+// 
+// FMD
+//
 #include "AliFMDSurveyToAlignObjs.h"
 #include "AliLog.h"
 #include "AliSurveyPoint.h"
@@ -64,9 +70,23 @@ AliFMDSurveyToAlignObjs::CalculatePlane(const     TVector3& a,
                                        Double_t* trans,
                                        Double_t* rot) const
 {
+  // 
+  // Calculate the plane translation and rotation from 3 survey points
+  // 
+  // Parameters:
+  //    a     1st Survey point 
+  //    b     2nd Survey point
+  //    c     3rd Survey point
+  //    trans Translation vector
+  //    rot   Rotation matrix (direction cosines)
+  // 
+  // Return:
+  //    
+  //
+
   // Vector a->b, b->c, and normal to plane defined by these two
   // vectors. 
-  TVector3 ab(b-a), bc(c-b);
+  TVector3 ab(b-a), bc(c-a);
   
   // Normal vector to the plane of the fiducial marks obtained
   // as cross product of the two vectors on the plane d0^d1
@@ -96,14 +116,22 @@ AliFMDSurveyToAlignObjs::CalculatePlane(const     TVector3& a,
   // The center of the square with the fiducial marks as the
   // corners.  The mid-point of one diagonal - md.  Used to get the
   // center of the surveyd box. 
-  TVector3 md(a + c);
-  md *= 1/2.;
+  // TVector3 md(a + c);
+  // md *= 1/2.;
+  //Info("CalculatePlane", "corner=(%8f,%8f,%8f)", c.X(),c.Y(),c.Z());
+  //Info("CalculatePlane", "corner=(%8f,%8f,%8f)", b.X(),b.Y(),b.Z());
+  TVector3 md(c + b);
+  md *= 1./2;
+  //Info("CalculatePlane", "mid=(%8f,%8f,%8f)", md.X(),md.Y(),md.Z());
+  //Info("CalculatePlane", "normal=(%8f,%8f,%8f)", n.X(),n.Y(),n.Z());
 
   // The center of the box. 
   TVector3 orig(md - depth * n);
+  // Info("CalculatePlane", "orig=(%8f,%8f,%8f)", orig.X(),orig.Y(),orig.Z());
   trans[0] = orig[0];
   trans[1] = orig[1];
   trans[2] = orig[2];
+  //Info("CalculatePlane", "trans=(%8f,%8f,%8f)", trans[0],trans[1],trans[2]);
   
   // Normalize the spanning vectors 
   TVector3 uab(ab.Unit());
@@ -112,6 +140,8 @@ AliFMDSurveyToAlignObjs::CalculatePlane(const     TVector3& a,
   for (size_t i = 0; i < 3; i++) { 
     rot[i * 3 + 0] = ubc[i];
     rot[i * 3 + 1] = uab[i];
+    // rot[i * 3 + 0] = uab[i];
+    // rot[i * 3 + 1] = ubc[i];
     rot[i * 3 + 2] = n[i];
   }
   return kTRUE;
@@ -125,6 +155,24 @@ AliFMDSurveyToAlignObjs::FitPlane(const TObjArray& points,
                                  Double_t*        trans,
                                  Double_t*        rot) const
 {
+  // 
+  // Calculate the plane rotation and translation by doing a fit of
+  // the plane equation to the surveyed points.  At least 4 points
+  // must be passed in the @a points array with corresponding errors
+  // in the array @a errors.  The arrays are assumed to contain
+  // TVector3 objects.
+  // 
+  // Parameters:
+  //    points Array surveyed positions
+  //    errors Array of errors corresponding to @a points
+  //    depth  Survey targets depth (perpendicular to the plane)
+  //    trans  On return, translation of the plane
+  //    rot    On return, rotation (direction cosines) of the plane
+  // 
+  // Return:
+  //    @c true on success, @c false otherwise
+  //
+
   Int_t nPoints = points.GetEntries();
   if (nPoints < 4) { 
     AliError(Form("Cannot fit a plane equation to less than 4 survey points, "
@@ -142,6 +190,7 @@ AliFMDSurveyToAlignObjs::FitPlane(const TObjArray& points,
     
     g.SetPoint(i, p->X(), p->Y(), p->Z());
     g.SetPointError(i, e->X(), e->Y(), e->Z());
+
   }
 
   // Check that we have enough points
@@ -197,11 +246,23 @@ AliFMDSurveyToAlignObjs::FitPlane(const TObjArray& points,
 //____________________________________________________________________
 Bool_t
 AliFMDSurveyToAlignObjs::MakeDelta(const char*  path, 
-                                  Double_t*    rot, 
-                                  Double_t*    trans, 
+                                  const Double_t*    rot, 
+                                  const Double_t*    trans, 
                                   TGeoHMatrix& delta) const
 {
-  // Make delta transformation
+  // 
+  // Create a delta transform from a global rotation matrix and
+  // translation. 
+  // 
+  // Parameters:
+  //    path   Path of element to transform.
+  //    rot    Rotation matrix (direction cosines)
+  //    trans  Translation 
+  //    delta  On return, the delta transform
+  // 
+  // Return:
+  //    Newly 
+  //
   if (!gGeoManager)           return kFALSE;
   if (!gGeoManager->cd(path)) return kFALSE;
   
@@ -217,34 +278,89 @@ AliFMDSurveyToAlignObjs::MakeDelta(const char*  path,
 
 //____________________________________________________________________
 Bool_t
-AliFMDSurveyToAlignObjs::MakeDelta(TGeoMatrix*  global,
-                                  Double_t*    rot, 
-                                  Double_t*    trans, 
+AliFMDSurveyToAlignObjs::MakeDelta(const TGeoMatrix*  global,
+                                  const Double_t*    rot, 
+                                  const Double_t*    trans, 
                                   TGeoHMatrix& delta) const
 {
-  // Make delta transformation
+  // 
+  // Create a delta transform from a global rotation matrix and
+  // translation. 
+  // 
+  // Parameters:
+  //    global Global matrix of element to transform.
+  //    rot    Rotation matrix (direction cosines)
+  //    trans  Translation 
+  //    delta  On return, the delta transform
+  // 
+  // Return:
+  //    Newly 
+  //
   TGeoHMatrix* geoM = new TGeoHMatrix;
   geoM->SetTranslation(trans);
   geoM->SetRotation(rot);
+  // Info("MakeDelta", "The HMatrix from survey");
+  // geoM->Print();
+  // Info("MakeDelta", "The global matrix");
+  // global->Print();
 
   delta = global->Inverse();
+  // Info("MakeDelta", "The inverse global matrix");
+  // delta.Print();
   delta.MultiplyLeft(geoM);
-  
+  // Info("MakeDelta", "The delta matrix");
+  // delta.Print();
   return true;
 }
 
+namespace {
+  Double_t getFMD1Offset()
+  {
+    static Double_t off = 0;
+    return off;
+
+#if 0
+    if (off != 0) return off;
+    
+    const char* lidN = "FMD1_lid_mat0";
+    TGeoMatrix* lidM = static_cast<TGeoMatrix*>(gGeoManager->GetListOfMatrices()
+                                               ->FindObject(lidN));
+    if (!lidM) { 
+      Error("getFMD1Offset", "Couldn't find FMD1 lid transformation %s", lidN);
+      return 0;
+    }
+
+    const Double_t* lidT = lidM->GetTranslation();
+    Double_t        lidZ = lidT[2];
+    off                  = lidZ-3.3;
+    
+    return off;
+#endif
+  }
+}
+
 //____________________________________________________________________
 Bool_t
 AliFMDSurveyToAlignObjs::GetFMD1Plane(Double_t* rot, Double_t* trans) const
 {
+  // 
+  // Get the FMD1 plane from the survey points
+  // 
+  // Parameters:
+  //    rot    Rotation matrix (direction cosines)
+  //    trans  Translation
+  // 
+  // Return:
+  //    @c true on success, @c false otherwise.
+  //
 
   // The possile survey points 
-  TVector3  icb, ict, ocb, oct, dummy;
+  TVector3  icb, ict, ocb, oct, eicb, eict, eocb, eoct;
   Int_t     missing = 0;
-  if (!GetPoint("V0L_ICB", icb, dummy)) missing++;
-  if (!GetPoint("V0L_ICT", ict, dummy)) missing++;
-  if (!GetPoint("V0L_OCB", ocb, dummy)) missing++;
-  if (!GetPoint("V0L_OCT", oct, dummy)) missing++;
+  if (!GetPoint("V0L_ICB", icb, eicb)) missing++;
+  if (!GetPoint("V0L_ICT", ict, eict)) missing++;
+  if (!GetPoint("V0L_OCB", ocb, eocb)) missing++;
+  if (!GetPoint("V0L_OCT", oct, eoct)) missing++;
 
   // Check that we have enough points
   if (missing > 1) { 
@@ -252,41 +368,86 @@ AliFMDSurveyToAlignObjs::GetFMD1Plane(Double_t* rot, Double_t* trans) const
                    4-missing));
     return kFALSE;
   }
-
 #if 0
-  const char* lidN = "FMD1_lid_mat0";
-  TGeoMatrix* lidM = static_cast<TGeoMatrix*>(gGeoManager->GetListOfMatrices()
-                                             ->FindObject(lidN));
-  if (!lidM) { 
-    AliError(Form("Couldn't find FMD1 lid transformation %s", lidN));
-    return kFALSE;
+  TObjArray points;
+  TObjArray errors;
+  points.Add(&icb); errors.Add(&eicb);
+  points.Add(&ict); errors.Add(&eict);
+  points.Add(&oct); errors.Add(&eoct);
+  points.Add(&ocb); errors.Add(&eocb);
+  
+  Bool_t ret = FitPlane(points, errors, 0, trans, rot);
+  if (!ret) { 
+    Warning("GetFMD1Plane", "fit to plane failed");
+  }
+  for (Int_t i = 0; i < 4; i++) { 
+    TVector3* v = static_cast<TVector3*>(points.At(i));
+    TVector3* e = static_cast<TVector3*>(errors.At(i));
+    Info("GetFMD1Plane", "p%d=(%8f,%8f,%8f)+/-(%8f,%8f,%8f)", 
+        i, v->X(), v->Y(), v->Z(), e->X(), e->Y(), e->Z());
   }
-
-  const Double_t* lidT = lidM->GetTranslation();
-  Double_t        lidZ = lidT[2];
-  Double_t        off  = lidZ-3.3;
 #else
-  Double_t        off  = 0;
-#endif
+  Double_t off = getFMD1Offset();
+  Info("GetFMD1Plane", "Lid offset is %f", off);
 
-  if (!CalculatePlane(ocb, icb, ict, off, trans, rot)) return kFALSE;
-  // PrintRotation("FMD1 rotation:",  rot);
-  // PrintVector("FMD1 translation:", trans);
+  // if (!CalculatePlane(ocb, icb, ict, off, trans, rot)) return kFALSE;
+  // Bool_t ret = CalculatePlane(ocb, icb, oct, off, trans, rot);
+  Bool_t ret = CalculatePlane(oct, ocb, ict, off, trans, rot);
+#endif
+  PrintRotation("FMD1 rotation:",  rot);
+  PrintVector("FMD1 translation:", trans);
 
-  return kTRUE;
+  return ret;
 }
 
 //____________________________________________________________________
 Bool_t
 AliFMDSurveyToAlignObjs::DoFMD1()
 {
+  // 
+  // Do the FMD1 analysis.  We have 4 survey targets on V0-A on the
+  // C-side.  These are 
+  //
+  //  - V0A_ICT  In-side, C-side, top.
+  //  - V0A_ICB  In-side, C-side, bottom.  
+  //  - V0A_OCT  Out-side, C-side, top.         
+  //  - V0A_OCB         Out-side, C-side, bottom.
+  // 
+  // These 4 survey targets sit 3.3mm over the V0-A C-side surface, or
+  // 3.3mm over the back surface of FMD1.  
+  //
+  // Since these are really sitting on a plane, we can use the method
+  // proposed by the CORE offline. 
+  // 
+  // Return:
+  //    @c true on success, @c false otherwise.
+  //
+
   // Do the FMD1 stuff
   Double_t rot[9], trans[3];
   if (!GetFMD1Plane(rot, trans)) return kFALSE;
   // const char* path = "/ALIC_1/F1MT_1/FMD1_lid_0";
-  
+
+#if 0  
   // TGeoHMatrix delta;
-  TGeoTranslation global(0,0,324.670);
+  Double_t gRot[9], gTrans[3];
+  TVector3 ocb(-127, -220, 324.67);
+  TVector3 oct(-127, +220, 324.67);
+  TVector3 icb(+127, -220, 324.67);
+  TVector3 ict(+127, +220, 324.67);
+  if (!CalculatePlane(ocb, icb, oct, 0, gTrans, gRot)) { 
+    Warning("DoFMD1", "Failed to make reference plane");
+    return kFALSE;
+  }
+  PrintRotation("FMD1 ref rotation:",  gRot);
+  PrintVector("FMD1 ref translation:", gTrans);
+  TGeoRotation ggRot; ggRot.SetMatrix(gRot);
+  TGeoCombiTrans global(gTrans[0], gTrans[1], gTrans[2], &ggRot);
+#endif
+  Double_t off = getFMD1Offset();
+  Info("DoFMD1", "Lid offset is %f", off);
+
+  TGeoTranslation global(0,0,324.670-off);
   if (!MakeDelta(&global, rot, trans, fFMD1Delta)) 
     return kFALSE;
   
@@ -300,6 +461,18 @@ AliFMDSurveyToAlignObjs::DoFMD1()
 Bool_t
 AliFMDSurveyToAlignObjs::GetFMD2Plane(Double_t* rot, Double_t* trans) const
 {
+  // 
+  // Get the surveyed plane corresponding to the backside of FMD2.
+  // The plane is done as a best fit of the plane equation to at least
+  // 4 of the available survey points.
+  // 
+  // Parameters:
+  //    rot    Rotation matrix (direction cosines)
+  //    trans  Translation vector.
+  // 
+  // Return:
+  //    @c true on success, @c false otherwise
+  //
 
   // The possible survey points 
   const char*    names[] = { "FMD2_ITOP",  "FMD2_OTOP", 
@@ -312,12 +485,23 @@ AliFMDSurveyToAlignObjs::GetFMD2Plane(Double_t* rot, Double_t* trans) const
   TObjArray errors;
   
   // Loop and fill graph 
+  int i = 0;
   while (*name) {
     TVector3 p, e;
-    if (!GetPoint(*name++, p, e)) continue;
+    if (!GetPoint(*name, p, e)) {
+      name++;
+      i++;
+      continue;
+    }
     
+    if (i == 5) {
+      Warning("GetFMD2plane", "Setting error on %d, %s to 0.4", i, *name);
+      e.SetXYZ(0.4, 0.4, 0.4); // OBOT
+    }
     points.Add(new TVector3(p));
     errors.Add(new TVector3(e));
+    name++;
+    i++;
   }
   if (points.GetEntries() < 4) { 
     AliWarning(Form("Only got %d survey points - no good for FMD2 plane",
@@ -333,17 +517,43 @@ AliFMDSurveyToAlignObjs::GetFMD2Plane(Double_t* rot, Double_t* trans) const
 Bool_t
 AliFMDSurveyToAlignObjs::DoFMD2()
 {
+  // 
+  // Do the FMD2 calculations.  We have 6 survey points of which only
+  // 5 are normally surveyed.  These are all sittings 
+  //
+  //  - FMD2_ITOP   - In-side, top
+  //  - FMD2_IBOTM  - In-side, middle bottom
+  //  - FMD2_IBOT   - In-side, bottom
+  //  - FMD2_OTOP   - Out-side, top
+  //  - FMD2_OBOTM  - Out-side, middle bottom
+  //  - FMD2_OBOT   - Out-side, bottom
+  //
+  // The nominal coordinates of these retro-fitted survey stickers
+  // isn't known.  Also, these stickers are put on a thin (0.3mm
+  // thick) carbon cover which flexes quite easily.  This means, that
+  // to rotations and xy-translation obtained from the survey data
+  // cannot be used, and left is only the z-translation.
+  //
+  // Further more, since FMD2 to is attached to the ITS SPD thermal
+  // screen, it is questionable if the FMD2 survey will ever be used. 
+  // 
+  // Return:
+  //    @c true on success, @c false otherwise.
+  //
+
   // Do the FMD2 stuff
   Double_t rot[9], trans[3];
   if (!GetFMD2Plane(rot, trans)) return kFALSE;
-  // PrintRotation("FMD2 rotation:",  rot);
-  // PrintVector("FMD2 translation:", trans);
+  PrintRotation("FMD2 rotation:",  rot);
+  PrintVector("FMD2 translation:", trans);
 
+#if 0
   for (int i = 0; i < 3; i++) { 
     for (int j = 0; j < 3; j++) { 
       rot[i*3+j] = (i == j ? 1 : 0);
     }
   }
+#endif
   trans[0] = trans[1] = 0;
   trans[2] += 0.015;
   // PrintRotation("FMD2 rotation:",  rot);
@@ -363,6 +573,11 @@ AliFMDSurveyToAlignObjs::DoFMD2()
 void
 AliFMDSurveyToAlignObjs::Run()
 {
+  // 
+  // Run the task.
+  // 
+  //  
+
   AliFMDGeometry* geom = AliFMDGeometry::Instance();
   geom->Init();
   geom->InitTransformations();
@@ -371,20 +586,136 @@ AliFMDSurveyToAlignObjs::Run()
   DoFMD2();
 }
 
+//____________________________________________________________________
+void
+AliFMDSurveyToAlignObjs::Run(const char** files)
+{
+  // 
+  // Run the task.
+  // 
+  //  
+
+  AliFMDGeometry* geom = AliFMDGeometry::Instance();
+  geom->Init();
+  geom->InitTransformations();
+
+  const char** file = files; 
+  while (*file) { 
+    if ((*file)[0] == '\0') { 
+      Warning("Run", "no file specified");
+      file++;
+      continue;
+    }
+    if (!LoadSurveyFromLocalFile(*file)) { 
+      Warning("Run", "Failed to load %s", *file);
+      file++;
+      continue;
+    }
+    TString sDet(fSurveyObj->GetDetector());
+    Int_t   d    = Int_t(sDet[sDet.Length()-1] - '0');
+    Info("Run", "Making alignment for %s (%d)", sDet.Data(), d);
+    Bool_t ret = true;
+    switch (d) { 
+    case 1: ret = DoFMD1(); break;
+    case 2: ret = DoFMD2(); break;
+    default: 
+      Warning("Run", "Do not know how to deal with %s", sDet.Data());
+      break;
+    }
+    if (!ret) { 
+      Warning("Run", "Calculation for %s failed", sDet.Data());
+    }
+    file++;
+  }
+  CreateAlignObjs();
+  GetAlignObjArray()->Print();
+  FillDefaultAlignObjs();
+}
+
+//____________________________________________________________________
+AliAlignObjParams*
+AliFMDSurveyToAlignObjs::CreateDefaultAlignObj(const TString& path, 
+                                              Int_t id)
+{
+  Int_t nAlign = fAlignObjArray->GetEntries();
+  AliAlignObjParams* obj = 
+    new ((*fAlignObjArray)[nAlign]) AliAlignObjParams(path.Data(),
+                                                     id,0,0,0,0,0,0,kTRUE);
+  if (!obj) {
+    AliError(Form("Failed to create alignment object for %s", path.Data()));
+    return 0;
+  }
+  if (!obj->SetLocalPars(0, 0, 0, 0, 0, 0)) {
+    AliError(Form("Failed to set local transforms on %s", path.Data()));
+    return obj;
+  }
+  return obj;
+}
+
+//____________________________________________________________________
+AliAlignObjParams*
+AliFMDSurveyToAlignObjs::FindAlignObj(const TString& path) const 
+{
+  AliAlignObjParams* p = 0;
+  for (int i = 0; i < fAlignObjArray->GetEntries(); i++) { 
+    p = static_cast<AliAlignObjParams*>(fAlignObjArray->At(i));
+    if (path.EqualTo(p->GetSymName())) return p;
+  }
+  return 0;
+}
+
+//____________________________________________________________________
+Bool_t 
+AliFMDSurveyToAlignObjs::FillDefaultAlignObjs()
+{
+  for (int d = 1; d <= 3; d++) { 
+    const char sides[] = { 'T', 'B', 0 };
+    const char* side   = sides;
+    while (*side) { 
+      TString path = TString::Format("FMD/FMD%d_%c", d, *side);
+      if (!FindAlignObj(path)) CreateDefaultAlignObj(path, 0);
+
+      const char halves[] = { 'I', d == 1 ? '\0' : 'O', 0 };
+      const char*  half = halves;
+      while (*half) { 
+       int nsec  = *half == 'I' ? 10 : 20;
+       int start = *side == 'T' ? 0      : nsec/2;
+       int end   = *side == 'T' ? nsec/2 : nsec;
+       for (int s=start; s < end; s++) {
+         path = TString::Format("FMD/FMD%d_%c/FMD%c_%02d", 
+                                d, *side, *half, s);
+         CreateDefaultAlignObj(path, 0);
+       }
+       half++;
+      }
+      side++;
+    }
+    
+  }
+  return true;
+}
+
 //____________________________________________________________________
 Bool_t 
 AliFMDSurveyToAlignObjs::CreateAlignObjs()
 {
+  // 
+  // 
+  // Method to create the alignment objects
+  // 
+  // Return:
+  //    @c true on success, @c false otherwise
+  //  
   TClonesArray& array = *fAlignObjArray;
   Int_t         n     = array.GetEntriesFast();
 
   if (!fFMD1Delta.IsIdentity()) { 
-    new (array[n++]) AliAlignObjParams("FMD1/FMD1_T", 0, fFMD1Delta, kTRUE);
-    new (array[n++]) AliAlignObjParams("FMD1/FMD1_B", 0, fFMD1Delta, kTRUE);
+    new (array[n++]) AliAlignObjParams("FMD/FMD1_T", 0, fFMD1Delta, kTRUE);
+    new (array[n++]) AliAlignObjParams("FMD/FMD1_B", 0, fFMD1Delta, kTRUE);
   }
   if (!fFMD2Delta.IsIdentity()) { 
-    new (array[n++]) AliAlignObjParams("FMD2/FMD2_T", 0, fFMD2Delta, kTRUE);
-    new (array[n++]) AliAlignObjParams("FMD2/FMD2_B", 0, fFMD2Delta, kTRUE);
+    new (array[n++]) AliAlignObjParams("FMD/FMD2_T", 0, fFMD2Delta, kTRUE);
+    new (array[n++]) AliAlignObjParams("FMD/FMD2_B", 0, fFMD2Delta, kTRUE);
   }
   // array.Print();
   
@@ -395,7 +726,13 @@ AliFMDSurveyToAlignObjs::CreateAlignObjs()
 void 
 AliFMDSurveyToAlignObjs::PrintVector(const char* text, const TVector3& v)
 {
-  // Print a vector
+  // 
+  // Service member function to print a vector
+  // 
+  // Parameters:
+  //    text Prefix text
+  //    v    Vector
+  //
   Double_t va[] = { v.X(), v.Y(), v.Z() };
   PrintVector(text, va);
 }
@@ -403,7 +740,13 @@ AliFMDSurveyToAlignObjs::PrintVector(const char* text, const TVector3& v)
 void 
 AliFMDSurveyToAlignObjs::PrintVector(const char* text, const Double_t* v)
 {
-  // Print a vector
+  // 
+  // Service member function to print a vector
+  // 
+  // Parameters:
+  //    text Prefix text
+  //    v    Vector (array of 3 doubles)
+  //
   std::cout << text 
            << std::setw(15) << v[0] 
            << std::setw(15) << v[1]
@@ -416,7 +759,14 @@ AliFMDSurveyToAlignObjs::PrintVector(const char* text, const Double_t* v)
 void 
 AliFMDSurveyToAlignObjs::PrintRotation(const char* text, const Double_t* rot)
 {
-  // Print a rotation matrix
+  // 
+  // Service member function to print a rotation matrix
+  // 
+  // Parameters:
+  //    text Prefix text
+  //    v    Matrix (array of 9 doubles)
+  //
+
   std::cout << text << std::endl;
   for (size_t i = 0; i < 3; i++) { 
     for (size_t j = 0; j < 3; j++)