Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 13:08:05

0001 #ifndef GEOMETRY_ECALGEOMETRYLOADER_ICC
0002 #define GEOMETRY_ECALGEOMETRYLOADER_ICC 1
0003 
0004 #include "Geometry/EcalTestBeam/test/ee/CaloGeometryLoaderTest.h"
0005 #include "Geometry/CaloGeometry/interface/CaloCellGeometry.h"
0006 #include "DetectorDescription/Core/interface/DDCompactView.h"
0007 #include "DetectorDescription/Core/interface/DDSolid.h"
0008 #include "DetectorDescription/Core/interface/DDSpecifics.h"
0009 
0010 #include "CLHEP/Units/GlobalSystemOfUnits.h"
0011 #include "CLHEP/Geometry/Plane3D.h"
0012 
0013 #include <vector>
0014 #include <fstream>
0015 
0016 typedef CaloCellGeometry::CCGFloat CCGFloat;
0017 
0018 template <class T>
0019 const double CaloGeometryLoaderTest<T>::k_ScaleFromDDDtoGeant(0.1);
0020 
0021 template <class T>
0022 CaloGeometryLoaderTest<T>::CaloGeometryLoaderTest() {
0023   m_filter.setCriteria(DDValue("Volume",
0024                                "EndcapSC",
0025                                //                 "EESCEnv1",
0026                                0),
0027                        DDCompOp::equals);
0028 }
0029 
0030 template <class T>
0031 typename CaloGeometryLoaderTest<T>::PtrType CaloGeometryLoaderTest<T>::load(const DDCompactView* cpv,
0032                                                                             const Alignments* alignments,
0033                                                                             const Alignments* globals) {
0034   PtrType geom = std::make_unique<T>();
0035 
0036   makeGeometry(cpv, dynamic_cast<T*>(geom.get()), alignments, globals);
0037 
0038   return geom;
0039 }
0040 
0041 template <class T>
0042 void CaloGeometryLoaderTest<T>::makeGeometry(const DDCompactView* cpv,
0043                                              T* geom,
0044                                              const Alignments* alignments,
0045                                              const Alignments* globals) {
0046   geom->allocateCorners(T::k_NumberOfCellsForCorners);
0047   geom->allocatePar(T::k_NumberOfParametersPerShape * T::k_NumberOfShapes, T::k_NumberOfParametersPerShape);
0048 
0049   DDFilteredView fv(*cpv, m_filter);
0050 
0051   unsigned int counter(0);
0052   for (bool doSubDets = fv.firstChild(); doSubDets; doSubDets = fv.nextSibling()) {
0053     ++counter;
0054     const DDSolid& solid(fv.logicalPart().solid());
0055 
0056     const ParmVec& parameters(solid.parameters());
0057 
0058     DD3Vector x, y, z;
0059     fv.rotation().GetComponents(x, y, z);
0060     const CLHEP::HepRep3x3 temp(x.X(), y.X(), z.X(), x.Y(), y.Y(), z.Y(), x.Z(), y.Z(), z.Z());
0061     const CLHEP::HepRotation hr(temp);
0062     const CLHEP::Hep3Vector h3v(fv.translation().X(), fv.translation().Y(), fv.translation().Z());
0063     const HepGeom::Transform3D ht3d(hr,  // only scale translation
0064                                     k_ScaleFromDDDtoGeant * h3v);
0065 
0066     const unsigned int id(getDetIdForDDDNode(fv));
0067 
0068     myFillGeom(geom, parameters, ht3d, id);
0069   }
0070 }
0071 
0072 template <class T>
0073 void CaloGeometryLoaderTest<T>::myFillGeom(T* geom,
0074                                            const ParmVec& vv,
0075                                            const HepGeom::Transform3D& tr,
0076                                            const unsigned int id) {
0077   static const double kRadToDeg(180 / M_PI);
0078   static const HepGeom::Point3D<double> PIVOT(13533.35, -000.02, 0);        //-152.23 ) ;
0079   static const HepGeom::Point3D<double> FLOOR(13533.35, -000.02, -152.23);  //-152.23 ) ;
0080 
0081   static const HepGeom::Point3D<double> IP(0, 0, 0);  //-152.23 ) ;
0082   HepGeom::Point3D<double> AA0;
0083   HepGeom::Point3D<double> BB0;
0084   HepGeom::Point3D<double> CC0;
0085   HepGeom::Point3D<double> DD0;
0086   HepGeom::Point3D<double> EE0;
0087   HepGeom::Point3D<double> FF0;
0088   HepGeom::Point3D<double> GG0;
0089   HepGeom::Point3D<double> HH0;
0090 
0091   std::vector<CCGFloat> pv;
0092   pv.reserve(vv.size());
0093   for (unsigned int i(0); i != vv.size(); ++i) {
0094     const double factor(1 == i || 2 == i || 6 == i || 10 == i ? 1 : k_ScaleFromDDDtoGeant);
0095     pv.push_back(factor * vv[i]);
0096   }
0097 
0098   std::vector<GlobalPoint> corners(8);
0099 
0100   TruncatedPyramid::createCorners(pv, tr, corners);
0101 
0102   if ((160 == id || 164 == id || 134 == id || 133 == id || 153 == id || 150 == id || 130 == id) &&
0103       (0 < corners[0].x() && 0 < corners[0].y() && 0 < corners[0].z())) {
0104     /*      std::cout<<"\n *************** id = "<<id
0105            <<"\n "<<corners[0]
0106            <<"\n "<<corners[1]
0107            <<"\n "<<corners[2]
0108            <<"\n "<<corners[3]
0109            <<"\n "<<corners[4]
0110            <<"\n "<<corners[5]
0111            <<"\n "<<corners[6]
0112            <<"\n "<<corners[7]
0113            <<std::endl ;*/
0114 
0115     if (160 == id)
0116       AA0 = HepGeom::Point3D<double>(corners[3].x(), corners[3].y(), corners[3].z());
0117     if (164 == id)
0118       BB0 = HepGeom::Point3D<double>(corners[2].x(), corners[2].y(), corners[2].z());
0119     if (134 == id)
0120       CC0 = HepGeom::Point3D<double>(corners[1].x(), corners[1].y(), corners[1].z());
0121     if (130 == id)
0122       DD0 = HepGeom::Point3D<double>(corners[0].x(), corners[0].y(), corners[0].z());
0123     if (150 == id)
0124       EE0 = HepGeom::Point3D<double>(corners[2].x(), corners[2].y(), corners[2].z());
0125     if (153 == id)
0126       FF0 = HepGeom::Point3D<double>(corners[2].x(), corners[2].y(), corners[2].z());
0127     if (133 == id)
0128       GG0 = HepGeom::Point3D<double>(corners[2].x(), corners[2].y(), corners[2].z());
0129     if (130 == id)
0130       HH0 = HepGeom::Point3D<double>(corners[2].x(), corners[2].y(), corners[2].z());
0131   }
0132   static bool done(false);
0133   if (0.001 < AA0.mag() && 0.001 < BB0.mag() && 0.001 < CC0.mag() && 0.001 < DD0.mag() && !done) {
0134     done = true;
0135 
0136     const HepGeom::Transform3D H4ToBeamLine(PIVOT + HepGeom::Point3D<double>(0, 1, 0),
0137                                             PIVOT + HepGeom::Point3D<double>(0, 0, 1),
0138                                             PIVOT + HepGeom::Point3D<double>(1, 0, 0),
0139                                             HepGeom::Point3D<double>(-1, 0, 0),
0140                                             HepGeom::Point3D<double>(0, -1, 0),
0141                                             HepGeom::Point3D<double>(0, 0, 1));
0142     /*                    HepGeom::Point3D<double> ( 0, -1, 0),
0143                       HepGeom::Point3D<double> ( 1,  0, 0),
0144                       HepGeom::Point3D<double> ( 0,  0, 1) ) ;*/
0145 
0146     std::cout << "\nH4ToBeamLine"
0147               << "   Theta_x = " << H4ToBeamLine.getRotation().thetaX() * kRadToDeg
0148               << ",  Theta_y = " << H4ToBeamLine.getRotation().thetaY() * kRadToDeg
0149               << ",  Theta_z = " << H4ToBeamLine.getRotation().thetaZ() * kRadToDeg
0150               << ",  eta = " << -log(tan(H4ToBeamLine.getRotation().thetaZ() / 2.)) << std::endl;
0151 
0152     std::cout << "\nH4ToBeamLine"
0153               << "   Phi_x = " << H4ToBeamLine.getRotation().phiX() * kRadToDeg
0154               << ",  Phi_y = " << H4ToBeamLine.getRotation().phiY() * kRadToDeg
0155               << ",  Phi_z = " << H4ToBeamLine.getRotation().phiZ() * kRadToDeg << std::endl;
0156 
0157     std::cout << "pivot to beamline=" << H4ToBeamLine * PIVOT << std::endl;
0158     std::cout << "endcap to beamline=" << H4ToBeamLine * HepGeom::Point3D<double>(PIVOT.x() + 316, 0, 50) << std::endl;
0159 
0160     std::cout << "\nA0 = " << AA0 << "\nB0 = " << BB0 << "\nC0 = " << CC0 << "\nD0 = " << DD0 << std::endl;
0161 
0162     const HepGeom::Point3D<double>& A0(AA0);
0163     const HepGeom::Point3D<double>& B0(BB0);
0164     const HepGeom::Point3D<double>& C0(CC0);
0165     const HepGeom::Point3D<double>& D0(DD0);
0166 
0167     const double angle((C0 - B0).angle(A0 - B0));
0168 
0169     static const std::vector<HepGeom::Point3D<double> > extra = {HepGeom::Point3D<double>(13859.05, 005.79, -002.66),
0170                                                                  HepGeom::Point3D<double>(13860.98, -038.39, -002.11),
0171                                                                  HepGeom::Point3D<double>(13853.19, -039.25, -030.77),
0172                                                                  HepGeom::Point3D<double>(13851.23, 004.99, -031.39)};
0173 
0174     static const std::vector<HepGeom::Point3D<double> > sur = {
0175         HepGeom::Point3D<double>(13862.30, 20.68, 11.60),    HepGeom::Point3D<double>(13865.69, -52.85, 12.64),
0176         HepGeom::Point3D<double>(13851.53, -54.27, -44.75),  HepGeom::Point3D<double>(13848.21, 19.18, -45.89),
0177 
0178         HepGeom::Point3D<double>(13861.72, 28.62, 11.61),    HepGeom::Point3D<double>(13866.88, -44.81, 12.64),
0179         HepGeom::Point3D<double>(13852.76, -46.56, -44.73),  HepGeom::Point3D<double>(13847.67, 26.79, -45.88),
0180         HepGeom::Point3D<double>(13856.82, 63.41, 11.60),    HepGeom::Point3D<double>(13869.76, -9.06, 12.64),
0181         HepGeom::Point3D<double>(13855.90, -12.31, -44.73),  HepGeom::Point3D<double>(13843.04, 60.08, -45.89),
0182         HepGeom::Point3D<double>(13860.86, -36.96, 11.60),   HepGeom::Point3D<double>(13851.38, -109.96, 12.63),
0183         HepGeom::Point3D<double>(13837.21, -108.89, -44.73), HepGeom::Point3D<double>(13846.73, -35.98, -45.88),
0184         HepGeom::Point3D<double>(13854.03, -76.21, 11.60),   HepGeom::Point3D<double>(13835.82, -147.53, 12.64),
0185         HepGeom::Point3D<double>(13821.88, -144.76, -44.73), HepGeom::Point3D<double>(13840.13, -73.54, -45.88),
0186         HepGeom::Point3D<double>(13861.10, 35.48, 09.88),    HepGeom::Point3D<double>(13867.78, -37.82, 11.07),
0187         HepGeom::Point3D<double>(13853.42, -40.01, -46.22),  HepGeom::Point3D<double>(13846.80, 33.21, -47.52),
0188         HepGeom::Point3D<double>(13860.68, 35.23, 19.89),    HepGeom::Point3D<double>(13867.35, -38.08, 20.23),
0189         HepGeom::Point3D<double>(13854.67, -39.43, -37.48),  HepGeom::Point3D<double>(13848.07, 33.79, -37.92),
0190         HepGeom::Point3D<double>(13859.95, 34.97, 29.87),    HepGeom::Point3D<double>(13866.62, -38.34, 29.33),
0191         HepGeom::Point3D<double>(13855.61, -38.88, -28.76),  HepGeom::Point3D<double>(13849.00, 34.34, -28.32),
0192         HepGeom::Point3D<double>(13858.93, 34.41, 39.81),    HepGeom::Point3D<double>(13865.65, -38.90, 38.40),
0193         HepGeom::Point3D<double>(13856.34, -38.58, -19.97),  HepGeom::Point3D<double>(13849.70, 34.65, -18.68),
0194         HepGeom::Point3D<double>(13857.67, 33.69, 49.76),    HepGeom::Point3D<double>(13864.43, -39.60, 47.48),
0195         HepGeom::Point3D<double>(13856.81, -38.43, -11.13),  HepGeom::Point3D<double>(13850.14, 34.78, -8.95)};
0196 
0197     double fsumx(0);
0198     double fsumy(0);
0199     double fsumz(0);
0200     double fsumx2(0);
0201     double fsumy2(0);
0202     double fsumz2(0);
0203 
0204     double psumx(0);
0205     double psumy(0);
0206     double psumz(0);
0207     double psumx2(0);
0208     double psumy2(0);
0209     double psumz2(0);
0210 
0211     double sumx(0);
0212     double sumy(0);
0213     double sumz(0);
0214     double sumx2(0);
0215     double sumy2(0);
0216     double sumz2(0);
0217     unsigned int nsur(sur.size() / 4);
0218 
0219     const std::string names[10] = {"H0", "H1", "H2", "H3", "H4", "V1", "V2", "V3", "V4", "V5"};
0220 
0221     HepGeom::Transform3D V3Transform;
0222 
0223     for (unsigned int i(0); i != nsur; ++i) {
0224       const HepGeom::Point3D<double> SA(sur[4 * i]);
0225       const HepGeom::Point3D<double> SB(sur[4 * i + 1]);
0226       const HepGeom::Point3D<double> SC(sur[4 * i + 2]);
0227       const HepGeom::Point3D<double> SD(sur[4 * i + 3]);
0228       /*     const HepGeom::Point3D<double>  SD ( sur[4*i] ) ;
0229      const HepGeom::Point3D<double>  SA ( sur[4*i+1] ) ;
0230      const HepGeom::Point3D<double>  SB ( sur[4*i+2] ) ;
0231      const HepGeom::Point3D<double>  SC ( sur[4*i+3] ) ;*/
0232 
0233       const HepGeom::Point3D<double>& aa(SA);
0234       const HepGeom::Point3D<double> bb(SA + (SB - SA).unit() * (B0 - A0).mag());
0235 
0236       const HepGeom::Plane3D<double> plane(aa, bb, SC);
0237 
0238       const HepGeom::Point3D<double> cc(
0239           bb +
0240           (HepGeom::Rotate3D(-angle, plane.normal()) * HepGeom::Vector3D<double>(aa - bb)).unit() * (C0 - B0).mag());
0241 
0242       const HepGeom::Transform3D trform(A0, B0, C0, aa, bb, cc);
0243 
0244       const HepGeom::Transform3D trforminv(aa, bb, cc, A0, B0, C0);
0245 
0246       const HepGeom::Point3D<double> dd(trform * D0);
0247 
0248       const HepGeom::Transform3D BeamLineToCMS(H4ToBeamLine * aa, H4ToBeamLine * bb, H4ToBeamLine * cc, A0, B0, C0);
0249 
0250       if (7 == i)
0251         V3Transform = BeamLineToCMS;
0252 
0253       const CLHEP::HepRotation blcmsrot(BeamLineToCMS.getRotation());
0254       const HepGeom::Vector3D<double> blcmstra(BeamLineToCMS.getTranslation());
0255 
0256       std::cout << "\n********* Name is " << names[i];
0257 
0258       std::cout << "\n\n B discrepancy is " << (bb - SB) << ", " << (bb - SB).mag() << "\n\n C discrepancy is "
0259                 << (cc - SC) << ", " << (cc - SC).mag() << "\n\n D discrepancy is " << (dd - SD) << ", "
0260                 << (dd - SD).mag() << std::endl;
0261       /*
0262      const double dxbb ( bb.x() - SB.x() ) ;
0263      const double dybb ( bb.y() - SB.y() ) ;
0264      const double dzbb ( bb.z() - SB.z() ) ;
0265 
0266      const double dxcc ( cc.x() - SC.x() ) ;
0267      const double dycc ( cc.y() - SC.y() ) ;
0268      const double dzcc ( cc.z() - SC.z() ) ;
0269 */
0270       const double dxdd(dd.x() - SD.x());
0271       const double dydd(dd.y() - SD.y());
0272       const double dzdd(dd.z() - SD.z());
0273 
0274       sumx += dxdd;          //dxbb + dxcc + dxdd ;
0275       sumy += dydd;          //dybb + dycc + dydd ;
0276       sumz += dzdd;          //dzbb + dzcc + dzdd ;
0277       sumx2 += dxdd * dxdd;  //dxbb*dxbb + dxcc*dxcc + dxdd*dxdd ;
0278       sumy2 += dydd * dydd;  //dybb*dybb + dycc*dycc + dydd*dydd ;
0279       sumz2 += dzdd * dzdd;  //dzbb*dzbb + dzcc*dzcc + dzdd*dzdd ;
0280                              //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
0281       if (0 == i) {
0282         const HepGeom::Point3D<double> ee(trform * EE0);
0283         const HepGeom::Point3D<double> ff(trform * FF0);
0284         const HepGeom::Point3D<double> gg(trform * GG0);
0285         const HepGeom::Point3D<double> hh(trform * HH0);
0286 
0287         std::cout << "\n\n E discrepancy is " << (ee - extra[0]) << ", " << (ee - extra[0]).mag()
0288                   << "\n\n F discrepancy is " << (ff - extra[1]) << ", " << (ff - extra[1]).mag()
0289                   << "\n\n G discrepancy is " << (gg - extra[2]) << ", " << (gg - extra[2]).mag()
0290                   << "\n\n H discrepancy is " << (hh - extra[3]) << ", " << (hh - extra[3]).mag() << std::endl;
0291 
0292         std::cout << "\n Point SA in CMS coord = " << trforminv * SA << std::endl;
0293       }
0294       const HepGeom::Point3D<double> cmsPiv(trforminv * PIVOT);
0295       psumx += cmsPiv.x();
0296       psumy += cmsPiv.y();
0297       psumz += cmsPiv.z();
0298       psumx2 += cmsPiv.x() * cmsPiv.x();
0299       psumy2 += cmsPiv.y() * cmsPiv.y();
0300       psumz2 += cmsPiv.z() * cmsPiv.z();
0301       const HepGeom::Point3D<double> cmsFlr(trforminv * FLOOR);
0302       fsumx += cmsFlr.x();
0303       fsumy += cmsFlr.y();
0304       fsumz += cmsFlr.z();
0305       fsumx2 += cmsFlr.x() * cmsFlr.x();
0306       fsumy2 += cmsFlr.y() * cmsFlr.y();
0307       fsumz2 += cmsFlr.z() * cmsFlr.z();
0308 
0309       const CLHEP::HepRotation rotInv(trforminv.getRotation());
0310 
0311       std::cout << "\n CMS IP in H4 coord = " << trform * IP << std::endl;
0312       std::cout << "\n Pivot point in CMS coord = " << cmsPiv << std::endl;
0313       std::cout << "\n Floor point in CMS coord = " << cmsFlr << std::endl;
0314       //     std::cout<<"\n\nTransform rotation is "<<rotInv
0315       //          <<"\n\n          translation="<<trforminv.getTranslation() <<std::endl ;
0316 
0317       /*     std::cout<<"\n"<<names[i]
0318           <<"   Theta_x == "<<(blcmsrot*HepGeom::Point3D<double> (1,0,0)).theta()*kRadToDeg
0319           <<",  Theta_y == "<<(blcmsrot*HepGeom::Point3D<double> (0,1,0)).theta()*kRadToDeg
0320           <<",  Theta_z == "<<(blcmsrot*HepGeom::Point3D<double> (0,0,1)).theta()*kRadToDeg<<std::endl ;
0321 
0322      std::cout<<"\n"<<names[i]
0323           <<"   Phi_x == "<<(blcmsrot*HepGeom::Point3D<double> (1,0,0)).phi()*kRadToDeg
0324           <<",  Phi_y == "<<(blcmsrot*HepGeom::Point3D<double> (0,1,0)).phi()*kRadToDeg
0325           <<",  Phi_z == "<<(blcmsrot*HepGeom::Point3D<double> (0,0,1)).phi()*kRadToDeg<<std::endl ;
0326 */
0327 
0328       std::cout << "\n"
0329                 << names[i] << "   Theta_x = " << blcmsrot.thetaX() * kRadToDeg
0330                 << ",  Theta_y = " << blcmsrot.thetaY() * kRadToDeg << ",  Theta_z = " << blcmsrot.thetaZ() * kRadToDeg
0331                 << ",  eta = " << -log(tan(blcmsrot.thetaZ() / 2.)) << std::endl;
0332 
0333       std::cout << "\n"
0334                 << names[i] << "   Phi_x = " << blcmsrot.phiX() * kRadToDeg
0335                 << ",  Phi_y = " << blcmsrot.phiY() * kRadToDeg << ",  Phi_z = " << blcmsrot.phiZ() * kRadToDeg
0336                 << std::endl;
0337 
0338       std::cout << "\n\n BeamLine to CMS translation" << blcmstra << std::endl;
0339     }
0340     const double ntot(1.0 * nsur);
0341     sumx /= ntot;
0342     sumy /= ntot;
0343     sumz /= ntot;
0344     sumx2 = sqrt(sumx2 / ntot - sumx * sumx);
0345     sumy2 = sqrt(sumy2 / ntot - sumy * sumy);
0346     sumz2 = sqrt(sumz2 / ntot - sumz * sumz);
0347 
0348     std::cout << "\n Number of survey points = " << ntot;
0349     std::cout << "\n====== Mean x deviation = " << sumx << " +- " << sumx2;
0350     std::cout << "\n====== Mean y deviation = " << sumy << " +- " << sumy2;
0351     std::cout << "\n====== Mean z deviation = " << sumz << " +- " << sumz2 << std::endl;
0352 
0353     psumx /= ntot;
0354     psumy /= ntot;
0355     psumz /= ntot;
0356     psumx2 = sqrt(psumx2 / ntot - psumx * psumx);
0357     psumy2 = sqrt(psumy2 / ntot - psumy * psumy);
0358     psumz2 = sqrt(psumz2 / ntot - psumz * psumz);
0359 
0360     std::cout << "\n====== Mean x of pivot = " << psumx << " +- " << psumx2;
0361     std::cout << "\n====== Mean y of pivot = " << psumy << " +- " << psumy2;
0362     std::cout << "\n====== Mean z of pivot = " << psumz << " +- " << psumz2 << std::endl;
0363 
0364     fsumx /= ntot;
0365     fsumy /= ntot;
0366     fsumz /= ntot;
0367     fsumx2 = sqrt(fsumx2 / ntot - fsumx * fsumx);
0368     fsumy2 = sqrt(fsumy2 / ntot - fsumy * fsumy);
0369     fsumz2 = sqrt(fsumz2 / ntot - fsumz * fsumz);
0370 
0371     std::cout << "\n====== Mean x of floor = " << fsumx << " +- " << fsumx2;
0372     std::cout << "\n====== Mean y of floor = " << fsumy << " +- " << fsumy2;
0373     std::cout << "\n====== Mean z of floor = " << fsumz << " +- " << fsumz2 << std::endl;
0374 
0375     const HepGeom::Vector3D<double> hpl[] = {HepGeom::Plane3D<double>(sur[16], sur[0], sur[8]).normal(),
0376                                              HepGeom::Plane3D<double>(sur[17], sur[1], sur[9]).normal(),
0377                                              HepGeom::Plane3D<double>(sur[18], sur[2], sur[10]).normal(),
0378                                              HepGeom::Plane3D<double>(sur[19], sur[3], sur[11]).normal(),
0379                                              HepGeom::Plane3D<double>(sur[16], sur[12], sur[8]).normal(),
0380                                              HepGeom::Plane3D<double>(sur[17], sur[13], sur[9]).normal(),
0381                                              HepGeom::Plane3D<double>(sur[18], sur[14], sur[10]).normal(),
0382                                              HepGeom::Plane3D<double>(sur[19], sur[15], sur[11]).normal()};
0383 
0384     const HepGeom::Vector3D<double> hplavr(0.125 *
0385                                            (hpl[0] + hpl[1] + hpl[2] + hpl[3] + hpl[4] + hpl[5] + hpl[6] + hpl[7]));
0386 
0387     std::cout << "\n\n+++ HNormal for A =" << hpl[0].theta() * kRadToDeg << ", " << hpl[0].phi() * kRadToDeg
0388               << "\n+++ HNormal for B =" << hpl[1].theta() * kRadToDeg << ", " << hpl[1].phi() * kRadToDeg
0389               << "\n+++ HNormal for C =" << hpl[2].theta() * kRadToDeg << ", " << hpl[2].phi() * kRadToDeg
0390               << "\n+++ HNormal for D =" << hpl[3].theta() * kRadToDeg << ", " << hpl[3].phi() * kRadToDeg
0391               << "\n+++ HNormal for A =" << hpl[4].theta() * kRadToDeg << ", " << hpl[4].phi() * kRadToDeg
0392               << "\n+++ HNormal for B =" << hpl[5].theta() * kRadToDeg << ", " << hpl[5].phi() * kRadToDeg
0393               << "\n+++ HNormal for C =" << hpl[6].theta() * kRadToDeg << ", " << hpl[6].phi() * kRadToDeg
0394               << "\n+++ HNormal for D =" << hpl[7].theta() * kRadToDeg << ", " << hpl[7].phi() * kRadToDeg
0395               << "\n+++ HAverage      =" << hplavr.theta() * kRadToDeg << ", " << hplavr.phi() * kRadToDeg << std::endl;
0396 
0397     const HepGeom::Vector3D<double> vpl[] = {HepGeom::Plane3D<double>(sur[20], sur[28], sur[36]).normal(),
0398                                              HepGeom::Plane3D<double>(sur[21], sur[29], sur[37]).normal(),
0399                                              HepGeom::Plane3D<double>(sur[22], sur[30], sur[38]).normal(),
0400                                              HepGeom::Plane3D<double>(sur[23], sur[31], sur[39]).normal(),
0401 
0402                                              HepGeom::Plane3D<double>(sur[20], sur[24], sur[36]).normal(),
0403                                              HepGeom::Plane3D<double>(sur[21], sur[25], sur[37]).normal(),
0404                                              HepGeom::Plane3D<double>(sur[22], sur[26], sur[38]).normal(),
0405                                              HepGeom::Plane3D<double>(sur[23], sur[27], sur[39]).normal(),
0406 
0407                                              HepGeom::Plane3D<double>(sur[20], sur[32], sur[36]).normal(),
0408                                              HepGeom::Plane3D<double>(sur[21], sur[33], sur[37]).normal(),
0409                                              HepGeom::Plane3D<double>(sur[22], sur[34], sur[38]).normal(),
0410                                              HepGeom::Plane3D<double>(sur[23], sur[35], sur[39]).normal(),
0411 
0412                                              HepGeom::Plane3D<double>(sur[20], sur[32], sur[36]).normal(),
0413                                              HepGeom::Plane3D<double>(sur[21], sur[33], sur[37]).normal(),
0414                                              HepGeom::Plane3D<double>(sur[22], sur[34], sur[38]).normal(),
0415                                              HepGeom::Plane3D<double>(sur[23], sur[35], sur[39]).normal()};
0416 
0417     const HepGeom::Vector3D<double> vplavr(0.0625 * (vpl[0] + vpl[1] + vpl[2] + vpl[3] + vpl[4] + vpl[5] + vpl[6] +
0418                                                      vpl[7] + vpl[8] + vpl[9] + vpl[10] + vpl[11] + vpl[12] + vpl[13] +
0419                                                      vpl[14] + vpl[15]));
0420 
0421     std::cout << "\n\n+++ VNormal for A =" << vpl[0].theta() * kRadToDeg << ", " << vpl[0].phi() * kRadToDeg
0422               << "\n+++ VNormal for B =" << vpl[1].theta() * kRadToDeg << ", " << vpl[1].phi() * kRadToDeg
0423               << "\n+++ VNormal for C =" << vpl[2].theta() * kRadToDeg << ", " << vpl[2].phi() * kRadToDeg
0424               << "\n+++ VNormal for D =" << vpl[3].theta() * kRadToDeg << ", " << vpl[3].phi() * kRadToDeg
0425               << "\n+++ VNormal for A =" << vpl[4].theta() * kRadToDeg << ", " << vpl[4].phi() * kRadToDeg
0426               << "\n+++ VNormal for B =" << vpl[5].theta() * kRadToDeg << ", " << vpl[5].phi() * kRadToDeg
0427               << "\n+++ VNormal for C =" << vpl[6].theta() * kRadToDeg << ", " << vpl[6].phi() * kRadToDeg
0428               << "\n+++ VNormal for D =" << vpl[7].theta() * kRadToDeg << ", " << vpl[7].phi() * kRadToDeg
0429               << "\n+++ VNormal for A =" << vpl[8].theta() * kRadToDeg << ", " << vpl[8].phi() * kRadToDeg
0430               << "\n+++ VNormal for B =" << vpl[9].theta() * kRadToDeg << ", " << vpl[9].phi() * kRadToDeg
0431               << "\n+++ VNormal for C =" << vpl[10].theta() * kRadToDeg << ", " << vpl[10].phi() * kRadToDeg
0432               << "\n+++ VNormal for D =" << vpl[11].theta() * kRadToDeg << ", " << vpl[11].phi() * kRadToDeg
0433               << "\n+++ VNormal for A =" << vpl[12].theta() * kRadToDeg << ", " << vpl[12].phi() * kRadToDeg
0434               << "\n+++ VNormal for B =" << vpl[13].theta() * kRadToDeg << ", " << vpl[13].phi() * kRadToDeg
0435               << "\n+++ VNormal for C =" << vpl[14].theta() * kRadToDeg << ", " << vpl[14].phi() * kRadToDeg
0436               << "\n+++ VNormal for D =" << vpl[15].theta() * kRadToDeg << ", " << vpl[15].phi() * kRadToDeg
0437               << "\n+++ VAverage      =" << vplavr.theta() * kRadToDeg << ", " << vplavr.phi() * kRadToDeg << std::endl;
0438 
0439     unsigned int Hinput(0);
0440     unsigned int Vinput(0);
0441 
0442     std::cout << "\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n"
0443               << "\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n"
0444               << "\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n"
0445               << "\n+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++"
0446               << "\n+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++"
0447               << "\n+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++"
0448               << "\n+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++"
0449               << "\n+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++"
0450               << "\n+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++"
0451               << "\n+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++" << std::endl;
0452 
0453     bool keepGoing(true);
0454     while (keepGoing) {
0455       std::cout << "\n\n+++ Preparing to make two files for the table position you specify";
0456 
0457       std::cout << "\n\n**Enter desired H position (integer, 25000-300000) [0 to quit]: " << std::ends;
0458       std::cin >> Hinput;
0459 
0460       keepGoing = (0 != Hinput);
0461 
0462       if (keepGoing) {
0463         std::cout << "\n\n**Enter desired V position (integer, 5000-30000): " << std::ends;
0464         std::cin >> Vinput;
0465 
0466         std::cout << "\n\nH input value =" << Hinput << ", V input value =" << Vinput << std::endl;
0467 
0468         if (25000 > Hinput || 300000 < Hinput || 5000 > Vinput || 30000 < Vinput) {
0469           std::cout << "!!!!!!!! Invalid input(s), try again !!!!!!!" << std::endl;
0470         } else {
0471           // calculate transform parameters ================================================
0472 
0473           const HepGeom::Plane3D<double> plaV(sur[20], sur[28], sur[36]);
0474 
0475           const HepGeom::Point3D<double> projPiv(plaV.point(PIVOT));
0476 
0477           const double angleVPerUnit((plaV.point(sur[36]) - projPiv).angle(plaV.point(sur[20]) - projPiv) / 16000.);
0478 
0479           const HepGeom::Plane3D<double> plaH(sur[16], sur[0], sur[8]);
0480 
0481           const HepGeom::Point3D<double> projPivH(plaH.point(PIVOT));
0482 
0483           const double angleHPerUnit((plaH.point(sur[8]) - projPivH).angle(plaH.point(sur[16]) - projPivH) / 176458);
0484           /*
0485 
0486 
0487            const double angleV1 ( ( plaV.point( sur[24] ) - projPiv ).angle( plaV.point( sur[20] ) - projPiv ) ) ;
0488            const double angleV2 ( ( plaV.point( sur[28] ) - projPiv ).angle( plaV.point( sur[24] ) - projPiv ) ) ;
0489            const double angleV3 ( ( plaV.point( sur[32] ) - projPiv ).angle( plaV.point( sur[28] ) - projPiv ) ) ;
0490            const double angleV4 ( ( plaV.point( sur[36] ) - projPiv ).angle( plaV.point( sur[32] ) - projPiv ) ) ;
0491 
0492            const double angleH1 ( ( plaH.point( sur[12] ) - projPivH ).angle( plaH.point( sur[16] ) - projPivH ) ) ;
0493            const double angleH2 ( ( plaH.point( sur[ 0] ) - projPivH ).angle( plaH.point( sur[12] ) - projPivH ) ) ;
0494            const double angleH3 ( ( plaH.point( sur[ 4] ) - projPivH ).angle( plaH.point( sur[ 0] ) - projPivH ) ) ;
0495            const double angleH4 ( ( plaH.point( sur[ 8] ) - projPivH ).angle( plaH.point( sur[ 4] ) - projPivH ) ) ;
0496 
0497            std::cout<<"\n V angles are: "
0498             <<angleV1*kRadToDeg/4000.<<", "
0499             <<angleV2*kRadToDeg/4000.<<", "
0500             <<angleV3*kRadToDeg/4000.<<", "
0501             <<angleV4*kRadToDeg/4000.<<std::endl ;
0502 
0503            std::cout<<"\n H angles are: "
0504             <<angleH1*kRadToDeg/50000.<<", "
0505             <<angleH2*kRadToDeg/72402.<<", "
0506             <<angleH3*kRadToDeg/9995.<<", "
0507             <<angleH4*kRadToDeg/44601.<<std::endl ;
0508 
0509            std::cout<<"angleHPerUnit="<<kRadToDeg*angleHPerUnit*1.e3<<std::endl ;
0510 
0511            std::cout<<"angleVPerUnit="<<kRadToDeg*angleVPerUnit*1.e3<<std::endl ;
0512 */
0513 
0514           const double hAngle((1. * Hinput - 191000.) * angleHPerUnit);
0515 
0516           const double vAngle((1. * Vinput - 16000.) * angleVPerUnit);
0517 
0518           const HepGeom::Rotate3D vRotate(-vAngle, H4ToBeamLine * vpl[0]);
0519           const HepGeom::RotateY3D hRotate(hAngle);
0520 
0521           const HepGeom::Transform3D myTrf(V3Transform * vRotate * hRotate);
0522 
0523           const CLHEP::HepRotation myRot(myTrf.getRotation());
0524           const HepGeom::Vector3D<double> myTran(myTrf.getTranslation());
0525 
0526           const HepGeom::RotateZ3D R1(myRot.phiZ() - M_PI);
0527           const HepGeom::Point3D<double> xUnit(0, 1, 0);
0528           const HepGeom::Point3D<double> yUnit(-1, 0, 0);
0529           const HepGeom::Point3D<double> zUnit(0, 0, 1);
0530 
0531           const HepGeom::Transform3D RXRZ(HepGeom::Rotate3D(-myRot.thetaZ(), R1 * xUnit) * R1);
0532 
0533           const HepGeom::Point3D<double> zzUnit(RXRZ * zUnit);
0534           const HepGeom::Point3D<double> yyUnit(RXRZ * yUnit);
0535           const HepGeom::Point3D<double> xxUnit(RXRZ * xUnit);
0536           /*
0537            std::cout<<"new z axis = theta "
0538             << zzUnit.theta()*kRadToDeg 
0539             << ", "
0540             << zzUnit.phi()*kRadToDeg 
0541             <<"\nnew y axis = theta "
0542             << yyUnit.theta()*kRadToDeg 
0543             << ", "
0544             << yyUnit.phi()*kRadToDeg 
0545             <<"\nnew x axis = theta "
0546             << xxUnit.theta()*kRadToDeg 
0547             << ", "
0548             << xxUnit.phi()*kRadToDeg 
0549             << std::endl ;
0550 */
0551           const HepGeom::Transform3D RPSI(HepGeom::Transform3D(myRot, HepGeom::Point3D<double>(0, 0, 0)) *
0552                                           (RXRZ.inverse()));
0553 
0554           CLHEP::Hep3Vector axis;
0555           double psi;
0556           RPSI.getRotation().getAngleAxis(psi, axis);
0557 
0558           const double zdiff(axis.angle(zzUnit - HepGeom::Point3D<double>(0, 0, 0)));
0559 
0560           if (0.95 * M_PI < fabs(zdiff))
0561             psi = -psi;
0562           /*
0563            if( 0.005 < fabs( zdiff ) &&
0564            0.005 < fabs( fabs( zdiff ) - M_PI ) )
0565            {
0566           std::cout<<"**************** BAD 3RD AXIS, zdiff="
0567                <<zdiff*kRadToDeg<<" degrees"<<std::endl ;
0568            }
0569            std::cout<<"\n3rd rotation = "
0570             << axis.angle(zzUnit-HepPoint3D(0,0,0))*kRadToDeg 
0571             <<", "<<angle*kRadToDeg
0572             << std::endl ;
0573 
0574            std::cout<<"PSI="<<(RPSI*xUnit).phi()*kRadToDeg<<std::endl;
0575            std::cout<<"xUnitT="<<(RPSI*xUnit)<<std::endl;
0576            std::cout<<"RPSI="<<RPSI.getRotation()<<", "<<RPSI.getTranslation()<<std::endl;
0577 
0578 //HepRotateZ3D(-0.75*M_PI)*
0579            const HepGeom::Transform3D TEMP1 ( HepGeom::RotateX3D(-myRot.thetaZ())*HepGeom::RotateZ3D(M_PI/2.-myRot.phiZ()) ) ;
0580            const HepGeom::Transform3D TEMP ( HepGeom::Transform3D(myTrf.getRotation(),HepGeom::Vector3D<double> (0,0,0) )*TEMP1.inverse() ) ;
0581            std::cout<<"TEMP="<<TEMP.getRotation()<<", "<<TEMP.getTranslation()<<std::endl;
0582 
0583            std::cout<<"Trot="<<myTrf.getRotation()<<", "<<myTrf.getTranslation()<<std::endl;
0584 */
0585           const double TthetaX(myRot.thetaX());
0586           const double TthetaY(myRot.thetaY());
0587           const double TthetaZ(myRot.thetaZ());
0588           const double TphiX(myRot.phiX());
0589           const double TphiY(myRot.phiY());
0590           const double TphiZ(myRot.phiZ());
0591           const double Txtra(myTran.x());
0592           const double Tytra(myTran.y());
0593           const double Tztra(myTran.z());
0594           const double TetaZ(-log(tan(TthetaZ / 2.)));
0595 
0596           const double Tpsi(psi);
0597 
0598           // write to files ================================================================
0599 
0600           std::stringstream ssTpsir;
0601           ssTpsir << std::setw(7) << setiosflags(std::ios::fixed) << std::setprecision(5) << Tpsi;
0602           const std::string strTpsir(ssTpsir.str());
0603 
0604           std::stringstream ssTxtra;
0605           ssTxtra << std::setw(7) << setiosflags(std::ios::fixed) << std::setprecision(5) << Txtra << "*cm";
0606           const std::string strTxtra(ssTxtra.str());
0607           std::stringstream ssTytra;
0608           ssTytra << std::setw(7) << setiosflags(std::ios::fixed) << std::setprecision(5) << Tytra << "*cm";
0609           const std::string strTytra(ssTytra.str());
0610           std::stringstream ssTztra;
0611           ssTztra << std::setw(7) << setiosflags(std::ios::fixed) << std::setprecision(5) << Tztra << "*cm";
0612           const std::string strTztra(ssTztra.str());
0613 
0614           std::stringstream ssTetaZ;
0615           ssTetaZ << std::setw(7) << setiosflags(std::ios::fixed) << std::setprecision(5) << TetaZ;
0616           const std::string strTetaZ(ssTetaZ.str());
0617 
0618           std::stringstream ssTthetaX;
0619           ssTthetaX << std::setw(7) << setiosflags(std::ios::fixed) << std::setprecision(5) << TthetaX * kRadToDeg
0620                     << "*deg";
0621           const std::string strTthetaX(ssTthetaX.str());
0622           std::stringstream ssTphiX;
0623           ssTphiX << std::setw(7) << setiosflags(std::ios::fixed) << std::setprecision(5) << TphiX * kRadToDeg
0624                   << "*deg";
0625           const std::string strTphiX(ssTphiX.str());
0626           /*           std::stringstream ssTthetaXr ;
0627            ssTthetaXr<<std::setw(7) << setiosflags( std::ios::fixed ) << std::setprecision(5)
0628                << TthetaX ;
0629            const std::string strTthetaXr ( ssTthetaXr.str() ) ;
0630            std::stringstream ssTphiXr ;
0631            ssTphiXr<<std::setw(7) << setiosflags( std::ios::fixed ) << std::setprecision(5)
0632                << TphiX ;
0633            const std::string strTphiXr ( ssTphiXr.str() ) ;
0634 */
0635           std::stringstream ssTthetaY;
0636           ssTthetaY << std::setw(7) << setiosflags(std::ios::fixed) << std::setprecision(5) << TthetaY * kRadToDeg
0637                     << "*deg";
0638           const std::string strTthetaY(ssTthetaY.str());
0639           std::stringstream ssTphiY;
0640           ssTphiY << std::setw(7) << setiosflags(std::ios::fixed) << std::setprecision(5) << TphiY * kRadToDeg
0641                   << "*deg";
0642           const std::string strTphiY(ssTphiY.str());
0643           /*           std::stringstream ssTthetaYr ;
0644            ssTthetaYr<<std::setw(7) << setiosflags( std::ios::fixed ) << std::setprecision(5)
0645              << TthetaY ;
0646            const std::string strTthetaYr ( ssTthetaYr.str() ) ;
0647            std::stringstream ssTphiYr ;
0648            ssTphiYr<<std::setw(7) << setiosflags( std::ios::fixed ) << std::setprecision(5)
0649                << TphiY ;
0650            const std::string strTphiYr ( ssTphiYr.str() ) ;
0651 */
0652           std::stringstream ssTthetaZ;
0653           ssTthetaZ << std::setw(7) << setiosflags(std::ios::fixed) << std::setprecision(5) << TthetaZ * kRadToDeg
0654                     << "*deg";
0655           const std::string strTthetaZ(ssTthetaZ.str());
0656           std::stringstream ssTphiZ;
0657           ssTphiZ << std::setw(7) << setiosflags(std::ios::fixed) << std::setprecision(5) << TphiZ * kRadToDeg
0658                   << "*deg";
0659           const std::string strTphiZ(ssTphiZ.str());
0660           std::stringstream ssTphiZr;
0661           ssTphiZr << std::setw(7) << setiosflags(std::ios::fixed) << std::setprecision(5) << TphiZ;
0662           const std::string strTphiZr(ssTphiZr.str());
0663 
0664           std::stringstream ssTxpos;
0665           ssTxpos << std::setw(7) << setiosflags(std::ios::fixed) << std::setprecision(5) << psumy;
0666           const std::string strTxpos(ssTxpos.str());
0667 
0668           std::stringstream ssTypos;
0669           ssTypos << std::setw(7) << setiosflags(std::ios::fixed) << std::setprecision(5) << -psumx;
0670           const std::string strTypos(ssTypos.str());
0671 
0672           std::stringstream ssTzpos;
0673           ssTzpos << std::setw(11) << setiosflags(std::ios::fixed) << std::setprecision(2) << -26733.5 + psumz;
0674           const std::string strTzpos(ssTzpos.str());
0675           //-------------- first tbrot file ---------------------------
0676           std::stringstream ss;
0677           ss << "tbrot_H" << std::setw(6) << std::setfill('0') << Hinput << "_V" << std::setw(5) << std::setfill('0')
0678              << Vinput << ".xml" << std::ends;
0679           const std::string myname(ss.str());
0680 
0681           std::ofstream tbrotfile(myname.c_str());
0682 
0683           tbrotfile << "<?xml version=\"1.0\"?>"
0684                     << "<DDDefinition xmlns=\"http://www.cern.ch/cms/DDL\""
0685                     << " xmlns:xsi=\"http://www.w3.org/2001/XMLSchema-instance\""
0686                     << " xsi:schemaLocation=\"http://www.cern.ch/cms/DDL ../../"
0687                     << "DetectorDescription/Schema/DDLSchema.xsd\">"
0688                     << "\n\n<ConstantsSection  label=\"tbrot.xml\" eval=\"true\" >"
0689                     << "\n\n<Constant name=\"eta\"        value=\"" << strTetaZ
0690                     << "\"/>"
0691                     //          <<"\n\n<Constant name=\"theta\"      value=\"2*atan(exp(-[tbrot:eta]))\"/>"
0692                     << "\n\n<Constant name=\"theta\"      value=\"" << strTthetaZ << "\"/>"
0693                     << "\n\n<Constant name=\"phi\"        value=\"" << strTphiZ << "\"/>"
0694                     << "\n\n<Constant name=\"theta_X\"    value=\"" << strTthetaX << "\"/>"
0695                     << "\n\n<Constant name=\"phi_X\"      value=\"" << strTphiX << "\"/>"
0696                     << "\n\n<Constant name=\"theta_Y\"    value=\"" << strTthetaY << "\"/>"
0697                     << "\n\n<Constant name=\"phi_Y\"      value=\"" << strTphiY << "\"/>"
0698                     << "\n\n<Constant name=\"xtran\"        value=\"" << strTxtra << "\"/>"
0699                     << "\n<Constant name=\"ytran\"        value=\"" << strTytra << "\"/>"
0700                     << "\n<Constant name=\"ztran\"        value=\"" << strTztra << "\"/>"
0701                     << "\n\n</ConstantsSection>"
0702                     << "\n\n</DDDefinition>" << std::endl;
0703 
0704           tbrotfile.close();
0705 
0706           //--------------- now beam paramter file ----------------------
0707 
0708           std::stringstream ss2;
0709           ss2 << "ee_beam_H" << std::setw(6) << std::setfill('0') << Hinput << "_V" << std::setw(5) << std::setfill('0')
0710               << Vinput << "_cff.py" << std::ends;
0711           const std::string myname2(ss2.str());
0712 
0713           std::ofstream beamfile(myname2.c_str());
0714 
0715           beamfile << "import FWCore.ParameterSet.Config as cms"
0716                    << "\n\ncommon_beam_direction_parameters = cms.PSet("
0717                    << "\n    MinEta = cms.double(" << strTetaZ << "),"
0718                    << "\n    MaxEta = cms.double(" << strTetaZ << "),"
0719                    << "\n    MinPhi = cms.double(" << strTphiZr << "),"
0720                    << "\n    MaxPhi = cms.double(" << strTphiZr << "),"
0721                    << "\n    Psi    = cms.double(" << strTpsir << "),"
0722                    << "\n    BeamMeanX = cms.double(" << strTxpos << "),"
0723                    << "\n    BeamMeanY = cms.double(" << strTypos << "),"
0724                    << "\n    BeamPosition = cms.double(" << strTzpos << ")\n    )" << std::endl;
0725 
0726           beamfile.close();
0727         }
0728       }
0729     }
0730   }
0731 }
0732 
0733 template <class T>
0734 unsigned int CaloGeometryLoaderTest<T>::getDetIdForDDDNode(const DDFilteredView& fv) {
0735   // perform some consistency checks
0736   // get the parents and grandparents of this node
0737 
0738   const DDGeoHistory& parents(fv.geoHistory());
0739   const DDGeoHistory::size_type psize(parents.size());
0740 
0741   EcalBaseNumber baseNumber;
0742   baseNumber.setSize(psize);
0743 
0744   for (unsigned int i = 1; i <= psize; ++i) {
0745     baseNumber.addLevel(parents[psize - i].logicalPart().name().name(), parents[psize - i].copyno());
0746   }
0747 
0748   return baseNumber.getCopyNumber(0);
0749 
0750   //  return m_scheme.getUnitID( baseNumber );
0751 }
0752 
0753 #endif