Skip to content

Commit

Permalink
refs #8589. Move transmission correction code into base class.
Browse files Browse the repository at this point in the history
More work will need to be done here because we yet need to separate out the making from the utilisation parts in order for this to be useful to both Correction and Reduction algorithms which use it.
  • Loading branch information
OwenArnold committed Dec 13, 2013
1 parent c7868e4 commit bcc2ae5
Show file tree
Hide file tree
Showing 4 changed files with 139 additions and 116 deletions.
Expand Up @@ -61,20 +61,6 @@ namespace Mantid

void exec();


/** Algorithm running methods **/

/// Perform a transmission correction on the input IvsLam workspace
API::MatrixWorkspace_sptr transmissonCorrection(API::MatrixWorkspace_sptr IvsLam,
const MinMax& wavelengthInterval, const MinMax& wavelengthMonitorBackgroundInterval,
const MinMax& wavelengthMonitorIntegrationInterval, const int& i0MonitorIndex,
API::MatrixWorkspace_sptr firstTransmissionRun,
OptionalMatrixWorkspace_sptr secondTransmissionRun, const OptionalDouble& stitchingStartQ,
const OptionalDouble& stitchingDeltaQ, const OptionalDouble& stitchingEndQ,
const OptionalDouble& stitchingStartOverlapQ, const OptionalDouble& stitchingEndOverlapQ,
const double& wavelengthStep
);

/// Get the surface sample component
Mantid::Geometry::IComponent_const_sptr getSurfaceSampleComponent(Mantid::Geometry::Instrument_const_sptr inst);

Expand Down
Expand Up @@ -84,6 +84,17 @@ namespace Mantid
/// Init common stitching inputs
void initStitchingInputs();


/// Perform a transmission correction on the input IvsLam workspace
API::MatrixWorkspace_sptr transmissonCorrection(API::MatrixWorkspace_sptr IvsLam,
const MinMax& wavelengthInterval, const MinMax& wavelengthMonitorBackgroundInterval,
const MinMax& wavelengthMonitorIntegrationInterval, const int& i0MonitorIndex,
API::MatrixWorkspace_sptr firstTransmissionRun,
OptionalMatrixWorkspace_sptr secondTransmissionRun, const OptionalDouble& stitchingStartQ,
const OptionalDouble& stitchingDeltaQ, const OptionalDouble& stitchingEndQ,
const OptionalDouble& stitchingStartOverlapQ, const OptionalDouble& stitchingEndOverlapQ,
const double& wavelengthStep );

private:

/// Validate the transmission correction property inputs
Expand Down
102 changes: 0 additions & 102 deletions Code/Mantid/Framework/Algorithms/src/ReflectometryReductionOne.cpp
Expand Up @@ -227,108 +227,6 @@ namespace Mantid
}



/**
* Perform Transmission Corrections.
* @param IvsLam : Run workspace which is to be normalized by the results of the transmission corrections.
* @param wavelengthInterval : Wavelength interval for the run workspace.
* @param wavelengthMonitorBackgroundInterval : Wavelength interval for the monitor background
* @param wavelengthMonitorIntegrationInterval : Wavelength interval for the monitor integration
* @param i0MonitorIndex : Monitor index for the I0 monitor
* @param firstTransmissionRun : The first transmission run
* @param secondTransmissionRun : The second transmission run (optional)
* @param stitchingStartQ : Stitching start Q (optional but dependent on secondTransmissionRun)
* @param stitchingDeltaQ : Stitching delta Q (optional but dependent on secondTransmissionRun)
* @param stitchingEndQ : Stitching end Q (optional but dependent on secondTransmissionRun)
* @param stitchingStartOverlapQ : Stitching start Q overlap (optional but dependent on secondTransmissionRun)
* @param stitchingEndOverlapQ : Stitching end Q overlap (optional but dependent on secondTransmissionRun)
* @param wavelengthStep : Step in angstroms for rebinning for workspaces converted into wavelength.
* @return Normalized run workspace by the transmission workspace, which have themselves been converted to Lam, normalized by monitors and possibly stitched together.
*/
MatrixWorkspace_sptr ReflectometryReductionOne::transmissonCorrection(MatrixWorkspace_sptr IvsLam,
const MinMax& wavelengthInterval,
const MinMax& wavelengthMonitorBackgroundInterval,
const MinMax& wavelengthMonitorIntegrationInterval,
const int& i0MonitorIndex,
MatrixWorkspace_sptr firstTransmissionRun,
OptionalMatrixWorkspace_sptr secondTransmissionRun,
const OptionalDouble& stitchingStartQ,
const OptionalDouble& stitchingDeltaQ,
const OptionalDouble& stitchingEndQ,
const OptionalDouble& stitchingStartOverlapQ,
const OptionalDouble& stitchingEndOverlapQ,
const double& wavelengthStep
)
{
g_log.debug("Extracting first transmission run workspace indexes from spectra");
const WorkspaceIndexList detectorIndexes = createWorkspaceIndexListFromDetectorWorkspace(IvsLam, firstTransmissionRun);
auto trans1InLam = toLam(firstTransmissionRun, detectorIndexes, i0MonitorIndex, wavelengthInterval, wavelengthMonitorBackgroundInterval, wavelengthStep);
MatrixWorkspace_sptr trans1Detector = trans1InLam.get<0>();
MatrixWorkspace_sptr trans1Monitor = trans1InLam.get<1>();

// Monitor integration ... can this happen inside the toLam routine?
auto integrationAlg = this->createChildAlgorithm("Integration");
integrationAlg->initialize();
integrationAlg->setProperty("InputWorkspace", trans1Monitor);
integrationAlg->setProperty("RangeLower", wavelengthMonitorIntegrationInterval.get<0>());
integrationAlg->setProperty("RangeUpper", wavelengthMonitorIntegrationInterval.get<1>());
integrationAlg->execute();
trans1Monitor = integrationAlg->getProperty("OutputWorkspace");

MatrixWorkspace_sptr denominator = trans1Detector / trans1Monitor;

if (secondTransmissionRun.is_initialized())
{
auto transRun2 = secondTransmissionRun.get();
g_log.debug("Extracting second transmission run workspace indexes from spectra");
const WorkspaceIndexList detectorIndexes = createWorkspaceIndexListFromDetectorWorkspace(IvsLam, transRun2);

auto trans2InLam = toLam(transRun2, detectorIndexes, i0MonitorIndex, wavelengthInterval,
wavelengthMonitorBackgroundInterval, wavelengthStep);

// Unpack the conversion results.
MatrixWorkspace_sptr trans2Detector = trans2InLam.get<0>();
MatrixWorkspace_sptr trans2Monitor = trans2InLam.get<1>();

// Monitor integration ... can this happen inside the toLam routine?
auto integrationAlg = this->createChildAlgorithm("Integration");
integrationAlg->initialize();
integrationAlg->setProperty("InputWorkspace", trans2Monitor);
integrationAlg->setProperty("RangeLower", wavelengthMonitorIntegrationInterval.get<0>());
integrationAlg->setProperty("RangeUpper", wavelengthMonitorIntegrationInterval.get<1>());
integrationAlg->execute();
trans2Monitor = integrationAlg->getProperty("OutputWorkspace");

MatrixWorkspace_sptr normalizedTrans2 = trans2Detector / trans2Monitor;

// Stitch the results.
auto stitch1DAlg = this->createChildAlgorithm("Stitch1D");
stitch1DAlg->initialize();
AnalysisDataService::Instance().addOrReplace("denominator", denominator);
AnalysisDataService::Instance().addOrReplace("normalizedTrans2", normalizedTrans2);
stitch1DAlg->setProperty("LHSWorkspace", denominator);
stitch1DAlg->setProperty("RHSWorkspace", normalizedTrans2);
stitch1DAlg->setProperty("StartOverlap", stitchingStartOverlapQ.get() );
stitch1DAlg->setProperty("EndOverlap", stitchingEndOverlapQ.get() );
const std::vector<double> params = boost::assign::list_of(stitchingStartQ.get())(stitchingDeltaQ.get())(stitchingEndQ.get()).convert_to_container<std::vector< double > >();
stitch1DAlg->setProperty("Params", params);
stitch1DAlg->execute();
denominator = stitch1DAlg->getProperty("OutputWorkspace");
AnalysisDataService::Instance().remove("denominator");
AnalysisDataService::Instance().remove("normalizedTrans2");
}

auto rebinToWorkspaceAlg = this->createChildAlgorithm("RebinToWorkspace");
rebinToWorkspaceAlg->initialize();
rebinToWorkspaceAlg->setProperty("WorkspaceToMatch", IvsLam);
rebinToWorkspaceAlg->setProperty("WorkspaceToRebin", denominator);
rebinToWorkspaceAlg->execute();
denominator = rebinToWorkspaceAlg->getProperty("OutputWorkspace");

MatrixWorkspace_sptr normalizedIvsLam = IvsLam / denominator;
return normalizedIvsLam;
}

/**
* Correct the position of the detectors based on the input theta value.
* @param toCorrect : Workspace to correct detector posisitions on.
Expand Down
128 changes: 128 additions & 0 deletions Code/Mantid/Framework/Algorithms/src/ReflectometryWorkflowBase.cpp
Expand Up @@ -16,6 +16,29 @@ namespace Mantid
{
namespace
{
/**
* Helper non-member function for translating all the workspace indexes in an origin workspace into workspace indexes
* of a host end-point workspace. This is done using spectrum numbers as the intermediate.
*
* This function will throw a runtime error if the specId are not found to exist on the host end-point workspace.
*
* @param originWS : Origin workspace, which provides the original workspace index to spectrum id mapping.
* @param hostWS : Workspace onto which the resulting workspace indexes will be hosted
* @return Remapped wokspace indexes applicable for the host workspace.
*/
ReflectometryWorkflowBase::WorkspaceIndexList createWorkspaceIndexListFromDetectorWorkspace(
MatrixWorkspace_const_sptr originWS, MatrixWorkspace_const_sptr hostWS)
{
auto spectrumMap = originWS->getSpectrumToWorkspaceIndexMap();
ReflectometryWorkflowBase::WorkspaceIndexList translatedIndexList;
for (auto it = spectrumMap.begin(); it != spectrumMap.end(); ++it)
{
specid_t specId = (*it).first;
translatedIndexList.push_back(static_cast<int>(hostWS->getIndexFromSpectrumNumber(specId))); // Could be slow to do it like this.
}
return translatedIndexList;
}

/**
* Helper method used with the stl to determine whether values are negative
* @param value : Value to check
Expand Down Expand Up @@ -444,5 +467,110 @@ namespace Mantid
return DetectorMonitorWorkspacePair(detectorWS, monitorWS);
}




/**
* Perform Transmission Corrections.
* @param IvsLam : Run workspace which is to be normalized by the results of the transmission corrections.
* @param wavelengthInterval : Wavelength interval for the run workspace.
* @param wavelengthMonitorBackgroundInterval : Wavelength interval for the monitor background
* @param wavelengthMonitorIntegrationInterval : Wavelength interval for the monitor integration
* @param i0MonitorIndex : Monitor index for the I0 monitor
* @param firstTransmissionRun : The first transmission run
* @param secondTransmissionRun : The second transmission run (optional)
* @param stitchingStartQ : Stitching start Q (optional but dependent on secondTransmissionRun)
* @param stitchingDeltaQ : Stitching delta Q (optional but dependent on secondTransmissionRun)
* @param stitchingEndQ : Stitching end Q (optional but dependent on secondTransmissionRun)
* @param stitchingStartOverlapQ : Stitching start Q overlap (optional but dependent on secondTransmissionRun)
* @param stitchingEndOverlapQ : Stitching end Q overlap (optional but dependent on secondTransmissionRun)
* @param wavelengthStep : Step in angstroms for rebinning for workspaces converted into wavelength.
* @return Normalized run workspace by the transmission workspace, which have themselves been converted to Lam, normalized by monitors and possibly stitched together.
*/
MatrixWorkspace_sptr ReflectometryWorkflowBase::transmissonCorrection(MatrixWorkspace_sptr IvsLam,
const MinMax& wavelengthInterval,
const MinMax& wavelengthMonitorBackgroundInterval,
const MinMax& wavelengthMonitorIntegrationInterval,
const int& i0MonitorIndex,
MatrixWorkspace_sptr firstTransmissionRun,
OptionalMatrixWorkspace_sptr secondTransmissionRun,
const OptionalDouble& stitchingStartQ,
const OptionalDouble& stitchingDeltaQ,
const OptionalDouble& stitchingEndQ,
const OptionalDouble& stitchingStartOverlapQ,
const OptionalDouble& stitchingEndOverlapQ,
const double& wavelengthStep
)
{
g_log.debug("Extracting first transmission run workspace indexes from spectra");
const WorkspaceIndexList detectorIndexes = createWorkspaceIndexListFromDetectorWorkspace(IvsLam, firstTransmissionRun);
auto trans1InLam = toLam(firstTransmissionRun, detectorIndexes, i0MonitorIndex, wavelengthInterval, wavelengthMonitorBackgroundInterval, wavelengthStep);
MatrixWorkspace_sptr trans1Detector = trans1InLam.get<0>();
MatrixWorkspace_sptr trans1Monitor = trans1InLam.get<1>();

// Monitor integration ... can this happen inside the toLam routine?
auto integrationAlg = this->createChildAlgorithm("Integration");
integrationAlg->initialize();
integrationAlg->setProperty("InputWorkspace", trans1Monitor);
integrationAlg->setProperty("RangeLower", wavelengthMonitorIntegrationInterval.get<0>());
integrationAlg->setProperty("RangeUpper", wavelengthMonitorIntegrationInterval.get<1>());
integrationAlg->execute();
trans1Monitor = integrationAlg->getProperty("OutputWorkspace");

MatrixWorkspace_sptr denominator = trans1Detector / trans1Monitor;

if (secondTransmissionRun.is_initialized())
{
auto transRun2 = secondTransmissionRun.get();
g_log.debug("Extracting second transmission run workspace indexes from spectra");
const WorkspaceIndexList detectorIndexes = createWorkspaceIndexListFromDetectorWorkspace(IvsLam, transRun2);

auto trans2InLam = toLam(transRun2, detectorIndexes, i0MonitorIndex, wavelengthInterval,
wavelengthMonitorBackgroundInterval, wavelengthStep);

// Unpack the conversion results.
MatrixWorkspace_sptr trans2Detector = trans2InLam.get<0>();
MatrixWorkspace_sptr trans2Monitor = trans2InLam.get<1>();

// Monitor integration ... can this happen inside the toLam routine?
auto integrationAlg = this->createChildAlgorithm("Integration");
integrationAlg->initialize();
integrationAlg->setProperty("InputWorkspace", trans2Monitor);
integrationAlg->setProperty("RangeLower", wavelengthMonitorIntegrationInterval.get<0>());
integrationAlg->setProperty("RangeUpper", wavelengthMonitorIntegrationInterval.get<1>());
integrationAlg->execute();
trans2Monitor = integrationAlg->getProperty("OutputWorkspace");

MatrixWorkspace_sptr normalizedTrans2 = trans2Detector / trans2Monitor;

// Stitch the results.
auto stitch1DAlg = this->createChildAlgorithm("Stitch1D");
stitch1DAlg->initialize();
AnalysisDataService::Instance().addOrReplace("denominator", denominator);
AnalysisDataService::Instance().addOrReplace("normalizedTrans2", normalizedTrans2);
stitch1DAlg->setProperty("LHSWorkspace", denominator);
stitch1DAlg->setProperty("RHSWorkspace", normalizedTrans2);
stitch1DAlg->setProperty("StartOverlap", stitchingStartOverlapQ.get() );
stitch1DAlg->setProperty("EndOverlap", stitchingEndOverlapQ.get() );
const std::vector<double> params = boost::assign::list_of(stitchingStartQ.get())(stitchingDeltaQ.get())(stitchingEndQ.get()).convert_to_container<std::vector< double > >();
stitch1DAlg->setProperty("Params", params);
stitch1DAlg->execute();
denominator = stitch1DAlg->getProperty("OutputWorkspace");
AnalysisDataService::Instance().remove("denominator");
AnalysisDataService::Instance().remove("normalizedTrans2");
}

auto rebinToWorkspaceAlg = this->createChildAlgorithm("RebinToWorkspace");
rebinToWorkspaceAlg->initialize();
rebinToWorkspaceAlg->setProperty("WorkspaceToMatch", IvsLam);
rebinToWorkspaceAlg->setProperty("WorkspaceToRebin", denominator);
rebinToWorkspaceAlg->execute();
denominator = rebinToWorkspaceAlg->getProperty("OutputWorkspace");

// Do normalization.
MatrixWorkspace_sptr normalizedIvsLam = IvsLam / denominator;
return normalizedIvsLam;
}

} // namespace Algorithms
} // namespace Mantid

0 comments on commit bcc2ae5

Please sign in to comment.