C++ Peak::getQSampleFrame方法代码示例

/**
 * Returns selected information for a "peak" at QLabFrame.
 *
 * @param qFrame      An arbitrary position in Q-space.  This does not have to
 *be the
 *                    position of a peak.
 * @param labCoords  Set true if the position is in the lab coordinate system,
 *false if
 *                    it is in the sample coordinate system.
 * @return a vector whose elements contain different information about the
 *"peak" at that position.
 *         each element is a pair of description of information and the string
 *form for the corresponding
 *         value.
 */
int PeaksWorkspace::peakInfoNumber(Kernel::V3D qFrame, bool labCoords) const {
  std::vector<std::pair<std::string, std::string>> Result;
  std::ostringstream oss;
  oss << std::setw(12) << std::fixed << std::setprecision(3) << (qFrame.norm());
  std::pair<std::string, std::string> QMag("|Q|", oss.str());
  Result.push_back(QMag);

  oss.str("");
  oss.clear();
  oss << std::setw(12) << std::fixed << std::setprecision(3)
      << (2.0 * M_PI / qFrame.norm());

  std::pair<std::string, std::string> dspc("d-spacing", oss.str());
  oss.str("");
  oss.clear();
  Result.push_back(dspc);

  int seqNum = -1;
  double minDist = 10000000;

  for (int i = 0; i < getNumberPeaks(); i++) {
    Peak pk = getPeak(i);
    V3D Q = pk.getQLabFrame();
    if (!labCoords)
      Q = pk.getQSampleFrame();
    double D = qFrame.distance(Q);
    if (D < minDist) {
      minDist = D;
      seqNum = i + 1;
    }
  }
  return seqNum;
}

/**
 * Calculates the h,k, and l offsets from an integer for (some of )the peaks,
 *given the parameter values.
 *
 * @param out  For each peak there are 3 consecutive elements in this array. The
 *first is for the h offset from an
 *             integer, the second is the k offset and the 3rd is the l offset
 * @param xValues  xValues give the index in the PeaksWorkspace for the peak.
 *For each peak considered there are
 *              three consecutive entries all with the same index
 * @param nData The size of the xValues and out arrays
 */
void PeakHKLErrors::function1D(double *out, const double *xValues,
                               const size_t nData) const {
  PeaksWorkspace_sptr Peaks =
      AnalysisDataService::Instance().retrieveWS<PeaksWorkspace>(
          PeakWorkspaceName);

  boost::shared_ptr<Geometry::Instrument> instNew = getNewInstrument(Peaks);

  if (!Peaks)
    throw std::invalid_argument("Peaks not stored under the name " +
                                PeakWorkspaceName);

  std::map<int, Mantid::Kernel::Matrix<double>> RunNum2GonMatrixMap;
  getRun2MatMap(Peaks, OptRuns, RunNum2GonMatrixMap);
  const DblMatrix &UBx = Peaks->sample().getOrientedLattice().getUB();

  DblMatrix UBinv(UBx);
  UBinv.Invert();
  UBinv /= (2 * M_PI);

  double GonRotx = getParameter("GonRotx");
  double GonRoty = getParameter("GonRoty");
  double GonRotz = getParameter("GonRotz");
  Matrix<double> GonRot = RotationMatrixAboutRegAxis(GonRotx, 'x') *
                          RotationMatrixAboutRegAxis(GonRoty, 'y') *
                          RotationMatrixAboutRegAxis(GonRotz, 'z');

  double ChiSqTot = 0.0;
  for (size_t i = 0; i < nData; i += 3) {
    int peakNum = boost::math::iround(xValues[i]);
    IPeak &peak_old = Peaks->getPeak(peakNum);

    int runNum = peak_old.getRunNumber();
    std::string runNumStr = std::to_string(runNum);
    Peak peak = createNewPeak(peak_old, instNew, 0, peak_old.getL1());

    size_t N = OptRuns.find("/" + runNumStr + "/");
    if (N < OptRuns.size()) {
      peak.setGoniometerMatrix(GonRot * RunNum2GonMatrixMap[runNum]);

    } else {
      peak.setGoniometerMatrix(GonRot * peak.getGoniometerMatrix());
    }
    V3D sampOffsets(getParameter("SampleXOffset"),
                    getParameter("SampleYOffset"),
                    getParameter("SampleZOffset"));
    peak.setSamplePos(peak.getSamplePos() + sampOffsets);

    V3D hkl = UBinv * peak.getQSampleFrame();

    for (int k = 0; k < 3; k++) {
      double d1 = hkl[k] - floor(hkl[k]);
      if (d1 > .5)
        d1 = d1 - 1;
      if (d1 < -.5)
        d1 = d1 + 1;

      out[i + k] = d1;
      ChiSqTot += d1 * d1;
    }
  }

  g_log.debug() << "------------------------Function---------------------------"
                   "--------------------\n";
  for (size_t p = 0; p < nParams(); p++) {
    g_log.debug() << parameterName(p) << "(" << getParameter(p) << "),";
    if ((p + 1) % 6 == 0)
      g_log.debug() << '\n';
  }
  g_log.debug() << '\n';
  g_log.debug() << "Off constraints=";
  for (size_t p = 0; p < nParams(); p++) {
    IConstraint *constr = getConstraint(p);
    if (constr)
      if ((constr->check() > 0))
        g_log.debug() << "(" << parameterName(p) << "=" << constr->check()
                      << ");";
  }
  g_log.debug() << '\n';

  g_log.debug() << "    Chi**2 = " << ChiSqTot << "     nData = " << nData
                << '\n';
}

     void PeakHKLErrors::functionDeriv1D ( Jacobian* out, const double *xValues, const size_t nData )
     {
      PeaksWorkspace_sptr Peaks =
          AnalysisDataService::Instance().retrieveWS<PeaksWorkspace> ( PeakWorkspaceName );
      boost::shared_ptr<Geometry::Instrument> instNew = getNewInstrument( Peaks );

      const DblMatrix & UB = Peaks->sample().getOrientedLattice().getUB();
      DblMatrix UBinv( UB );
      UBinv.Invert();
      UBinv /= 2 * M_PI;

      double GonRotx = getParameter( "GonRotx" );
      double GonRoty = getParameter( "GonRoty" );
      double GonRotz = getParameter( "GonRotz" );
      Matrix<double> InvGonRotxMat = RotationMatrixAboutRegAxis( GonRotx,'x');
      Matrix<double> InvGonRotyMat = RotationMatrixAboutRegAxis( GonRoty,'y');
      Matrix<double> InvGonRotzMat = RotationMatrixAboutRegAxis( GonRotz,'z');
      Matrix<double> GonRot= InvGonRotxMat * InvGonRotyMat * InvGonRotzMat;

      InvGonRotxMat.Invert();
      InvGonRotyMat.Invert();
      InvGonRotzMat.Invert();

      std::map<int, Kernel::Matrix<double> > RunNums2GonMatrix;
      getRun2MatMap(  Peaks, OptRuns,RunNums2GonMatrix);

      g_log.debug()<<"----------------------------Derivative------------------------" << std::endl;

      V3D samplePosition = instNew->getSample()->getPos();
      IPeak& ppeak = Peaks->getPeak( 0 );
      double L0 = ppeak.getL1();
      double velocity = ( L0 + ppeak.getL2() ) / ppeak.getTOF();

       double K = 2 * M_PI / ppeak.getWavelength() / velocity;//2pi/lambda = K* velocity
       V3D beamDir = instNew->getBeamDirection();

       size_t paramNums[] = { parameterIndex( std::string( "SampleXOffset" ) ),
                              parameterIndex( std::string( "SampleYOffset" ) ),
                              parameterIndex( std::string( "SampleZOffset" ) ) };

      for ( size_t i = 0; i < nData; i += 3 )
      {
        int peakNum = (int)( .5 + xValues[i] );
        IPeak & peak_old = Peaks->getPeak( peakNum );
        Peak peak = SCDPanelErrors::createNewPeak( peak_old, instNew, 0, peak_old.getL1() );

        int runNum = peak_old.getRunNumber();
        std::string runNumStr = boost::lexical_cast<std::string>( runNum );

        for ( int kk = 0; kk <(int) nParams(); kk++ )
        {
          out->set( i, kk , 0.0 );
          out->set( i+1, kk , 0.0 );
          out->set( i+1, kk , 0.0 );
        }

        double chi, phi, omega;
        size_t chiParamNum, phiParamNum, omegaParamNum;

        size_t N = OptRuns.find( "/" +  runNumStr );
        if ( N < OptRuns.size() )
        {
          chi = getParameter( "chi" + ( runNumStr ) );
          phi = getParameter( "phi" + ( runNumStr ) );
          omega = getParameter( "omega" + ( runNumStr ) );

          peak.setGoniometerMatrix(GonRot*RunNums2GonMatrix[runNum] );

          chiParamNum = parameterIndex( "chi" + ( runNumStr ) );
          phiParamNum = parameterIndex( "phi" + ( runNumStr ) );
          omegaParamNum = parameterIndex( "omega" + ( runNumStr ) );
        }
        else
        {

          Geometry::Goniometer Gon( peak.getGoniometerMatrix() );
          std::vector<double> phichiOmega = Gon.getEulerAngles( "YZY" );
          chi = phichiOmega[1];
          phi = phichiOmega[2];
          omega = phichiOmega[0];
         // peak.setGoniometerMatrix( GonRot*Gon.getR());
          chiParamNum = phiParamNum =  omegaParamNum = nParams() + 10;
          peak.setGoniometerMatrix( GonRot*peak.getGoniometerMatrix());
        }
        //NOTE:Use getQLabFrame except for below.
        // For parameters the getGoniometerMatrix should remove GonRot, for derivs wrt GonRot*, wrt chi*,phi*,etc.


        V3D hkl = UBinv * peak.getQSampleFrame();

        //Deriv wrt chi phi and omega
        if ( phiParamNum < nParams() )
        {
          Matrix<double> chiMatrix= RotationMatrixAboutRegAxis( chi,'z');
          Matrix<double> phiMatrix= RotationMatrixAboutRegAxis( phi,'y');
          Matrix<double> omegaMatrix= RotationMatrixAboutRegAxis( omega,'y');

          Matrix<double> dchiMatrix= DerivRotationMatrixAboutRegAxis( chi,'z');
          Matrix<double> dphiMatrix= DerivRotationMatrixAboutRegAxis( phi,'y');
          Matrix<double> domegaMatrix= DerivRotationMatrixAboutRegAxis( omega,'y');
//.........这里部分代码省略.........