DDCore/src/AlignmentsCalculator.cpp

0001 //==========================================================================
0002 //  AIDA Detector description implementation
0003 //--------------------------------------------------------------------------
0004 // Copyright (C) Organisation europeenne pour la Recherche nucleaire (CERN)
0005 // All rights reserved.
0006 //
0007 // For the licensing terms see $DD4hepINSTALL/LICENSE.
0008 // For the list of contributors see $DD4hepINSTALL/doc/CREDITS.
0009 //
0010 // Author     : M.Frank
0011 //
0012 //==========================================================================
0013
0014 // Framework include files
0015 #include <DD4hep/Printout.h>
0016 #include <DD4hep/Conditions.h>
0017 #include <DD4hep/ConditionsMap.h>
0018 #include <DD4hep/InstanceCount.h>
0019 #include <DD4hep/MatrixHelpers.h>
0020 #include <DD4hep/ConditionDerived.h>
0021 #include <DD4hep/DetectorProcessor.h>
0022 #include <DD4hep/AlignmentsProcessor.h>
0023 #include <DD4hep/AlignmentsCalculator.h>
0024 #include <DD4hep/detail/AlignmentsInterna.h>
0025
0026 using namespace dd4hep;
0027 using namespace dd4hep::align;
0028 using Result = AlignmentsCalculator::Result;
0029
0030 /// Namespace for the AIDA detector description toolkit
0031 namespace dd4hep {
0032
0033   /// Namespace for the alignment part of the AIDA detector description toolkit
0034   namespace align {
0035
0036     /// Anonymous implementation classes
0037     namespace {
0038       static Delta        identity_delta;
0039
0040       /// Alignment calculator.
0041       /**
0042        *  Uses internally the conditions mechanism to calculator the alignment conditions.
0043        *
0044        *  \author  M.Frank
0045        *  \version 2.0
0046        *  \ingroup DD4HEP_ALIGNMENTS
0047        */
0048       class Calculator {
0049       public:
0050         class Entry;
0051         class Context;
0052
0053       public:
0054         /// Initializing constructor
0055         Calculator() = default;
0056         /// Default destructor
0057         ~Calculator() = default;
0058         /// Compute all alignment conditions of the lower levels
0059         Result compute(Context& context, Entry& entry) const;
0060         /// Resolve child dependencies for a given context
0061         void resolve(Context& context, DetElement child) const;
0062       };
0063
0064       class Calculator::Entry  {
0065       public:
0066         DetElement::Object*         det   = 0;
0067         const Delta*                delta = 0;
0068         AlignmentCondition::Object* cond  = 0;
0069         unsigned char               key   = 0, valid = 0, created = 0, _pad[1] { 0 };
0070         Entry(DetElement d, const Delta* del) : det(d.ptr()), delta(del), key(d.key())  {}
0071       };
0072
0073       class Calculator::Context  {
0074       public:
0075         typedef std::map<DetElement,std::size_t,AlignmentsCalculator::PathOrdering>  DetectorMap;
0076         typedef std::map<unsigned int,std::size_t>             Keys;
0077         typedef std::vector<Entry>                        Entries;
0078
0079         DetectorMap    detectors;
0080         Keys           keys;
0081         Entries        entries;
0082         ConditionsMap& mapping;
0083         Context(ConditionsMap& m) : mapping(m)  {
0084           InstanceCount::increment(this);
0085         }
0086         ~Context()  {
0087           InstanceCount::decrement(this);
0088         }
0089         void insert(DetElement det, const Delta* delta)   {
0090           if ( det.isValid() )  {
0091             Entry entry(det,delta);
0092             detectors.emplace(det, entries.size());
0093             keys.emplace(entry.key, entries.size());
0094             entries.emplace_back(entry);
0095             return;
0096           }
0097           except("AlignContext","Failed to add entry: invalid detector handle!");
0098         }
0099       };
0100     }
0101   }       /* End namespace align */
0102 }         /* End namespace dd4hep     */
0103
0104 //static PrintLevel s_PRINT = DEBUG;
0105 //static PrintLevel s_PRINT = INFO;
0106 static PrintLevel s_PRINT = WARNING;
0107
0108 /// Callback to output alignments information
0109 int AlignmentsCalculator::Scanner::operator()(DetElement de, int)  const  {
0110   if ( de.isValid() )  {
0111     Condition::key_type key(ConditionKey::KeyMaker(de.key(),align::Keys::deltaKey).hash);
0112     Condition c = context.condition(key, false);
0113     if ( c.isValid() )  {
0114       const Delta& d = c.get<Delta>();
0115       deltas.emplace(de,&d);
0116       if ( iov ) iov->iov_intersection(c->iov->key());
0117     }
0118     return 1;
0119   }
0120   return 0;
0121 }
0122
0123 /// Compute all alignment conditions of the lower levels
0124 Result Calculator::compute(Context& context, Entry& e)   const  {
0125   Result result;
0126   DetElement det = e.det;
0127
0128   if ( e.valid == 1 )  {
0129     printout(DEBUG,"ComputeAlignment","================ IGNORE %s (already valid)",det.path().c_str());
0130     return result;
0131   }
0132   AlignmentCondition c = context.mapping.get(det, Keys::alignmentKey);
0133   AlignmentCondition cond = c.isValid() ? c : AlignmentCondition(det.path()+"#alignment");
0134   AlignmentData&     align = cond.data();
0135   const Delta*       delta = e.delta ? e.delta : &identity_delta;
0136   TGeoHMatrix        transform_for_delta;
0137
0138   printout(DEBUG,"ComputeAlignment",
0139            "============================== Compute transformation of %s",det.path().c_str());
0140   e.valid     = 1;
0141   e.cond      = cond.ptr();
0142   align.delta = *delta;
0143   delta->computeMatrix(transform_for_delta);
0144   result.multiply += 2;
0145
0146   DetElement parent_det = det.parent();
0147   AlignmentCondition parent_cond = context.mapping.get(parent_det, Keys::alignmentKey);
0148   TGeoHMatrix parent_transform;
0149   if (parent_cond.isValid()) {
0150     AlignmentData&     parent_align = parent_cond.data();
0151     parent_transform = parent_align.worldTrafo;
0152   }
0153   else if ( det.parent().isValid() )   {
0154     parent_transform = parent_det.nominal().worldTransformation();
0155   }
0156   else {
0157     // The tranformation from the "world" to its parent is non-existing i.e. unity
0158   }
0159
0160   align.detectorTrafo = det.nominal().detectorTransformation() * transform_for_delta;
0161   align.worldTrafo    = parent_transform * align.detectorTrafo;
0162   align.trToWorld     = detail::matrix::_transform(&align.worldTrafo);
0163   ++result.computed;
0164   result.multiply += 3;
0165   // Update mapping if the condition is freshly created
0166   if ( !c.isValid() )  {
0167     e.created = 1;
0168     cond->flags |= Condition::ALIGNMENT_DERIVED;
0169     cond->hash = ConditionKey(e.det,Keys::alignmentKey).hash;
0170     context.mapping.insert(e.det, Keys::alignmentKey, cond);
0171   }
0172   if ( s_PRINT <= INFO )  {
0173     printout(INFO,"ComputeAlignment","Level:%d Path:%s DetKey:%08X: Cond:%s key:%16llX",
0174              det.level(), det.path().c_str(), det.key(),
0175              yes_no(e.delta != 0), (long long int)cond.key());
0176     if ( s_PRINT <= DEBUG )  {
0177       ::printf("Nominal:     '%s' ", det.path().c_str());
0178       det.nominal().worldTransformation().Print();
0179       ::printf("Parent: '%s' -> '%s' ", det.path().c_str(), parent_det.path().c_str());
0180       parent_transform.Print();
0181       ::printf("DetectorTrafo: '%s' -> '%s' ", det.path().c_str(), det.parent().path().c_str());
0182       det.nominal().detectorTransformation().Print();
0183       ::printf("Delta:       '%s' ", det.path().c_str());
0184       transform_for_delta.Print();
0185       ::printf("Result:      '%s' ", det.path().c_str());
0186       align.worldTrafo.Print();
0187     }
0188   }
0189   return result;
0190 }
0191
0192 /// Resolve child dependencies for a given context
0193 void Calculator::resolve(Context& context, DetElement detector) const   {
0194   auto children = detector.children();
0195   auto item = context.detectors.find(detector);
0196   if ( item == context.detectors.end() ) context.insert(detector,0);
0197   for(const auto& c : children )
0198     resolve(context, c.second);
0199 }
0200
0201 /// Optimized call using already properly ordered Deltas
0202 Result AlignmentsCalculator::compute(const OrderedDeltas& deltas,
0203                                      ConditionsMap& alignments)  const
0204 {
0205   Result  result;
0206   Calculator obj;
0207   Calculator::Context context(alignments);
0208   for( const auto& i : deltas )
0209     context.insert(i.first, i.second);
0210   for( const auto& i : deltas )
0211     obj.resolve(context,i.first);
0212   for( auto& i : context.entries )
0213     result += obj.compute(context, i);
0214   return result;
0215 }
0216
0217 /// Compute all alignment conditions of the internal dependency list
0218 Result AlignmentsCalculator::compute(const std::map<DetElement, Delta>& deltas,
0219                                      ConditionsMap& alignments)  const
0220 {
0221   Calculator::Context context(alignments);
0222   // This is a tricky one. We absolutely need the detector elements ordered
0223   // by their depth aka. the distance to /world.
0224   // Unfortunately one cannot use the raw pointer of the DetElement here,
0225   // But has to insert them in a map which is ordered by the DetElement path.
0226   //
0227   // Otherwise memory randomization gives us the wrong order and the
0228   // corrections are calculated in the wrong order ie. not top -> down the
0229   // hierarchy, but in "some" order depending on the pointer values!
0230   //
0231   OrderedDeltas ordered_deltas;
0232   for( const auto& i : deltas )
0233     ordered_deltas.emplace(i.first, &i.second);
0234   return compute(ordered_deltas, alignments);
0235 }
0236
0237 /// Compute all alignment conditions of the internal dependency list
0238 Result AlignmentsCalculator::compute(const std::map<DetElement, const Delta*>& deltas,
0239                                      ConditionsMap& alignments)  const
0240 {
0241   Calculator::Context context(alignments);
0242   // This is a tricky one. We absolutely need the detector elements ordered
0243   // by their depth aka. the distance to /world.
0244   // Unfortunately one cannot use the raw pointer of the DetElement here,
0245   // But has to insert them in a map which is ordered by the DetElement path.
0246   //
0247   // Otherwise memory randomization gives us the wrong order and the
0248   // corrections are calculated in the wrong order ie. not top -> down the
0249   // hierarchy, but in "some" order depending on the pointer values!
0250   //
0251   OrderedDeltas ordered_deltas;
0252   for( const auto& i : deltas )
0253     ordered_deltas.insert(i);
0254   return compute(ordered_deltas, alignments);
0255 }
0256
0257 /// Helper: Extract all Delta-conditions from the conditions map
0258 size_t AlignmentsCalculator::extract_deltas(cond::ConditionUpdateContext& ctxt,
0259                                             ExtractContext& extract_context,
0260                                             OrderedDeltas& deltas,
0261                                             IOV* effective_iov)   const
0262 {
0263   return extract_deltas(ctxt.world(), ctxt, extract_context, deltas, effective_iov);
0264 }
0265
0266 /// Helper: Extract all Delta-conditions from the conditions map starting at a certain sub-tree
0267 size_t AlignmentsCalculator::extract_deltas(DetElement start,
0268                                             cond::ConditionUpdateContext& ctxt,
0269                                             ExtractContext& extract_context,
0270                                             OrderedDeltas& deltas,
0271                                             IOV* effective_iov)   const
0272 {
0273   if ( !extract_context.empty() )   {
0274     struct DeltaScanner : public Condition::Processor   {
0275       OrderedDeltas& delta_conditions;
0276       ExtractContext& extract_context;
0277       IOV* effective_iov = 0;
0278       /// Constructor
0279       DeltaScanner(OrderedDeltas& d, ExtractContext& e, IOV* eff_iov)
0280         : delta_conditions(d), extract_context(e), effective_iov(eff_iov) {}
0281       /// Conditions callback for object processing
0282       virtual int process(Condition c)  const override  {
0283         ConditionKey::KeyMaker key_maker(c->hash);
0284         if ( key_maker.values.item_key == align::Keys::deltaKey )   {
0285           auto idd = extract_context.find(key_maker.values.det_key);
0286           if ( idd != extract_context.end() )   {
0287             const Delta& d  = c.get<Delta>();
0288             DetElement   de = idd->second;
0289             delta_conditions.emplace(de,&d);
0290             if (effective_iov) effective_iov->iov_intersection(c->iov->key());
0291             return 1;
0292           }
0293           // Fatal or not? Depends if the context should be re-usable or the projection
0294           //except("extract_deltas","++ Inconsistency between extraction context and conditions content!!");
0295         }
0296         return 0;
0297       }
0298     };
0299     DeltaScanner scanner(deltas,extract_context,effective_iov);
0300     ctxt.resolver->conditionsMap().scan(scanner);
0301     return deltas.size();
0302   }
0303   DetectorScanner().scan(AlignmentsCalculator::Scanner(ctxt,deltas,effective_iov),start);
0304   for( const auto& d : deltas ) extract_context.emplace(d.first.key(), d.first);
0305   return deltas.size();
0306 }
0307
0308 /// Helper: Extract all Delta-conditions from the conditions map
0309 size_t AlignmentsCalculator::extract_deltas(cond::ConditionUpdateContext& ctxt,
0310                                             OrderedDeltas& deltas,
0311                                             IOV* effective_iov)   const
0312 {
0313   return extract_deltas(ctxt.world(), ctxt, deltas, effective_iov);
0314 }
0315
0316
0317 /// Helper: Extract all Delta-conditions from the conditions map
0318 size_t AlignmentsCalculator::extract_deltas(DetElement start,
0319                                             cond::ConditionUpdateContext& ctxt,
0320                                             OrderedDeltas& deltas,
0321                                             IOV* effective_iov)   const
0322 {
0323   DetectorScanner().scan(AlignmentsCalculator::Scanner(ctxt,deltas,effective_iov),start);
0324   return deltas.size();
0325 }
0326
0327 #include <DD4hep/GrammarUnparsed.h>
0328 static auto s_registry = GrammarRegistry::pre_note<AlignmentsCalculator::OrderedDeltas>(1);