DetectorDescription/Parser/src/DDDividedPolyhedra.cc - master - githubmirror/cmssw

DetectorDescription/Parser/src/DDDividedPolyhedra.cc
291 строка12 KB
Cms Build
Clang-Format
11 май 2019, 12:09
11 май 2019, 12:0902cea70
#include "DetectorDescription/Parser/src/DDDividedPolyhedra.h"
#include "DetectorDescription/Core/interface/DDRotationMatrix.h"
#include "DetectorDescription/Core/interface/DDLogicalPart.h"
#include "DetectorDescription/Core/interface/DDMaterial.h"
#include "DetectorDescription/Core/interface/DDName.h"
#include "DetectorDescription/Core/interface/DDSolid.h"
#include "DetectorDescription/Core/interface/DDTransform.h"
#include "DataFormats/Math/interface/GeantUnits.h"
#include "DetectorDescription/Parser/src/DDDividedGeometryObject.h"
#include "FWCore/MessageLogger/interface/MessageLogger.h"
#include "FWCore/Utilities/interface/Exception.h"

#include <cstddef>
#include <iostream>
#include <string>
#include <utility>
#include <vector>

class DDCompactView;

using namespace geant_units::operators;

DDDividedPolyhedraRho::DDDividedPolyhedraRho(const DDDivision& div, DDCompactView* cpv)
    : DDDividedGeometryObject(div, cpv) {
  checkParametersValidity();
  setType("DivisionPolyhedraRho");

  DDPolyhedra msol = (DDPolyhedra)(div_.parent().solid());

  if (divisionType_ == DivWIDTH) {
    compNDiv_ = calculateNDiv(msol.rMaxVec()[0] - msol.rMinVec()[0], div_.width(), div_.offset());
  } else if (divisionType_ == DivNDIV) {
    compWidth_ = calculateWidth(msol.rMaxVec()[0] - msol.rMinVec()[0], div_.nReplicas(), div_.offset());
  }
}

void DDDividedPolyhedraRho::checkParametersValidity(void) {
  DDDividedGeometryObject::checkParametersValidity();

  DDPolyhedra msol = (DDPolyhedra)(div_.parent().solid());

  if (divisionType_ == DivNDIVandWIDTH || divisionType_ == DivWIDTH) {
    std::cout << "WARNING - "
              << "DDDividedPolyhedraRho::checkParametersValidity()" << std::endl
              << "          Solid " << msol << std::endl
              << "          Division along R will be done with a width "
              << "different for each solid section." << std::endl
              << "          WIDTH will not be used !" << std::endl;
  }
  if (div_.offset() != 0.) {
    std::cout << "WARNING - "
              << "DDDividedPolyhedraRho::checkParametersValidity()" << std::endl
              << "          Solid " << msol << std::endl
              << "          Division along  R will be done with a width "
              << "different for each solid section." << std::endl
              << "          OFFSET will not be used !" << std::endl;
  }
}

double DDDividedPolyhedraRho::getMaxParameter(void) const {
  DDPolyhedra msol = (DDPolyhedra)(div_.parent().solid());
  return msol.rMaxVec()[0] - msol.rMinVec()[0];
}

DDTranslation DDDividedPolyhedraRho::makeDDTranslation(const int copyNo) const { return DDTranslation(); }

DDRotation DDDividedPolyhedraRho::makeDDRotation(const int copyNo) const { return DDRotation(); }

DDLogicalPart DDDividedPolyhedraRho::makeDDLogicalPart(const int copyNo) const {
  DDPolyhedra msol = (DDPolyhedra)(div_.parent().solid());
  DDMaterial usemat = div_.parent().material();

  std::vector<double> localrMaxVec = msol.rMaxVec();
  std::vector<double> localrMinVec = msol.rMinVec();
  std::vector<double> localzVec = msol.zVec();
  std::vector<double> newrMinVec;
  std::vector<double> newrMaxVec;
  int nZplanes = localzVec.size();

  double width = 0.;
  for (int ii = 0; ii < nZplanes; ++ii) {
    //     width = CalculateWidth( origparamMother->Rmax[ii]
    //                           - origparamMother->Rmin[ii], compNDiv_, foffset );
    //     origparam.Rmin[ii] = origparamMother->Rmin[ii]+foffset+width*copyNo;
    //     origparam.Rmax[ii] = origparamMother->Rmin[ii]+foffset+width*(copyNo+1);
    width = calculateWidth(localrMaxVec[ii] - localrMinVec[ii], compNDiv_, div_.offset());
    newrMinVec[ii] = localrMinVec[ii] + div_.offset() + width * copyNo;
    newrMaxVec[ii] = localrMaxVec[ii] + div_.offset() + width * (copyNo + 1);
  }

  //   phedra.SetOriginalParameters(&origparam); // copy values & transfer pointers
  //   phedra.Reset();                           // reset to new solid parameters

  DDName solname(div_.parent().ddname().name() + "_DIVCHILD" + std::to_string(copyNo), div_.parent().ddname().ns());

  DDSolid dsol = DDSolidFactory::polyhedra(
      solname, msol.sides(), msol.startPhi(), msol.deltaPhi(), localzVec, newrMinVec, newrMaxVec);
  DDLogicalPart ddlp = DDLogicalPart(solname, usemat, dsol);
  return ddlp;
}

DDDividedPolyhedraPhi::DDDividedPolyhedraPhi(const DDDivision& div, DDCompactView* cpv)
    : DDDividedGeometryObject(div, cpv) {
  checkParametersValidity();
  setType("DivisionPolyhedraPhi");

  DDPolyhedra msol = (DDPolyhedra)(div_.parent().solid());

  if (divisionType_ == DivWIDTH) {
    //If you divide a tube of 360 degrees the offset displaces the starting angle, but you still fill the 360 degrees
    if (msol.deltaPhi() == 360._deg) {
      compNDiv_ = calculateNDiv(msol.deltaPhi(), div_.width(), 0.);
    } else {
      compNDiv_ = calculateNDiv(msol.deltaPhi(), div_.width(), div_.offset());
    }
  } else if (divisionType_ == DivNDIV) {
    if (msol.deltaPhi() == 360._deg) {
      compWidth_ = calculateWidth(msol.deltaPhi(), div_.nReplicas(), 0.);
    } else {
      // original line looks wrong!
      compWidth_ = calculateWidth(msol.deltaPhi(), div_.nReplicas(), div_.offset());
    }
  }
}

double DDDividedPolyhedraPhi::getMaxParameter(void) const {
  DDPolyhedra msol = (DDPolyhedra)(div_.parent().solid());
  return msol.deltaPhi();  //msol->GetEndPhi() - msol->GetStartPhi();
}

void DDDividedPolyhedraPhi::checkParametersValidity(void) {
  DDDividedGeometryObject::checkParametersValidity();

  DDPolyhedra msol = (DDPolyhedra)(div_.parent().solid());

  if (divisionType_ == DivNDIVandWIDTH || divisionType_ == DivWIDTH) {
    std::cout << "WARNING - "
              << "DDDividedPolyhedraPhi::checkParametersValidity()" << std::endl
              << "          Solid " << msol << std::endl
              << "          Division along PHI will be done splitting "
              << "in the defined numSide." << std::endl
              << "          WIDTH will not be used !" << std::endl;
  }
  if (div_.offset() != 0.) {
    std::cout << "WARNING - "
              << "DDDividedPolyhedraPhi::checkParametersValidity()" << std::endl
              << "          Solid " << msol << std::endl
              << "          Division along PHI will be done splitting "
              << "in the defined numSide." << std::endl
              << "          OFFSET will not be used !" << std::endl;
  }

  if (msol.sides() != compNDiv_) {
    std::cout << "ERROR - "
              << "DDDividedPolyhedraPhi::checkParametersValidity()" << std::endl
              << "        Division along PHI will be done splitting in the defined" << std::endl
              << "        numSide, i.e, the number of division would be :"
              << "        " << msol.sides() << " instead of " << compNDiv_ << " !" << std::endl;
    std::string s = "DDDividedPolyhedraPhi::checkParametersValidity() Not supported configuration.";
    throw cms::Exception("DDException") << s;
  }
}

DDTranslation DDDividedPolyhedraPhi::makeDDTranslation(const int copyNo) const { return DDTranslation(); }

DDRotation DDDividedPolyhedraPhi::makeDDRotation(const int copyNo) const {
  DDRotation myddrot;  // sets to identity.
  double posi = (copyNo - 1) * compWidth_;
  DDName ddrotname(div_.parent().ddname().name() + "_DIVCHILD_ROT" + std::to_string(copyNo),
                   div_.parent().ddname().ns());
  myddrot = DDrot(ddrotname, changeRotMatrix(posi));

  return myddrot;
}

DDLogicalPart DDDividedPolyhedraPhi::makeDDLogicalPart(const int copyNo) const {
  DDPolyhedra msol = (DDPolyhedra)(div_.parent().solid());
  DDMaterial usemat = div_.parent().material();

  DDName solname(div_.parent().ddname().name() + "_DIVCHILD", div_.parent().ddname().ns());
  DDSolid dsol(solname);
  if (!dsol.isDefined().second) {
    dsol = DDSolidFactory::polyhedra(
        solname, msol.sides(), msol.startPhi() + div_.offset(), compWidth_, msol.zVec(), msol.rMinVec(), msol.rMaxVec());
  }
  DDLogicalPart ddlp(solname);
  if (!ddlp.isDefined().second)
    DDLogicalPart ddlp2 = DDLogicalPart(solname, usemat, dsol);
  return ddlp;
}

DDDividedPolyhedraZ::DDDividedPolyhedraZ(const DDDivision& div, DDCompactView* cpv)
    : DDDividedGeometryObject(div, cpv) {
  checkParametersValidity();
  setType("DivisionPolyhedraZ");

  DDPolyhedra msol = (DDPolyhedra)(div_.parent().solid());

  std::vector<double> zvec = msol.zVec();

  if (divisionType_ == DivWIDTH) {
    compNDiv_ = calculateNDiv(zvec[zvec.size() - 1] - zvec[0], div_.width(), div_.offset());
  } else if (divisionType_ == DivNDIV) {
    compWidth_ = calculateWidth(zvec[zvec.size() - 1] - zvec[0], div_.nReplicas(), div_.offset());
  }
}

double DDDividedPolyhedraZ::getMaxParameter(void) const {
  DDPolyhedra msol = (DDPolyhedra)(div_.parent().solid());

  std::vector<double> zvec = msol.zVec();
  return (zvec[zvec.size() - 1] - zvec[0]);
}

void DDDividedPolyhedraZ::checkParametersValidity(void) {
  DDDividedGeometryObject::checkParametersValidity();

  DDPolyhedra msol = (DDPolyhedra)(div_.parent().solid());

  if (divisionType_ == DivNDIVandWIDTH || divisionType_ == DivWIDTH) {
    std::cout << "WARNING - "
              << "DDDividedPolyhedraZ::checkParametersValidity()" << std::endl
              << "          Solid " << msol << std::endl
              << "          Division along Z will be done splitting "
              << "in the defined z_planes." << std::endl
              << "          WIDTH will not be used !" << std::endl;
  }

  if (div_.offset() != 0.) {
    std::cout << "WARNING - "
              << "DDDividedPolyhedraZ::checkParametersValidity()" << std::endl
              << "          Solid " << msol << std::endl
              << "          Division along Z will be done splitting "
              << "in the defined z_planes." << std::endl
              << "          OFFSET will not be used !" << std::endl;
  }

  std::vector<double> zvec = msol.zVec();

  if (zvec.size() - 1 != size_t(compNDiv_)) {
    std::cout << "ERROR - "
              << "DDDividedPolyhedraZ::checkParametersValidity()" << std::endl
              << "        Division along Z can only be done by splitting in the defined" << std::endl
              << "        z_planes, i.e, the number of division would be :"
              << "        " << zvec.size() - 1 << " instead of " << compNDiv_ << " !" << std::endl;
    std::string s = "DDDividedPolyhedraZ::checkParametersValidity()";
    s += "Illegal Construct. Not a supported configuration.";
    throw cms::Exception("DDException") << s;
  }
}

DDTranslation DDDividedPolyhedraZ::makeDDTranslation(const int copyNo) const {
  DDPolyhedra msol = (DDPolyhedra)(div_.parent().solid());
  std::vector<double> zvec = msol.zVec();

  //----- set translation: along Z axis
  double posi = (zvec[copyNo] + zvec[copyNo + 1]) / 2;

  DDTranslation tr(0, 0, posi);
  //----- calculate rotation matrix: unit
  return tr;
}

DDRotation DDDividedPolyhedraZ::makeDDRotation(const int copyNo) const { return DDRotation(); }

DDLogicalPart DDDividedPolyhedraZ::makeDDLogicalPart(const int copyNo) const {
  // only for mother number of planes = 2!!
  // mec: what?  why?  comment above and = 2 below straight from G4 impl.
  DDPolyhedra msol = (DDPolyhedra)(div_.parent().solid());
  DDMaterial usemat = div_.parent().material();

  std::vector<double> zvec = msol.zVec();
  std::vector<double> rminvec = msol.rMinVec();
  std::vector<double> rmaxvec = msol.rMaxVec();

  double posi = (zvec[copyNo] + zvec[copyNo + 1]) / 2.0;

  DDName solname(div_.parent().ddname().name() + "_DIVCHILD" + std::to_string(copyNo), div_.parent().ddname().ns());
  std::vector<double> newRmin, newRmax, newZ;
  newZ.emplace_back(zvec[copyNo] - posi);
  newZ.emplace_back(zvec[copyNo + 1] - posi);
  newRmin.emplace_back(rminvec[copyNo]);
  newRmin.emplace_back(rminvec[copyNo + 1]);
  newRmax.emplace_back(rmaxvec[copyNo]);
  newRmax.emplace_back(rmaxvec[copyNo + 1]);

  DDSolid dsol =
      DDSolidFactory::polyhedra(solname, msol.sides(), msol.startPhi(), msol.deltaPhi(), newZ, newRmin, newRmax);
  DDLogicalPart lp(solname, usemat, dsol);
  return lp;
}