Skip to content

Commit

Permalink
Revert "Implementing generating masking workspace. Refs #4370."
Browse files Browse the repository at this point in the history 
This broke the direct inelastic system tests and will need to be
re-worked before being committed again.
This reverts commit c2e67af.
  • Loading branch information
@RussellTaylor
RussellTaylor committed Dec 30, 2011
1 parent b902933 commit 1fb1684
Show file tree
Hide file tree
Showing 5 changed files with 87 additions and 269 deletions.
1 change: 0 additions & 1 deletion Code/Mantid/Framework/Algorithms/inc/MantidAlgorithms/FindDetectorsOutsideLimits.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -5,7 +5,6 @@
// Includes
//----------------------------------------------------------------------
#include "MantidAlgorithms/DetectorDiagnostic.h"
#include "MantidDataObjects/SpecialWorkspace2D.h"

namespace Mantid
{
Expand Down
7 changes: 3 additions & 4 deletions Code/Mantid/Framework/Algorithms/inc/MantidAlgorithms/MedianDetectorTest.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -5,7 +5,6 @@
// Includes
//----------------------------------------------------------------------
#include "MantidAlgorithms/DetectorDiagnostic.h"
#include "MantidDataObjects/SpecialWorkspace2D.h"

namespace Mantid
{
Expand Down Expand Up @@ -84,9 +83,9 @@ namespace Mantid
API::MatrixWorkspace_sptr getSolidAngles(int firstSpec, int lastSpec);
/// Do the tests and mask those that fail
int doDetectorTests(const API::MatrixWorkspace_sptr countWorkspace,
DataObjects::SpecialWorkspace2D_sptr maskWS,
const double average,
const std::set<int> & badIndices);
API::MatrixWorkspace_sptr maskWS,
const double average,
const std::set<int> & badIndices);

/// Input workspace
API::MatrixWorkspace_sptr m_inputWS;
Expand Down
78 changes: 19 additions & 59 deletions Code/Mantid/Framework/Algorithms/src/FindDetectorsOutsideLimits.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -44,14 +44,9 @@ namespace Mantid
new WorkspaceProperty<MatrixWorkspace>("InputWorkspace","",Direction::Input),
"Name of the input workspace2D" );
declareProperty(
new WorkspaceProperty<DataObjects::SpecialWorkspace2D>("OutputWorkspace","",Direction::Output),
new WorkspaceProperty<MatrixWorkspace>("OutputWorkspace","",Direction::Output),
"Each histogram from the input workspace maps to a histogram in this\n"
"SpecialWorkspace2D (Masking workdpace) with one value that indicates if there was a dead detector" );
/*
declareProperty(
new WorkspaceProperty<MatrixWorkspace>("OutputWorkspace2", "", Direction::Output),
"Same as OutputWorksapce, but in Workspace2D");
*/
"workspace with one value that indicates if there was a dead detector" );
declareProperty("HighThreshold", 0.0,
"Spectra whose total number of counts are equal to or above this value\n"
"will be marked bad (default 0)" );
Expand All @@ -76,7 +71,6 @@ namespace Mantid
*/
void FindDetectorsOutsideLimits::exec()
{
// 1. Get properties (inputs)
double lowThreshold = getProperty("LowThreshold");
double highThreshold = getProperty("HighThreshold");
if ( highThreshold <= lowThreshold )
Expand All @@ -91,84 +85,53 @@ namespace Mantid
const double rangeLower = getProperty("RangeLower");
const double rangeUpper = getProperty("RangeUpper");

// 2. Create output workspace and set it up & check input workspace size
Geometry::Instrument_const_sptr instrument = inputWS->getInstrument();
DataObjects::SpecialWorkspace2D_sptr maskWS(new DataObjects::SpecialWorkspace2D(instrument));
this->setProperty("OutputWorkspace", maskWS);

if (maskWS->getNumberHistograms() > inputWS->getNumberHistograms()){
g_log.warning() << "Input Workspace " << inputWS->getName() << " has different number of spectra (" <<
inputWS->getNumberHistograms() << ") to new generated masking workspace (" <<
maskWS->getNumberHistograms() << "). " << " Input workspace might be focused" << std::endl;
}
else if (maskWS->getNumberHistograms() < inputWS->getNumberHistograms()){
g_log.error() << "Input Workspace " << inputWS->getName() << " has different number of spectra (" <<
inputWS->getNumberHistograms() << ") to new generated masking workspace (" <<
maskWS->getNumberHistograms() << "). " << " Input workspace might be focused" << std::endl;
throw std::runtime_error("Error! This cannot be right. Instrument must be wrong.");
}

// 3. Get the integrated input workspace; converting to a Workspace2D
// Get the integrated input workspace; converting to a Workspace2D
MatrixWorkspace_sptr integratedWorkspace =
integrateSpectra(inputWS, 0, EMPTY_INT(), rangeLower, rangeUpper, true);

EventWorkspace_sptr ew = boost::dynamic_pointer_cast<EventWorkspace>(integratedWorkspace);
EventWorkspace_sptr ew = boost::dynamic_pointer_cast<EventWorkspace>( integratedWorkspace);
if (ew)
throw std::runtime_error("Error! Integration output is not a Workspace2D.");

// i. set the workspace to have no units
// TODO remove the setup on integratedWorkspace after testing
// integratedWorkspace->isDistribution(false);
// integratedWorkspace->setYUnit("");
//set the workspace to have no units
integratedWorkspace->isDistribution(false);
integratedWorkspace->setYUnit("");

// 4. Read integrated and set up the Masking workspace
const double liveValue(1.0);
const double deadValue = 0.0;

// a) iterate over the data values setting the live and dead values
const double liveValue(1.0);
// iterate over the data values setting the live and dead values
const int numSpec = static_cast<int>(integratedWorkspace->getNumberHistograms());
const int progStep = static_cast<int>(ceil(numSpec / 100.0));

// b) Keep track of those reading low/high
// Keep track of those reading low/high
int numLow(0), numHigh(0);

// c) Loop
for (int i = 0; i < numSpec; ++i)
{
// i) Progress
if (i % progStep == 0)
{
progress(static_cast<double>(i)/numSpec);
interruption_point();
}

// ii) Set up value
const double & yValue = integratedWorkspace->readY(i)[0];
if ( yValue <= lowThreshold )
{
// Value failed
integratedWorkspace->maskWorkspaceIndex(i);

maskWS->dataY(i)[0] = deadValue;

++numLow;
continue;
}

} else if ( yValue >= highThreshold )
if ( yValue >= highThreshold )
{
// Value failed
integratedWorkspace->maskWorkspaceIndex(i);

maskWS->dataY(i)[0] = deadValue;

++numHigh;
continue;
}

} else {
// Value passed the tests, flag it.
integratedWorkspace->dataY(i)[0] = liveValue;

maskWS->dataY(i)[0] = liveValue;
} // ENDIFELSE
} // ENDFOR
// Value passed the tests, flag it.
integratedWorkspace->dataY(i)[0] = liveValue;
}

const int totalNumFailed(numLow+numHigh);
g_log.notice() << "Found a total of " << totalNumFailed << " failing "
Expand All @@ -177,10 +140,7 @@ namespace Mantid

setProperty("NumberOfFailures", totalNumFailed);
// Assign it to the output workspace property
// setProperty("OutputWorkspace2", integratedWorkspace);
setProperty("OutputWorkspace", maskWS);

return;
setProperty("OutputWorkspace", integratedWorkspace);
}

}
Expand Down
130 changes: 65 additions & 65 deletions Code/Mantid/Framework/Algorithms/src/MedianDetectorTest.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -46,7 +46,7 @@ namespace Mantid
new HistogramValidator<>),
"Name of the input workspace" );
declareProperty(
new WorkspaceProperty<DataObjects::SpecialWorkspace2D>("OutputWorkspace","",Direction::Output),
new WorkspaceProperty<>("OutputWorkspace","",Direction::Output),
"A MaskWorkspace where 0 denotes a masked spectra. Any spectra containing"
"a zero is also masked on the output");

Expand Down Expand Up @@ -90,39 +90,43 @@ namespace Mantid
*/
void MedianDetectorTest::exec()
{
// 1. Get input arguments
retrieveProperties();

// 2. Adds the counts from all the bins and puts them in one total bin


//Adds the counts from all the bins and puts them in one total bin
MatrixWorkspace_sptr counts = integrateSpectra(m_inputWS, m_minSpec, m_maxSpec,
m_rangeLower, m_rangeUpper, true);
counts = convertToRate(counts);

// FIXME: The next section that calculates the solid angle is commented out until
// the SolidAngle algorithm is corrected to return the correct number.
// (see http://trac.mantidproject.org/mantid/ticket/2596)

// FIXME: The next section that calculates the solid angle is commented out until
// the SolidAngle algorithm is corrected to return the correct number.
// (see http://trac.mantidproject.org/mantid/ticket/2596)
//MatrixWorkspace_sptr angles = getSolidAngles(m_minSpec, m_maxSpec);

//Gets the count rate per solid angle (in steradians), if it exists, for each spectrum
//this calculation is optional, it depends on angle information existing
//if ( angles.use_count() == 1 )
//{
//if some numbers in angles are zero we will get the infinity flag value
// in the output work space which needs to be dealt with later
//counts = counts/angles;
//if some numbers in angles are zero we will get the infinity flag value
// in the output work space which needs to be dealt with later
//counts = counts/angles;
//}

// End of Solid Angle commented out section


// 3. An average of the data, the median is less influenced by a small number of huge values than the mean
// An average of the data, the median is less influenced by a small number of huge values than the mean
std::set<int> badIndices;
double average = calculateMedian(counts, badIndices);

// 4. Initialize output masking workspace
Geometry::Instrument_const_sptr myinstrument = m_inputWS->getInstrument();
DataObjects::SpecialWorkspace2D_sptr maskWS(new DataObjects::SpecialWorkspace2D(myinstrument));
// Create a workspace for the output
MatrixWorkspace_sptr maskWS = WorkspaceFactory::Instance().create(counts);
// Make sure the output is simple
maskWS->isDistribution(false);
maskWS->setYUnit("");
maskWS->instrumentParameters().clear();

int numFailed = doDetectorTests(counts, maskWS, average, badIndices);

Expand Down Expand Up @@ -223,8 +227,8 @@ namespace Mantid
* @return The number of detectors that failed the tests, not including those skipped
*/
int MedianDetectorTest::doDetectorTests(const API::MatrixWorkspace_sptr countWorkspace,
DataObjects::SpecialWorkspace2D_sptr maskWS,
const double average, const std::set<int> & badIndices)
API::MatrixWorkspace_sptr maskWS,
const double average, const std::set<int> & badIndices)
{
g_log.information("Applying the criteria to find failing detectors");

Expand All @@ -238,87 +242,83 @@ namespace Mantid
const int numSpec(m_maxSpec - m_minSpec);
const int progStep = static_cast<int>(ceil(numSpec/30.0));

const double live_value(1.0);
int numLow(0), numHigh(0);
size_t numBadPixels = 0;

//set the workspace to have no units
countWorkspace->isDistribution(false);
countWorkspace->setYUnit("");

double liveValue = 1.0;
double deadValue = 0.0;

PARALLEL_FOR1(countWorkspace)
for (int i = 0; i <= numSpec; ++i)
{
PARALLEL_START_INTERUPT_REGION

// i. update the progressbar information
// update the progressbar information
if (i % progStep == 0)
{
progress(advanceProgress(progStep*static_cast<double>(RTMarkDetects)/numSpec));
}


// ii. If the value is not in the badIndices set then assume it is good
// else skip tests for it
// If the value is not in the badIndices set then assume it is good
// else skip tests for it
if( badIndices.count(i) == 1 )
{
// Bad pixel
maskWS->dataY(i)[0] = deadValue;
numBadPixels ++;

} else {
//Good pixel: do tests

const double sig = minSigma*countWorkspace->readE(i)[0];
maskWS->maskWorkspaceIndex(i);
continue;
}
//Do tests
const double sig = minSigma*countWorkspace->readE(i)[0];
// Check the significance value is okay
if( boost::math::isinf(std::abs(sig)) || boost::math::isinf(sig) )
{
PARALLEL_CRITICAL(MedianDetectorTest_failed_a)
{
maskWS->maskWorkspaceIndex(i);
++numLow;
}
continue;
}

// (a) Check the significance value is okay
if( boost::math::isinf(std::abs(sig)) || boost::math::isinf(sig) )
const double yIn = countWorkspace->dataY(i)[0];
if ( yIn <= lowLim )
{
// compare the difference against the size of the errorbar -statistical significance check
if(average - yIn > sig)
{
// sig is infinity
PARALLEL_CRITICAL(MedianDetectorTest_failed_a)
PARALLEL_CRITICAL(MedianDetectorTest_failed_b)
{
maskWS->dataY(i)[0] = deadValue;
maskWS->maskWorkspaceIndex(i);
++numLow;
}
} else {
// sig is not infinity
const double yIn = countWorkspace->dataY(i)[0];
if ( (yIn <= lowLim) && (average - yIn > sig))
continue;
}
}
if (yIn >= highLim)
{
// compare the difference against the size of the errorbar -statistical significance check
if(yIn - average > sig)
{
PARALLEL_CRITICAL(MedianDetectorTest_failed_c)
{
// compare the difference against the size of the errorbar -statistical significance check
PARALLEL_CRITICAL(MedianDetectorTest_failed_b)
{
maskWS->dataY(i)[0] = deadValue;
++numLow;
}
maskWS->maskWorkspaceIndex(i);
++numHigh;
}
else if ( (yIn >= highLim) && (yIn - average > sig) )
{
// compare the difference against the size of the errorbar -statistical significance check
PARALLEL_CRITICAL(MedianDetectorTest_failed_c)
{
maskWS->dataY(i)[0] = deadValue;
++numHigh;
}
} else {
// Reaching here passes the tests
maskWS->dataY(i)[0] = liveValue;

} // IF-ELSE: test
} // IF-ELSE: Sig is infinity?
} // IF-ELSE: Bad/good pixels
continue;
}
}
// Reaching here passes the tests
maskWS->dataY(i)[0] = live_value;

PARALLEL_END_INTERUPT_REGION
} // ENDFOR
}
PARALLEL_CHECK_INTERUPT_REGION

// Log finds
g_log.information() << "-- Detector tests --\n "
<< "Number recording low: " << numLow << "\n"
<< "Number recording high: " << numHigh << "\n"
<< "Number bad pixels: << " << numBadPixels << std::endl;
<< "Number recording high: " << numHigh << "\n";

return (numLow + numHigh);
}
Expand Down

0 comments on commit 1fb1684

Please sign in to comment.