C++ MatrixWorkspace_const_sptr::dataY方法代码示例

/** Do the initial copy of the data from the input to the output workspace for histogram workspaces.
 *  Takes out the bin width if necessary.
 *  @param inputWS  The input workspace
 *  @param outputWS The output workspace
 */
void ConvertUnits::fillOutputHist(const API::MatrixWorkspace_const_sptr inputWS, const API::MatrixWorkspace_sptr outputWS)
{
  const int size = static_cast<int>(inputWS->blocksize());

  // Loop over the histograms (detector spectra)
  Progress prog(this,0.0,0.2,m_numberOfSpectra);
  int64_t numberOfSpectra_i = static_cast<int64_t>(m_numberOfSpectra); // cast to make openmp happy
  PARALLEL_FOR2(inputWS,outputWS)
  for (int64_t i = 0; i < numberOfSpectra_i; ++i)
  {
    PARALLEL_START_INTERUPT_REGION
    // Take the bin width dependency out of the Y & E data
    if (m_distribution)
    {
      for (int j = 0; j < size; ++j)
      {
        const double width = std::abs( inputWS->dataX(i)[j+1] - inputWS->dataX(i)[j] );
        outputWS->dataY(i)[j] = inputWS->dataY(i)[j]*width;
        outputWS->dataE(i)[j] = inputWS->dataE(i)[j]*width;
      }
    }
    else
    {
      // Just copy over
      outputWS->dataY(i) = inputWS->dataY(i);
      outputWS->dataE(i) = inputWS->dataE(i);
    }
    // Copy over the X data
    outputWS->setX( i, inputWS->refX(i) );

    prog.report("Convert to " + m_outputUnit->unitID());
    PARALLEL_END_INTERUPT_REGION
  }
  PARALLEL_CHECK_INTERUPT_REGION
}

/** Execute the algorithm.
 */
void DampSq::exec()
{
    // TODO Auto-generated execute stub

    // 1. Generate new workspace
    API::MatrixWorkspace_const_sptr isqspace = getProperty("InputWorkspace");
    API::MatrixWorkspace_sptr osqspace = WorkspaceFactory::Instance().create(isqspace, 1, isqspace->size(), isqspace->size());

    int mode = getProperty("Mode");
    double qmax = getProperty("QMax");

    if (mode < 1 || mode > 4) {
        g_log.error("Damp mode can only be 1, 2, 3, or 4");
        return;
    }

    // 2. Get access to all
    const MantidVec& iQVec = isqspace->dataX(0);
    const MantidVec& iSVec = isqspace->dataY(0);
    const MantidVec& iEVec = isqspace->dataE(0);

    MantidVec& oQVec = osqspace->dataX(0);
    MantidVec& oSVec = osqspace->dataY(0);
    MantidVec& oEVec = osqspace->dataE(0);

    // 3. Calculation
    double dqmax = qmax - iQVec[0];

    double damp;
    for (unsigned int i = 0; i < iQVec.size(); i ++) {
        // a) calculate damp coefficient
        switch (mode) {
        case 1:
            damp = dampcoeff1(iQVec[i], qmax, dqmax);
            break;
        case 2:
            damp = dampcoeff2(iQVec[i], qmax, dqmax);;
            break;
        case 3:
            damp = dampcoeff3(iQVec[i], qmax, dqmax);;
            break;
        case 4:
            damp = dampcoeff4(iQVec[i], qmax, dqmax);;
            break;
        default:
            damp = 0;
            break;
        }
        // b) calculate new S(q)
        oQVec[i] = iQVec[i];
        oSVec[i] = 1 + damp*(iSVec[i]-1);
        oEVec[i] = damp*iEVec[i];
    }  // i

    // 4. Over
    setProperty("OutputWorkspace", osqspace);

    return;
}

/** Initialization method:
@param bkgWS    -- shared pointer to the workspace which contains background
@param sourceWS -- shared pointer to the workspace to remove background from
@param emode    -- energy conversion mode used during internal units conversion
(0 -- elastic, 1-direct, 2 indirect, as defined in Units conversion
@param pLog     -- pointer to the logger class which would report errors
@param nThreads -- number of threads to be used for background removal
@param inPlace  -- if the background removal occurs from the existing workspace
or target workspace has to be cloned.
*/
void BackgroundHelper::initialize(const API::MatrixWorkspace_const_sptr &bkgWS,
                                  const API::MatrixWorkspace_sptr &sourceWS,
                                  int emode, Kernel::Logger *pLog, int nThreads,
                                  bool inPlace) {
  m_bgWs = bkgWS;
  m_wkWS = sourceWS;
  m_Emode = emode;
  m_pgLog = pLog;
  m_inPlace = inPlace;

  std::string bgUnits = bkgWS->getAxis(0)->unit()->unitID();
  if (bgUnits != "TOF")
    throw std::invalid_argument(" Background Workspace: " + bkgWS->getName() +
                                " should be in the units of TOF");

  if (!(bkgWS->getNumberHistograms() == 1 ||
        sourceWS->getNumberHistograms() == bkgWS->getNumberHistograms()))
    throw std::invalid_argument(" Background Workspace: " + bkgWS->getName() +
                                " should have the same number of spectra as "
                                "source workspace or be a single histogram "
                                "workspace");

  auto WSUnit = sourceWS->getAxis(0)->unit();
  if (!WSUnit)
    throw std::invalid_argument(" Source Workspace: " + sourceWS->getName() +
                                " should have units");

  Geometry::IComponent_const_sptr source =
      sourceWS->getInstrument()->getSource();
  m_Sample = sourceWS->getInstrument()->getSample();
  if ((!source) || (!m_Sample))
    throw std::invalid_argument(
        "Instrument on Source workspace:" + sourceWS->getName() +
        "is not sufficiently defined: failed to get source and/or sample");
  m_L1 = source->getDistance(*m_Sample);

  // just in case.
  this->deleteUnitsConverters();
  // allocate the array of units converters to avoid units reallocation within a
  // loop
  m_WSUnit.assign(nThreads, NULL);
  for (int i = 0; i < nThreads; i++) {
    m_WSUnit[i] = WSUnit->clone();
  }

  m_singleValueBackground = false;
  if (bkgWS->getNumberHistograms() == 0)
    m_singleValueBackground = true;
  const MantidVec &dataX = bkgWS->dataX(0);
  const MantidVec &dataY = bkgWS->dataY(0);
  // const MantidVec& dataE = bkgWS->dataE(0);
  m_NBg = dataY[0];
  m_dtBg = dataX[1] - dataX[0];
  // m_ErrSq  = dataE[0]*dataE[0]; // needs further clarification

  m_Efix = this->getEi(sourceWS);
}

double BivariateNormal::initCommon() {

    double penalty = 0;
    bool ParamsOK = true;
    bool CommonsOK = true;
    if (!expVals)
        CommonsOK = false;

    API::MatrixWorkspace_const_sptr ws = getMatrixWorkspace();
    MantidVec D = ws->dataY(0);
    MantidVec X = ws->dataY(1);
    MantidVec Y = ws->dataY(2);

    if (NCells < 0) {
        NCells =
            (int)std::min<size_t>(D.size(), std::min<size_t>(X.size(), Y.size()));
        CommonsOK = false;
    }

    double Attrib[12] = {0.0};

    double MinX, MinY, MaxX, MaxY, MaxD, MinD;
    MinX = MaxX = X[0];
    MinY = MaxY = Y[0];
    MaxD = MinD = D[0];

    if (!CommonsOK) {

        for (int i = 0; i < NCells; i++) {
            Attrib[S_int] += D[i];
            Attrib[S_xint] += D[i] * X[i];
            Attrib[S_yint] += D[i] * Y[i];
            Attrib[S_x2int] += D[i] * X[i] * X[i];
            Attrib[S_y2int] += D[i] * Y[i] * Y[i];
            Attrib[S_xyint] += D[i] * X[i] * Y[i];

            Attrib[S_y] += Y[i];
            Attrib[S_x] += X[i];
            Attrib[S_x2] += X[i] * X[i];
            Attrib[S_y2] += Y[i] * Y[i];
            Attrib[S_xy] += X[i] * Y[i];
            Attrib[S_1] += 1.0;

            if (X[i] < MinX)
                MinX = X[i];
            if (X[i] > MaxX)
                MaxX = X[i];
            if (Y[i] < MinY)
                MinY = Y[i];
            if (Y[i] > MaxY)
                MaxY = Y[i];
            if (D[i] < MinD)
                MinD = D[i];
            if (D[i] > MaxD)
                MaxD = D[i];
        }

        mIx = Attrib[S_xint] / Attrib[S_int];
        mIy = Attrib[S_yint] / Attrib[S_int];
        mx = Attrib[S_x] / Attrib[S_1];
        my = Attrib[S_y] / Attrib[S_1];

        SIxx = Attrib[S_x2int] - (Attrib[S_xint] * Attrib[S_xint]) / Attrib[S_int];
        SIyy =
            Attrib[S_y2int] - (Attrib[S_yint]) * (Attrib[S_yint]) / Attrib[S_int];
        SIxy =
            Attrib[S_xyint] - (Attrib[S_xint]) * (Attrib[S_yint]) / Attrib[S_int];

        Sxx = Attrib[S_x2] - (Attrib[S_x]) * (Attrib[S_x]) / Attrib[S_1];
        Syy = Attrib[S_y2] - (Attrib[S_y]) * (Attrib[S_y]) / Attrib[S_1];
        Sxy = Attrib[S_xy] - (Attrib[S_x]) * (Attrib[S_y]) / Attrib[S_1];

        // CommonsOK = false;

        TotI = Attrib[S_int];
        TotN = Attrib[S_1];

        // CommonsOK = false;

        if (getConstraint(0) == NULL) {

            addConstraint((new BoundaryConstraint(this, "Background", 0,
                                                  Attrib[S_int] / Attrib[S_1])));
        }

        double maxIntensity = Attrib[S_int] + 3 * sqrt(Attrib[S_int]);

        if (maxIntensity < 100)
            maxIntensity = 100;

        if (getConstraint(1) == NULL) {
            addConstraint(new BoundaryConstraint(this, "Intensity", 0, maxIntensity));
        }

        double minMeany = MinY * .9 + .1 * MaxY;
        double maxMeany = MinY * .1 + .9 * MaxY;

        if (getConstraint(3) == NULL) {
            addConstraint(new BoundaryConstraint(this, "Mrow", minMeany, maxMeany));
        }
//.........这里部分代码省略.........