File indexing completed on 2023-03-17 11:14:30
0001
0002
0003
0004
0005
0006
0007 inline void volumeHandle::buildTrap(double x1,
0008 double x2,
0009 double x3,
0010 double x4,
0011 double y1,
0012 double y2,
0013 double theta,
0014 double phi,
0015 double halfZ,
0016 double alpha1,
0017 double alpha2) {
0018 LogTrace("MagGeoBuilder") << "Building trapezoid surfaces...: ";
0019 LogTrace("MagGeoBuilder") << "x1 " << x1 << newln << "x2 " << x2 << newln << "x3 " << x3 << newln << "x4 " << x4
0020 << newln << "y1 " << y1 << newln << "y2 " << y2 << newln << "theta " << theta << newln
0021 << "phi " << phi << newln << "halfZ " << halfZ << newln << "alpha1 " << alpha1 << newln
0022 << "alpha2 " << alpha2;
0023
0024
0025 if (debug) {
0026 const double epsilon = 1e-5;
0027 if (theta > epsilon || phi > epsilon || y1 - y2 > epsilon || x1 - x3 > epsilon || x2 - x4 > epsilon ||
0028 alpha1 - alpha2 > epsilon) {
0029 LogTrace("MagGeoBuilder") << "*** WARNING: unexpected trapezoid parameters.";
0030 }
0031 }
0032
0033
0034 GlobalVector planeXAxis = refPlane->toGlobal(LocalVector(1, 0, 0));
0035 GlobalVector planeYAxis = refPlane->toGlobal(LocalVector(0, 1, 0));
0036 GlobalVector planeZAxis = refPlane->toGlobal(LocalVector(0, 0, 1));
0037
0038
0039
0040
0041
0042
0043
0044
0045
0046
0047 GlobalVector Rvol(refPlane->position().x(), refPlane->position().y(), refPlane->position().z());
0048 theRN = std::abs(planeYAxis.dot(Rvol));
0049
0050 if (planeZAxis.z() < 0.) {
0051 x1 = -x1;
0052 x2 = -x2;
0053 halfZ = -halfZ;
0054 }
0055
0056 double halfX((x1 + x2) / 2.);
0057
0058 GlobalPoint pos_outer(refPlane->toGlobal(LocalPoint(0., y1, 0.)));
0059 GlobalPoint pos_inner(refPlane->toGlobal(LocalPoint(0., -y1, 0.)));
0060 GlobalPoint pos_zplus(refPlane->toGlobal(LocalPoint(0., 0., halfZ)));
0061 GlobalPoint pos_zminus(refPlane->toGlobal(LocalPoint(0., 0., -halfZ)));
0062 GlobalPoint pos_phiplus(refPlane->toGlobal(LocalPoint(-halfX, 0., 0.)));
0063 GlobalPoint pos_phiminus(refPlane->toGlobal(LocalPoint(halfX, 0., 0.)));
0064
0065 LogTrace("MagGeoBuilder") << "RN " << theRN << newln
0066 << "pos_outer " << pos_outer << " " << pos_outer.perp() << " " << pos_outer.phi()
0067 << newln << "pos_inner " << pos_inner << " " << pos_inner.perp() << " "
0068 << pos_inner.phi() << newln << "pos_zplus " << pos_zplus << " " << pos_zplus.perp()
0069 << " " << pos_zplus.phi() << newln << "pos_zminus " << pos_zminus << " "
0070 << pos_zminus.perp() << " " << pos_zminus.phi() << newln << "pos_phiplus " << pos_phiplus
0071 << " " << pos_phiplus.perp() << " " << pos_phiplus.phi() << newln << "pos_phiminus "
0072 << pos_phiminus << " " << pos_phiminus.perp() << " " << pos_phiminus.phi();
0073
0074
0075 if (debug) {
0076 if (pos_outer.perp() < pos_inner.perp()) {
0077 LogTrace("MagGeoBuilder") << "*** WARNING: pos_outer < pos_inner for trapezoid";
0078 }
0079 if (pos_zplus.z() < pos_zminus.z()) {
0080 LogTrace("MagGeoBuilder") << "*** WARNING: pos_zplus < pos_zminus for trapezoid";
0081 }
0082 if (Geom::Phi<float>(pos_phiplus.phi() - pos_phiminus.phi()) < 0.) {
0083 LogTrace("MagGeoBuilder") << "*** WARNING: pos_phiplus < pos_phiminus for trapezoid";
0084 }
0085 }
0086
0087
0088 GlobalVector y_phiplus = (refPlane->toGlobal(LocalVector((tan(alpha1) * y1 - (x2 - x1) / 2.), y1, 0.))).unit();
0089 GlobalVector y_phiminusV = (refPlane->toGlobal(LocalVector((tan(alpha1) * y1 + (x2 - x1) / 2.), y1, 0.)));
0090 GlobalVector y_phiminus = y_phiminusV.unit();
0091
0092 LogTrace("MagGeoBuilder") << "y_phiplus " << y_phiplus << newln << "y_phiminus " << y_phiminus;
0093
0094 Surface::RotationType rot_R(planeZAxis, planeXAxis);
0095 Surface::RotationType rot_Z(planeXAxis, planeYAxis);
0096 Surface::RotationType rot_phiplus(planeZAxis, y_phiplus);
0097 Surface::RotationType rot_phiminus(planeZAxis, y_phiminus);
0098
0099
0100 surfaces[outer] = new Plane(pos_outer, rot_R);
0101 surfaces[inner] = new Plane(pos_inner, rot_R);
0102 surfaces[zplus] = new Plane(pos_zplus, rot_Z);
0103 surfaces[zminus] = new Plane(pos_zminus, rot_Z);
0104 surfaces[phiplus] = new Plane(pos_phiplus, rot_phiplus);
0105 surfaces[phiminus] = new Plane(pos_phiminus, rot_phiminus);
0106
0107
0108 theRMin = std::abs(surfaces[inner]->toLocal(GlobalPoint(0, 0, 0)).z());
0109 theRMax = std::abs(surfaces[outer]->toLocal(GlobalPoint(0, 0, 0)).z());
0110
0111
0112
0113
0114
0115 if (name.substr(0, 3) == "V_Z") {
0116 thePhiMin = surfaces[phiminus]->position().phi();
0117 } else {
0118 thePhiMin = min((pos_phiminus + y_phiminusV).phi(), (pos_phiminus - y_phiminusV).phi());
0119 }
0120 LogTrace("MagGeoBuilder") << "rot_R " << surfaces[outer]->toGlobal(LocalVector(0., 0., 1.)) << newln << "rot_Z "
0121 << surfaces[zplus]->toGlobal(LocalVector(0., 0., 1.)) << newln << "rot_phi+ "
0122 << surfaces[phiplus]->toGlobal(LocalVector(0., 0., 1.)) << " phi "
0123 << surfaces[phiplus]->toGlobal(LocalVector(0., 0., 1.)).phi() << newln << "rot_phi- "
0124 << surfaces[phiminus]->toGlobal(LocalVector(0., 0., 1.)) << " phi "
0125 << surfaces[phiminus]->toGlobal(LocalVector(0., 0., 1.)).phi();
0126 if (debug) {
0127 if (std::abs(surfaces[phiminus]->position().phi() - thePhiMin) > 1e-5) {
0128 LogTrace("MagGeoBuilder") << " Trapez. phiMin = " << thePhiMin;
0129 }
0130 }
0131 }