C++ V3D::spherical方法代码示例

    /**
     * Set the new detector position given the r,theta and phi.
     * @param det :: A pointer to the detector
     * @param l2 :: A single l2
     * @param theta :: A single theta
     * @param phi :: A single phi
     */
    void UpdateInstrumentFromFile::setDetectorPosition(const Geometry::IDetector_const_sptr & det, const float l2,
                                                       const float theta, const float phi)
    {
      if( m_ignoreMonitors && det->isMonitor() ) return;

      Geometry::ParameterMap & pmap = m_workspace->instrumentParameters();
      Kernel::V3D pos;
      if (!m_ignorePhi)
      {
        pos.spherical(l2, theta, phi);
      }
      else
      {
        double r,t,p;
        det->getPos().getSpherical(r,t,p);
        pos.spherical(l2, theta, p);
      }
      Geometry::ComponentHelper::moveComponent(*det, pmap, pos, Geometry::ComponentHelper::Absolute);
    }

/**
 * Set the new detector position given the r,theta and phi.
 * @param detectorInfo :: Reference to the DetectorInfo
 * @param index :: Index into detectorInfo
 * @param l2 :: A single l2
 * @param theta :: A single theta
 * @param phi :: A single phi
 */
void UpdateInstrumentFromFile::setDetectorPosition(
    Geometry::DetectorInfo &detectorInfo, const size_t index, const float l2,
    const float theta, const float phi) {
  if (m_ignoreMonitors && detectorInfo.isMonitor(index))
    return;

  Kernel::V3D pos;
  pos.spherical(l2, theta, phi);
  detectorInfo.setPosition(index, pos);
}

//.........这里部分代码省略.........
  std::vector<double> polar_width;
  file.getData(polar_width);
  file.closeData();

  // distance might not have been saved in all NXSPE files
  std::vector<double> distance;
  if (entries.count("distance")) {
    file.openData("distance");
    file.getData(distance);
    file.closeData();
  }

  file.closeGroup(); // data group
  file.closeGroup(); // Main entry
  file.close();

  // check if dimensions of the vectors are correct
  if ((error.size() != data.size()) || (azimuthal.size() != numSpectra) ||
      (azimuthal_width.size() != numSpectra) || (polar.size() != numSpectra) ||
      (polar_width.size() != numSpectra) ||
      ((energies.size() != numBins) && (energies.size() != numBins + 1))) {
    throw std::invalid_argument(
        "incompatible sizes of fields in the NXSPE file");
  }

  MatrixWorkspace_sptr outputWS = boost::dynamic_pointer_cast<MatrixWorkspace>(
      WorkspaceFactory::Instance().create("Workspace2D", numSpectra,
                                          energies.size(), numBins));
  // Need to get hold of the parameter map
  outputWS->getAxis(0)->unit() = UnitFactory::Instance().create("DeltaE");
  outputWS->setYUnit("SpectraNumber");

  // add logs
  outputWS->mutableRun().addLogData(
      new PropertyWithValue<double>("Ei", fixed_energy));
  outputWS->mutableRun().addLogData(new PropertyWithValue<double>("psi", psi));
  outputWS->mutableRun().addLogData(new PropertyWithValue<std::string>(
      "ki_over_kf_scaling", kikfscaling == 1 ? "true" : "false"));

  // Set Goniometer
  Geometry::Goniometer gm;
  gm.pushAxis("psi", 0, 1, 0, psi);
  outputWS->mutableRun().setGoniometer(gm, true);

  // generate instrument
  Geometry::Instrument_sptr instrument(new Geometry::Instrument("NXSPE"));
  outputWS->setInstrument(instrument);

  Geometry::ObjComponent *source = new Geometry::ObjComponent("source");
  source->setPos(0.0, 0.0, -10.0);
  instrument->add(source);
  instrument->markAsSource(source);
  Geometry::ObjComponent *sample = new Geometry::ObjComponent("sample");
  instrument->add(sample);
  instrument->markAsSamplePos(sample);

  Geometry::Object_const_sptr cuboid(
      createCuboid(0.1, 0.1, 0.1)); // FIXME: memory hog on rendering. Also,
                                    // make each detector separate size
  for (std::size_t i = 0; i < numSpectra; ++i) {
    double r = 1.0;
    if (!distance.empty()) {
      r = distance.at(i);
    }

    Kernel::V3D pos;
    pos.spherical(r, polar.at(i), azimuthal.at(i));

    Geometry::Detector *det =
        new Geometry::Detector("pixel", static_cast<int>(i + 1), sample);
    det->setPos(pos);
    det->setShape(cuboid);
    instrument->add(det);
    instrument->markAsDetector(det);
  }

  Geometry::ParameterMap &pmap = outputWS->instrumentParameters();
  std::vector<double>::iterator itdata = data.begin(), iterror = error.begin(),
                                itdataend, iterrorend;
  API::Progress prog = API::Progress(this, 0.0, 0.9, numSpectra);
  for (std::size_t i = 0; i < numSpectra; ++i) {
    itdataend = itdata + numBins;
    iterrorend = iterror + numBins;
    outputWS->dataX(i) = energies;
    if ((!boost::math::isfinite(*itdata)) || (*itdata <= -1e10)) // masked bin
    {
      outputWS->dataY(i) = std::vector<double>(numBins, 0);
      outputWS->dataE(i) = std::vector<double>(numBins, 0);
      pmap.addBool(outputWS->getDetector(i)->getComponentID(), "masked", true);
    } else {
      outputWS->dataY(i) = std::vector<double>(itdata, itdataend);
      outputWS->dataE(i) = std::vector<double>(iterror, iterrorend);
    }
    itdata = (itdataend);
    iterror = (iterrorend);
    prog.report();
  }

  setProperty("OutputWorkspace", outputWS);
}

/** Read the scaling information from a file (e.g. merlin_detector.sca) or from
 * the RAW file (.raw)
 *  @param scalingFile :: Name of scaling file .sca
 *  @param truepos :: V3D vector of actual positions as read from the file
 *  @return False if unable to open file, True otherwise
 */
bool SetScalingPSD::processScalingFile(const std::string &scalingFile,
                                       std::vector<Kernel::V3D> &truepos) {
  // Read the scaling information from a text file (.sca extension) or from a
  // raw file (.raw)
  // This is really corrected positions as (r,theta,phi) for each detector
  // Compare these with the instrument values to determine the change in
  // position and the scaling
  // which may be necessary for each pixel if in a tube.
  // movePos is used to updated positions
  std::map<int, Kernel::V3D> posMap;
  std::map<int, double> scaleMap;
  std::map<int, double>::iterator its;

  Instrument_const_sptr instrument = m_workspace->getInstrument();
  if (scalingFile.find(".sca") != std::string::npos ||
      scalingFile.find(".SCA") != std::string::npos) {
    // read a .sca text format file
    // format consists of a short header followed by one line per detector

    std::ifstream sFile(scalingFile.c_str());
    if (!sFile) {
      g_log.error() << "Unable to open scaling file " << scalingFile
                    << std::endl;
      return false;
    }
    std::string str;
    getline(sFile,
            str); // skip header line should be <filename> generated by <prog>
    int detectorCount;
    getline(sFile, str); // get detector count line
    std::istringstream istr(str);
    istr >> detectorCount;
    if (detectorCount < 1) {
      g_log.error("Bad detector count in scaling file");
      throw std::runtime_error("Bad detector count in scaling file");
    }
    truepos.reserve(detectorCount);
    getline(sFile, str); // skip title line
    int detIdLast = -10;
    Kernel::V3D truPosLast, detPosLast;

    Progress prog(this, 0.0, 0.5, detectorCount);
    // Now loop through lines, one for each detector/monitor. The latter are
    // ignored.

    while (getline(sFile, str)) {
      if (str.empty() || str[0] == '#')
        continue;
      std::istringstream istr(str);

      // read 6 values from the line to get the 3 (l2,theta,phi) of interest
      int detIndex, code;
      double l2, theta, phi, offset;
      istr >> detIndex >> offset >> l2 >> code >> theta >> phi;

      // sanity check on angles - l2 should be +ve but sample file has a few -ve
      // values
      // on monitors
      if (theta > 181.0 || theta < -1 || phi < -181 || phi > 181) {
        g_log.error("Position angle data out of range in .sca file");
        throw std::runtime_error(
            "Position angle data out of range in .sca file");
      }
      Kernel::V3D truPos;
      // use abs as correction file has -ve l2 for first few detectors
      truPos.spherical(fabs(l2), theta, phi);
      truepos.push_back(truPos);
      //
      Geometry::IDetector_const_sptr det;
      try {
        det = instrument->getDetector(detIndex);
      } catch (Kernel::Exception::NotFoundError &) {
        continue;
      }
      Kernel::V3D detPos = det->getPos();
      Kernel::V3D shift = truPos - detPos;

      // scaling applied to dets that are not monitors and have sequential IDs
      if (detIdLast == detIndex - 1 && !det->isMonitor()) {
        Kernel::V3D diffI = detPos - detPosLast;
        Kernel::V3D diffT = truPos - truPosLast;
        double scale = diffT.norm() / diffI.norm();
        Kernel::V3D scaleDir = diffT / diffT.norm();
        // Wish to store the scaling in a map, if we already have a scaling
        // for this detector (i.e. from the other side) we average the two
        // values. End of tube detectors only have one scaling estimate.
        scaleMap[detIndex] = scale;
        its = scaleMap.find(detIndex - 1);
        if (its == scaleMap.end())
          scaleMap[detIndex - 1] = scale;
        else
          its->second = 0.5 * (its->second + scale);
        // std::cout << detIndex << scale << scaleDir << std::endl;
      }
//.........这里部分代码省略.........