use of uk.ac.sussex.gdsc.smlm.results.count.FailCounter in project GDSC-SMLM by aherbert.
the class FailCountManager method analyseData.
private void analyseData() {
final LocalList<FailCountData> failCountData = failCountDataRef.get();
if (failCountData.isEmpty()) {
IJ.error(TITLE, "No fail count data in memory");
return;
}
if (!showAnalysisDialog()) {
return;
}
final int maxCons = getMaxConsecutiveFailCount(failCountData);
final int maxFail = getMaxFailCount(failCountData);
// Create a set of fail counters
final LocalList<FailCounter> counters = new LocalList<>();
final TByteArrayList type = new TByteArrayList();
for (int i = 0; i <= maxCons; i++) {
counters.add(ConsecutiveFailCounter.create(i));
}
fill(type, counters, 0);
// The other counters are user configured.
// Ideally this would be a search to optimise the best parameters
// for each counter as any enumeration may be way off the mark.
// Note that 0 failures in a window can be scored using the consecutive fail counter.
int max = Math.min(maxFail, settings.getRollingCounterMaxAllowedFailures());
for (int fail = MathUtils.clip(1, maxFail, settings.getRollingCounterMinAllowedFailures()); fail <= max; fail++) {
// Note that n-1 failures in window n can be scored using the consecutive fail counter.
for (int window = Math.max(fail + 2, settings.getRollingCounterMinWindow()); window <= settings.getRollingCounterMaxWindow(); window++) {
counters.add(RollingWindowFailCounter.create(fail, window));
}
switch(checkCounters(counters)) {
case ANALYSE:
break;
case CONTINUE:
break;
case RETURN:
return;
default:
throw new IllegalStateException();
}
}
fill(type, counters, 1);
max = Math.min(maxFail, settings.getWeightedCounterMaxAllowedFailures());
for (int fail = MathUtils.min(maxFail, settings.getWeightedCounterMinAllowedFailures()); fail <= max; fail++) {
for (int w = settings.getWeightedCounterMinPassDecrement(); w <= settings.getWeightedCounterMaxPassDecrement(); w++) {
counters.add(WeightedFailCounter.create(fail, 1, w));
}
switch(checkCounters(counters)) {
case ANALYSE:
break;
case CONTINUE:
break;
case RETURN:
return;
default:
throw new IllegalStateException();
}
}
fill(type, counters, 2);
max = Math.min(maxFail, settings.getResettingCounterMaxAllowedFailures());
for (int fail = MathUtils.min(maxFail, settings.getResettingCounterMinAllowedFailures()); fail <= max; fail++) {
for (double f = settings.getResettingCounterMinResetFraction(); f <= settings.getResettingCounterMaxResetFraction(); f += settings.getResettingCounterIncResetFraction()) {
counters.add(ResettingFailCounter.create(fail, f));
}
switch(checkCounters(counters)) {
case ANALYSE:
break;
case CONTINUE:
break;
case RETURN:
return;
default:
throw new IllegalStateException();
}
}
fill(type, counters, 3);
for (int count = settings.getPassRateCounterMinAllowedCounts(); count <= settings.getPassRateCounterMaxAllowedCounts(); count++) {
for (double f = settings.getPassRateCounterMinPassRate(); f <= settings.getPassRateCounterMaxPassRate(); f += settings.getPassRateCounterIncPassRate()) {
counters.add(PassRateFailCounter.create(count, f));
}
switch(checkCounters(counters)) {
case ANALYSE:
break;
case CONTINUE:
break;
case RETURN:
return;
default:
throw new IllegalStateException();
}
}
fill(type, counters, 4);
counters.trimToSize();
// Score each of a set of standard fail counters against each frame using how
// close they are to the target.
final double[] score = new double[counters.size()];
final double targetPassFraction = settings.getTargetPassFraction();
final int nThreads = Prefs.getThreads();
final ExecutorService executor = Executors.newFixedThreadPool(nThreads);
final LocalList<Future<?>> futures = new LocalList<>(nThreads);
final Ticker ticker = ImageJUtils.createTicker(failCountData.size(), nThreads);
IJ.showStatus("Analysing " + TextUtils.pleural(counters.size(), "counter"));
for (int i = 0; i < failCountData.size(); i++) {
final FailCountData data = failCountData.unsafeGet(i);
futures.add(executor.submit(() -> {
if (IJ.escapePressed()) {
return;
}
// TODO - Ideally this plugin should be run on benchmark data with ground truth.
// The target could be to ensure all the correct results are fit
// and false positives are excluded from incrementing the pass counter.
// This could be done by saving the results from a benchmarking scoring
// plugin to memory as the current dataset.
data.initialiseAnalysis(targetPassFraction);
// Score in blocks and then do a synchronized write to the combined score
final Thread t = Thread.currentThread();
final double[] s = new double[8192];
int index = 0;
while (index < counters.size()) {
if (t.isInterrupted()) {
break;
}
final int block = Math.min(8192, counters.size() - index);
for (int j = 0; j < block; j++) {
final FailCounter counter = counters.unsafeGet(index + j).newCounter();
s[j] = data.score(counter);
}
// Write to the combined score
synchronized (score) {
for (int j = 0; j < block; j++) {
score[index + j] += s[j];
}
}
index += block;
}
ticker.tick();
}));
}
executor.shutdown();
ConcurrencyUtils.waitForCompletionUnchecked(futures);
IJ.showProgress(1);
if (IJ.escapePressed()) {
IJ.showStatus("");
IJ.error(TITLE, "Cancelled analysis");
return;
}
IJ.showStatus("Summarising results ...");
// TODO - check if the top filter is at the bounds of the range
final int minIndex = SimpleArrayUtils.findMinIndex(score);
ImageJUtils.log(TITLE + " Analysis : Best counter = %s (Score = %f)", counters.unsafeGet(minIndex).getDescription(), score[minIndex]);
// Show a table of results for the top N for each type
final int topN = Math.min(settings.getTableTopN(), score.length);
if (topN > 0) {
final byte[] types = type.toArray();
final byte maxType = types[types.length - 1];
final TextWindow resultsWindow = createTable();
for (byte b = 0; b <= maxType; b++) {
int[] indices;
// Use a heap to avoid a full sort
final IntDoubleMinHeap heap = new IntDoubleMinHeap(topN);
for (int i = 0; i < score.length; i++) {
if (types[i] == b) {
heap.offer(score[i], i);
}
}
if (heap.getSize() == 0) {
continue;
}
indices = heap.getItems();
// Ensure sorted
SortUtils.sortIndices(indices, score, false);
final StringBuilder sb = new StringBuilder();
try (final BufferedTextWindow tw = new BufferedTextWindow(resultsWindow)) {
for (int i = 0; i < topN; i++) {
sb.setLength(0);
final int j = indices[i];
sb.append(i + 1).append('\t');
sb.append(counters.unsafeGet(j).getDescription()).append('\t');
sb.append(score[j]);
tw.append(sb.toString());
}
}
}
}
// TODO - Save the best fail counter to the current fit configuration.
IJ.showStatus("");
}
use of uk.ac.sussex.gdsc.smlm.results.count.FailCounter in project GDSC-SMLM by aherbert.
the class FitWorker method run.
/**
* Locate all the peaks in the image specified by the fit job.
*
* <p>WARNING: The FitWorker fits a sub-region of the data for each maxima. It then updates the
* FitResult parameters with an offset reflecting the position. The initialParameters are not
* updated with this offset unless configured.
*
* @param job The fit job
*/
public void run(FitJob job) {
final long start = System.nanoTime();
job.start();
this.job = job;
benchmarking = false;
this.slice = job.slice;
// Used for debugging
// if (logger == null) logger = new gdsc.fitting.logging.ConsoleLogger();
// Crop to the ROI
cc = new CoordinateConverter(job.bounds);
// Note if the bounds change for efficient caching.
newBounds = !cc.dataBounds.equals(lastBounds);
if (newBounds) {
lastBounds = cc.dataBounds;
}
final int width = cc.dataBounds.width;
final int height = cc.dataBounds.height;
borderLimitX = width - border;
borderLimitY = height - border;
data = job.data;
// This is tied to the input data
dataEstimator = null;
// relative to the global origin.
if (isFitCameraCounts) {
cameraModel.removeBias(cc.dataBounds, data);
} else {
cameraModel.removeBiasAndGain(cc.dataBounds, data);
}
final FitParameters params = job.getFitParameters();
this.endT = (params != null) ? params.endT : -1;
candidates = indentifySpots(job, width, height, params);
if (candidates.getSize() == 0) {
finishJob(job, start);
return;
}
fittedBackground = new Statistics();
// Always get the noise and store it with the results.
if (params != null && !Float.isNaN(params.noise)) {
noise = params.noise;
fitConfig.setNoise(noise);
} else if (calculateNoise) {
noise = estimateNoise();
fitConfig.setNoise(noise);
}
// System.out.printf("Slice %d : Noise = %g\n", slice, noise);
if (logger != null) {
LoggerUtils.log(logger, Level.INFO, "Slice %d: Noise = %f", slice, noise);
}
final ImageExtractor ie = ImageExtractor.wrap(data, width, height);
double[] region = null;
final float offsetx = cc.dataBounds.x;
final float offsety = cc.dataBounds.y;
if (params != null && params.fitTask == FitTask.MAXIMA_IDENITIFICATION) {
final float sd0 = (float) xsd;
final float sd1 = (float) ysd;
for (int n = 0; n < candidates.getSize(); n++) {
// Find the background using the perimeter of the data.
// TODO - Perhaps the Gaussian Fitter should be used to produce the initial estimates but no
// actual fit done.
// This would produce coords using the centre-of-mass.
final Candidate candidate = candidates.get(n);
int x = candidate.x;
int y = candidate.y;
final Rectangle regionBounds = ie.getBoxRegionBounds(x, y, fitting);
region = ie.crop(regionBounds, region);
final float b = (float) Gaussian2DFitter.getBackground(region, regionBounds.width, regionBounds.height, 1);
// Offset the coords to the centre of the pixel. Note the bounds will be added later.
// Subtract the background to get the amplitude estimate then convert to signal.
final float amplitude = candidate.intensity - ((relativeIntensity) ? 0 : b);
final float signal = (float) (amplitude * 2.0 * Math.PI * sd0 * sd1);
final int index = y * width + x;
x += offsetx;
y += offsety;
final float[] peakParams = new float[1 + Gaussian2DFunction.PARAMETERS_PER_PEAK];
peakParams[Gaussian2DFunction.BACKGROUND] = b;
peakParams[Gaussian2DFunction.SIGNAL] = signal;
peakParams[Gaussian2DFunction.X_POSITION] = x + 0.5f;
peakParams[Gaussian2DFunction.Y_POSITION] = y + 0.5f;
// peakParams[Gaussian2DFunction.Z_POSITION] = 0;
peakParams[Gaussian2DFunction.X_SD] = sd0;
peakParams[Gaussian2DFunction.Y_SD] = sd1;
// peakParams[Gaussian2DFunction.ANGLE] = 0;
final float u = (float) Gaussian2DPeakResultHelper.getMeanSignalUsingP05(signal, sd0, sd1);
sliceResults.add(createResult(x, y, data[index], 0, noise, u, peakParams, null, n, 0));
}
} else {
initialiseFitting();
// Smooth the data to provide initial background estimates
final float[] smoothedData = backgroundSmoothing.process(data, width, height);
final ImageExtractor ie2 = ImageExtractor.wrap(smoothedData, width, height);
// Perform the Gaussian fit
// The SpotFitter is used to create a dynamic MultiPathFitResult object.
// This is then passed to a multi-path filter. Thus the same fitting decision process
// is used when benchmarking and when running on actual data.
// Note: The SpotFitter labels each PreprocessedFitResult using the offset in the FitResult
// object.
// The initial params and deviations can then be extracted for the results that pass the
// filter.
MultiPathFilter filter;
final IMultiPathFitResults multiPathResults = this;
final SelectedResultStore store = this;
coordinateStore = coordinateStore.resize(cc.dataBounds.x, cc.dataBounds.y, width, height);
if (params != null && params.fitTask == FitTask.BENCHMARKING) {
// Run filtering as normal. However in the event that a candidate is missed or some
// results are not generated we must generate them. This is done in the complete(int)
// method if we set the benchmarking flag.
benchmarking = true;
// Filter using the benchmark filter
filter = params.benchmarkFilter;
if (filter == null) {
// Create a default filter using the standard FitConfiguration to ensure sensible fits
// are stored as the current slice results.
// Note the current fit configuration for benchmarking may have minimal filtering settings
// so we do not use that object.
final FitConfiguration tmp = new FitConfiguration();
final double residualsThreshold = 0.4;
filter = new MultiPathFilter(tmp, createMinimalFilter(PrecisionMethod.POISSON_CRLB), residualsThreshold);
}
} else {
// Filter using the configuration.
if (this.filter == null) {
// This can be cached. Q. Clone the config?
this.filter = new MultiPathFilter(fitConfig, createMinimalFilter(fitConfig.getPrecisionMethod()), config.getResidualsThreshold());
}
filter = this.filter;
}
// If we are benchmarking then do not generate results dynamically since we will store all
// results in the fit job.
dynamicMultiPathFitResult = new DynamicMultiPathFitResult(ie, ie2, !benchmarking);
// dynamicMultiPathFitResult = new DynamicMultiPathFitResult(ie, false);
// The local background computation is only required for the precision method.
// Also compute it when benchmarking.
localBackground = benchmarking || fitConfig.getPrecisionMethodValue() == PrecisionMethod.MORTENSEN_LOCAL_BACKGROUND_VALUE;
// Debug where the fit config may be different between benchmarking and fitting
if (slice == -1) {
fitConfig.initialise(1, 1, 1);
final String newLine = System.lineSeparator();
final String tmpdir = System.getProperty("java.io.tmpdir");
try (BufferedWriter writer = Files.newBufferedWriter(Paths.get(tmpdir, String.format("config.%d.txt", slice)))) {
JsonFormat.printer().appendTo(config.getFitEngineSettings(), writer);
} catch (final IOException ex) {
logger.log(Level.SEVERE, "Unable to write message", ex);
}
FileUtils.save(Paths.get(tmpdir, String.format("filter.%d.xml", slice)).toString(), filter.toXml());
// filter.setDebugFile(String.format("/tmp/fitWorker.%b.txt", benchmarking));
final StringBuilder sb = new StringBuilder();
sb.append((benchmarking) ? ((uk.ac.sussex.gdsc.smlm.results.filter.Filter) filter.getFilter()).toXml() : fitConfig.getSmartFilterString()).append(newLine);
sb.append(((uk.ac.sussex.gdsc.smlm.results.filter.Filter) filter.getMinimalFilter()).toXml()).append(newLine);
sb.append(filter.residualsThreshold).append(newLine);
sb.append(config.getFailuresLimit()).append(newLine);
sb.append(config.getDuplicateDistance()).append(":");
sb.append(config.getDuplicateDistanceAbsolute()).append(newLine);
if (spotFilter != null) {
sb.append(spotFilter.getDescription()).append(newLine);
}
sb.append("MaxCandidate = ").append(candidates.getSize()).append(newLine);
for (int i = 0, len = candidates.getLength(); i < len; i++) {
TextUtils.formatTo(sb, "Fit %d [%d,%d = %.1f]%n", i, candidates.get(i).x, candidates.get(i).y, candidates.get(i).intensity);
}
FileUtils.save(Paths.get(tmpdir, String.format("candidates.%d.xml", slice)).toString(), sb.toString());
}
FailCounter failCounter = config.getFailCounter();
if (!benchmarking && params != null && params.pass != null) {
// We want to store the pass/fail for consecutive candidates
params.pass = new boolean[candidates.getLength()];
failCounter = new RecordingFailCounter(params.pass, failCounter);
filter.select(multiPathResults, failCounter, true, store, coordinateStore);
} else {
filter.select(multiPathResults, failCounter, true, store, coordinateStore);
}
// Note: We go deeper into the candidate list than max candidate
// for any candidate where we have a good fit result as an estimate.
// Q. Should this only be for benchmarking?
// if (benchmarking)
// System.out.printf("Slice %d: %d + %d\n", slice, dynamicMultiPathFitResult.extra,
// candidates.getSize());
// Create the slice results
final CandidateList fitted = gridManager.getFittedCandidates();
sliceResults.ensureCapacity(fitted.getSize());
for (int i = 0; i < fitted.getSize(); i++) {
if (fitted.get(i).fit) {
sliceResults.push(createResult(offsetx, offsety, fitted.get(i)));
}
}
if (logger != null) {
LoggerUtils.log(logger, Level.INFO, "Slice %d: %d / %d = %s", slice, success, candidates.getSize(), TextUtils.pleural(fitted.getSize(), "result"));
}
}
this.results.addAll(sliceResults);
finishJob(job, start);
}
use of uk.ac.sussex.gdsc.smlm.results.count.FailCounter in project GDSC-SMLM by aherbert.
the class FitEngineConfiguration method getFailCounter.
/**
* Gets the fail counter.
*
* @return the fail counter
*/
public FailCounter getFailCounter() {
// manually passed in. Either would allow more complex fail counters to be used.
if (failCounter == null) {
final int failuresLimit = getFailuresLimit();
final FailCounter f1 = // Negative will disable
(failuresLimit >= 0) ? ConsecutiveFailCounter.create(failuresLimit) : null;
final double passRate = getPassRate();
// TODO - the allowed counts could be an input
final FailCounter f2 = (passRate > 0) ? PassRateFailCounter.create(5, passRate) : null;
// All fail counters must pass to continue fitting
failCounter = CombinedAndFailCounter.join(f1, f2);
if (failCounter == null) {
failCounter = NullFailCounter.INSTANCE;
}
}
return failCounter.newCounter();
}
Aggregations