Search in sources :

Example 11 with Calibration

use of uk.ac.sussex.gdsc.smlm.data.config.CalibrationProtos.Calibration in project GDSC-SMLM by aherbert.

the class TrackPopulationAnalysis method run.

@Override
public void run(String arg) {
    SmlmUsageTracker.recordPlugin(this.getClass(), arg);
    if (MemoryPeakResults.isMemoryEmpty()) {
        IJ.error(TITLE, "No localisations in memory");
        return;
    }
    settings = Settings.load();
    // Saved by reference so just save now
    settings.save();
    // Read in multiple traced datasets
    // All datasets must have the same pixel pitch and exposure time
    // Get parameters
    // Convert datasets to tracks
    // For each track compute the 4 local track features using the configured window
    // 
    // Optional:
    // Fit a multi-variate Gaussian mixture model to the data
    // (using the configured number of components/populations)
    // Assign each point in the track using the model.
    // Smooth the assignments.
    // 
    // The alternative is to use the localisation category to assign populations.
    // 
    // Plot histograms of each track parameter, coloured by component
    final List<MemoryPeakResults> combinedResults = new LocalList<>();
    if (!showInputDialog(combinedResults)) {
        return;
    }
    final boolean hasCategory = showHasCategoryDialog(combinedResults);
    if (!showDialog(hasCategory)) {
        return;
    }
    ImageJUtils.log(TITLE + "...");
    final List<Trace> tracks = getTracks(combinedResults, settings.window, settings.minTrackLength);
    if (tracks.isEmpty()) {
        IJ.error(TITLE, "No tracks. Please check the input data and min track length setting.");
        return;
    }
    final Calibration cal = combinedResults.get(0).getCalibration();
    final CalibrationReader cr = new CalibrationReader(cal);
    // Use micrometer / second
    final TypeConverter<DistanceUnit> distanceConverter = cr.getDistanceConverter(DistanceUnit.UM);
    final double exposureTime = cr.getExposureTime() / 1000.0;
    final Pair<int[], double[][]> trackData = extractTrackData(tracks, distanceConverter, exposureTime, hasCategory);
    final double[][] data = trackData.getValue();
    // Histogram the raw data.
    final Array2DRowRealMatrix raw = new Array2DRowRealMatrix(data, false);
    final WindowOrganiser wo = new WindowOrganiser();
    // Store the histogram data for plotting the components
    final double[][] columns = new double[FEATURE_NAMES.length][];
    final double[][] limits = new double[FEATURE_NAMES.length][];
    // Get column data
    for (int i = 0; i < FEATURE_NAMES.length; i++) {
        columns[i] = raw.getColumn(i);
        if (i == FEATURE_D) {
            // Plot using a logarithmic scale
            SimpleArrayUtils.apply(columns[i], Math::log10);
        }
        limits[i] = MathUtils.limits(columns[i]);
    }
    // Compute histogram bins
    final int[] bins = new int[FEATURE_NAMES.length];
    if (settings.histogramBins > 0) {
        Arrays.fill(bins, settings.histogramBins);
    } else {
        for (int i = 0; i < FEATURE_NAMES.length; i++) {
            bins[i] = HistogramPlot.getBins(StoredData.create(columns[i]), BinMethod.FD);
        }
        // Use the maximum so all histograms look the same
        Arrays.fill(bins, MathUtils.max(bins));
    }
    // Compute plots
    final Plot[] plots = new Plot[FEATURE_NAMES.length];
    for (int i = 0; i < FEATURE_NAMES.length; i++) {
        final double[][] hist = HistogramPlot.calcHistogram(columns[i], limits[i][0], limits[i][1], bins[i]);
        plots[i] = new Plot(TITLE + " " + FEATURE_NAMES[i], getFeatureLabel(i, i == FEATURE_D), "Frequency");
        plots[i].addPoints(hist[0], hist[1], Plot.BAR);
        ImageJUtils.display(plots[i].getTitle(), plots[i], ImageJUtils.NO_TO_FRONT, wo);
    }
    wo.tile();
    // The component for each data point
    int[] component;
    // The number of components
    int numComponents;
    // Data used to fit the Gaussian mixture model
    double[][] fitData;
    // The fitted model
    MixtureMultivariateGaussianDistribution model;
    if (hasCategory) {
        // Use the category as the component.
        // No fit data and no output model
        fitData = null;
        model = null;
        // The component is stored at the end of the raw track data.
        final int end = data[0].length - 1;
        component = Arrays.stream(data).mapToInt(d -> (int) d[end]).toArray();
        numComponents = MathUtils.max(component) + 1;
        // In the EM algorithm the probability of each data point is computed and normalised to
        // sum to 1. The normalised probabilities are averaged to create the weights.
        // Note the probability of each data point uses the previous weight and the algorithm
        // iterates.
        // This is not a fitted model but the input model so use
        // zero weights to indicate no fitting was performed.
        final double[] weights = new double[numComponents];
        // Remove the trailing component to show the 'model' in a table.
        createModelTable(Arrays.stream(data).map(d -> Arrays.copyOf(d, end)).toArray(double[][]::new), weights, component);
    } else {
        // Multivariate Gaussian mixture EM
        // Provide option to not use the anomalous exponent in the population mix.
        int sortDimension = SORT_DIMENSION;
        if (settings.ignoreAlpha) {
            // Remove index 0. This shifts the sort dimension.
            sortDimension--;
            fitData = Arrays.stream(data).map(d -> Arrays.copyOfRange(d, 1, d.length)).toArray(double[][]::new);
        } else {
            fitData = SimpleArrayUtils.deepCopy(data);
        }
        final MultivariateGaussianMixtureExpectationMaximization mixed = fitGaussianMixture(fitData, sortDimension);
        if (mixed == null) {
            IJ.error(TITLE, "Failed to fit a mixture model");
            return;
        }
        model = sortComponents(mixed.getFittedModel(), sortDimension);
        // For the best model, assign to the most likely population.
        component = assignData(fitData, model);
        // Table of the final model using the original data (i.e. not normalised)
        final double[] weights = model.getWeights();
        numComponents = weights.length;
        createModelTable(data, weights, component);
    }
    // Output coloured histograms of the populations.
    final LUT lut = LutHelper.createLut(settings.lutIndex);
    IntFunction<Color> colourMap;
    if (LutHelper.getColour(lut, 0).equals(Color.BLACK)) {
        colourMap = i -> LutHelper.getNonZeroColour(lut, i, 0, numComponents - 1);
    } else {
        colourMap = i -> LutHelper.getColour(lut, i, 0, numComponents - 1);
    }
    for (int i = 0; i < FEATURE_NAMES.length; i++) {
        // Extract the data for each component
        final double[] col = columns[i];
        final Plot plot = plots[i];
        for (int n = 0; n < numComponents; n++) {
            final StoredData feature = new StoredData();
            for (int j = 0; j < component.length; j++) {
                if (component[j] == n) {
                    feature.add(col[j]);
                }
            }
            if (feature.size() == 0) {
                continue;
            }
            final double[][] hist = HistogramPlot.calcHistogram(feature.values(), limits[i][0], limits[i][1], bins[i]);
            // Colour the points
            plot.setColor(colourMap.apply(n));
            plot.addPoints(hist[0], hist[1], Plot.BAR);
        }
        plot.updateImage();
    }
    createTrackDataTable(tracks, trackData, fitData, model, component, cal, colourMap);
// Analysis.
// Assign the original localisations to their track component.
// Q. What about the start/end not covered by the window?
// Save tracks as a dataset labelled with the sub-track ID?
// Output for the bound component and free components track parameters.
// Compute dwell times.
// Other ...
// Track analysis plugin:
// Extract all continuous segments of the same component.
// Produce MSD plot with error bars.
// Fit using FBM model.
}
Also used : LocalList(uk.ac.sussex.gdsc.core.utils.LocalList) Array2DRowRealMatrix(org.apache.commons.math3.linear.Array2DRowRealMatrix) MixtureMultivariateGaussianDistribution(uk.ac.sussex.gdsc.smlm.math3.distribution.fitting.MultivariateGaussianMixtureExpectationMaximization.MixtureMultivariateGaussianDistribution) MemoryPeakResults(uk.ac.sussex.gdsc.smlm.results.MemoryPeakResults) DistanceUnit(uk.ac.sussex.gdsc.smlm.data.config.UnitProtos.DistanceUnit) Plot(ij.gui.Plot) HistogramPlot(uk.ac.sussex.gdsc.core.ij.HistogramPlot) Color(java.awt.Color) LUT(ij.process.LUT) Calibration(uk.ac.sussex.gdsc.smlm.data.config.CalibrationProtos.Calibration) WindowOrganiser(uk.ac.sussex.gdsc.core.ij.plugin.WindowOrganiser) CalibrationReader(uk.ac.sussex.gdsc.smlm.data.config.CalibrationReader) Trace(uk.ac.sussex.gdsc.smlm.results.Trace) MultivariateGaussianMixtureExpectationMaximization(uk.ac.sussex.gdsc.smlm.math3.distribution.fitting.MultivariateGaussianMixtureExpectationMaximization) StoredData(uk.ac.sussex.gdsc.core.utils.StoredData)

Example 12 with Calibration

use of uk.ac.sussex.gdsc.smlm.data.config.CalibrationProtos.Calibration in project GDSC-SMLM by aherbert.

the class ConvertResults method showDialog.

private static boolean showDialog(MemoryPeakResults results) {
    final ExtendedGenericDialog gd = new ExtendedGenericDialog(TITLE);
    gd.addMessage("Convert the current units for the results");
    gd.addHelp(HelpUrls.getUrl("convert-results"));
    final CalibrationReader cr = CalibrationWriter.create(results.getCalibration());
    gd.addChoice("Distance_unit", SettingsManager.getDistanceUnitNames(), cr.getDistanceUnitValue());
    gd.addNumericField("Calibration (nm/px)", cr.getNmPerPixel(), 2);
    gd.addChoice("Intensity_unit", SettingsManager.getIntensityUnitNames(), cr.getIntensityUnitValue());
    gd.addNumericField("Gain (Count/photon)", cr.getCountPerPhoton(), 2);
    gd.addChoice("Angle_unit", SettingsManager.getAngleUnitNames(), cr.getAngleUnitValue());
    gd.showDialog();
    if (gd.wasCanceled()) {
        return false;
    }
    final CalibrationWriter cw = results.getCalibrationWriterSafe();
    final DistanceUnit distanceUnit = SettingsManager.getDistanceUnitValues()[gd.getNextChoiceIndex()];
    cw.setNmPerPixel(Math.abs(gd.getNextNumber()));
    final IntensityUnit intensityUnit = SettingsManager.getIntensityUnitValues()[gd.getNextChoiceIndex()];
    cw.setCountPerPhoton(Math.abs(gd.getNextNumber()));
    final AngleUnit angleUnit = SettingsManager.getAngleUnitValues()[gd.getNextChoiceIndex()];
    // Don't set the calibration with bad values
    if (distanceUnit.getNumber() > 0 && !(cw.getNmPerPixel() > 0)) {
        IJ.error(TITLE, "Require positive nm/pixel for conversion");
        return false;
    }
    if (intensityUnit.getNumber() > 0 && !(cw.getCountPerPhoton() > 0)) {
        IJ.error(TITLE, "Require positive Count/photon for conversion");
        return false;
    }
    final Calibration newCalibration = cw.getCalibration();
    results.setCalibration(newCalibration);
    if (!results.convertToUnits(distanceUnit, intensityUnit, angleUnit)) {
        IJ.error(TITLE, "Conversion failed");
        return false;
    }
    return true;
}
Also used : AngleUnit(uk.ac.sussex.gdsc.smlm.data.config.UnitProtos.AngleUnit) CalibrationWriter(uk.ac.sussex.gdsc.smlm.data.config.CalibrationWriter) IntensityUnit(uk.ac.sussex.gdsc.smlm.data.config.UnitProtos.IntensityUnit) ExtendedGenericDialog(uk.ac.sussex.gdsc.core.ij.gui.ExtendedGenericDialog) Calibration(uk.ac.sussex.gdsc.smlm.data.config.CalibrationProtos.Calibration) CalibrationReader(uk.ac.sussex.gdsc.smlm.data.config.CalibrationReader) DistanceUnit(uk.ac.sussex.gdsc.smlm.data.config.UnitProtos.DistanceUnit)

Example 13 with Calibration

use of uk.ac.sussex.gdsc.smlm.data.config.CalibrationProtos.Calibration in project GDSC-SMLM by aherbert.

the class CalibrateResults method showDialog.

private boolean showDialog(MemoryPeakResults results) {
    final ExtendedGenericDialog gd = new ExtendedGenericDialog(TITLE);
    gd.addHelp(HelpUrls.getUrl("calibrate-results"));
    final Calibration oldCalibration = results.getCalibration();
    final boolean existingCalibration = oldCalibration != null;
    if (!existingCalibration) {
        gd.addMessage("No calibration found, using defaults");
    }
    final CalibrationWriter cw = results.getCalibrationWriterSafe();
    gd.addStringField("Name", results.getName(), Math.max(Math.min(results.getName().length(), 60), 20));
    if (existingCalibration) {
        gd.addCheckbox("Update_all_linked_results", settings.updateAll);
    }
    PeakFit.addCameraOptions(gd, PeakFit.FLAG_QUANTUM_EFFICIENCY, cw);
    gd.addNumericField("Calibration (nm/px)", cw.getNmPerPixel(), 2);
    gd.addNumericField("Exposure_time (ms)", cw.getExposureTime(), 2);
    gd.showDialog();
    if (gd.wasCanceled()) {
        return false;
    }
    final String name = gd.getNextString();
    if (!results.getName().equals(name)) {
        // If the name is changed then remove and add back to memory
        final MemoryPeakResults existingResults = MemoryPeakResults.removeResults(results.getName());
        if (existingResults != null) {
            results = existingResults;
            results.setName(name);
            MemoryPeakResults.addResults(results);
        }
    }
    if (existingCalibration) {
        settings.updateAll = gd.getNextBoolean();
    }
    cw.setCameraType(SettingsManager.getCameraTypeValues()[gd.getNextChoiceIndex()]);
    cw.setNmPerPixel(Math.abs(gd.getNextNumber()));
    cw.setExposureTime(Math.abs(gd.getNextNumber()));
    gd.collectOptions();
    final Calibration newCalibration = cw.getCalibration();
    results.setCalibration(newCalibration);
    if (settings.updateAll) {
        // be missed.
        for (final MemoryPeakResults r : MemoryPeakResults.getAllResults()) {
            if (r.getCalibration() == oldCalibration) {
                r.setCalibration(newCalibration);
            }
        }
    }
    return true;
}
Also used : CalibrationWriter(uk.ac.sussex.gdsc.smlm.data.config.CalibrationWriter) MemoryPeakResults(uk.ac.sussex.gdsc.smlm.results.MemoryPeakResults) ExtendedGenericDialog(uk.ac.sussex.gdsc.core.ij.gui.ExtendedGenericDialog) Calibration(uk.ac.sussex.gdsc.smlm.data.config.CalibrationProtos.Calibration)

Example 14 with Calibration

use of uk.ac.sussex.gdsc.smlm.data.config.CalibrationProtos.Calibration in project GDSC-SMLM by aherbert.

the class PeakFit method addCameraOptions.

/**
 * Adds the camera options.
 *
 * @param gd the dialog
 * @param options the options
 * @param calibrationProvider the calibration provider
 */
public static void addCameraOptions(final ExtendedGenericDialog gd, final int options, final CalibrationProvider calibrationProvider) {
    final CalibrationReader calibration = new CalibrationReader(calibrationProvider.getCalibration());
    gd.addChoice("Camera_type", SettingsManager.getCameraTypeNames(), CalibrationProtosHelper.getName(calibration.getCameraType()), new OptionListener<Integer>() {

        @Override
        public boolean collectOptions(Integer field) {
            final CalibrationWriter calibration = new CalibrationWriter(calibrationProvider.getCalibration());
            final CameraType t = SettingsManager.getCameraTypeValues()[field];
            if (calibration.getCameraType() != t) {
                calibration.setCameraType(t);
                calibrationProvider.saveCalibration(calibration.getCalibration());
            }
            return collectOptions(false);
        }

        @Override
        public boolean collectOptions() {
            return collectOptions(true);
        }

        private boolean collectOptions(boolean silent) {
            final CalibrationWriter calibration = new CalibrationWriter(calibrationProvider.getCalibration());
            final ExtendedGenericDialog egd = new ExtendedGenericDialog("Camera type options", null);
            if (calibration.isCcdCamera()) {
                egd.addNumericField("Camera_bias", calibration.getBias(), 2, 6, "Count");
                if (BitFlagUtils.anyNotSet(options, FLAG_NO_GAIN)) {
                    egd.addNumericField("Gain", calibration.getCountPerPhoton(), 4, 6, "Count/photon");
                }
                if (BitFlagUtils.anyNotSet(options, FLAG_NO_READ_NOISE)) {
                    egd.addNumericField("Read_noise", calibration.getReadNoise(), 4, 6, "Count");
                }
                if (BitFlagUtils.areSet(options, FLAG_QUANTUM_EFFICIENCY)) {
                    egd.addNumericField("Quantum_efficiency", calibration.getQuantumEfficiency(), 4, 6, "electron/photon");
                }
            } else if (calibration.isScmos()) {
                final String[] models = CameraModelManager.listCameraModels(true);
                egd.addChoice("Camera_model_name", models, calibration.getCameraModelName());
                if (BitFlagUtils.areSet(options, FLAG_QUANTUM_EFFICIENCY)) {
                    egd.addNumericField("Quantum_efficiency", calibration.getQuantumEfficiency(), 4, 6, "electron/photon");
                }
            } else {
                IJ.error("Unsupported camera type " + CalibrationProtosHelper.getName(calibration.getCameraType()));
                return false;
            }
            egd.setSilent(silent);
            egd.showDialog(true, gd);
            if (egd.wasCanceled()) {
                return false;
            }
            final Calibration old = calibration.getCalibration();
            if (calibration.isCcdCamera()) {
                calibration.setBias(Math.abs(egd.getNextNumber()));
                if (BitFlagUtils.anyNotSet(options, FLAG_NO_GAIN)) {
                    calibration.setCountPerPhoton(Math.abs(egd.getNextNumber()));
                }
                if (BitFlagUtils.anyNotSet(options, FLAG_NO_READ_NOISE)) {
                    calibration.setReadNoise(Math.abs(egd.getNextNumber()));
                }
                if (BitFlagUtils.areSet(options, FLAG_QUANTUM_EFFICIENCY)) {
                    calibration.setQuantumEfficiency(Math.abs(egd.getNextNumber()));
                }
            } else if (calibration.isScmos()) {
                // Note: Since this does not go through the FitConfiguration object the
                // camera model is not invalidated. However any code using this function
                // should later call configureFitSolver(...) which will set the camera model
                // using the camera model name.
                calibration.setCameraModelName(egd.getNextChoice());
                if (BitFlagUtils.areSet(options, FLAG_QUANTUM_EFFICIENCY)) {
                    calibration.setQuantumEfficiency(Math.abs(egd.getNextNumber()));
                }
            }
            final Calibration current = calibration.getCalibration();
            final boolean changed = !old.equals(current);
            if (changed) {
                calibrationProvider.saveCalibration(current);
            }
            return changed;
        }
    });
}
Also used : CalibrationWriter(uk.ac.sussex.gdsc.smlm.data.config.CalibrationWriter) ExtendedGenericDialog(uk.ac.sussex.gdsc.core.ij.gui.ExtendedGenericDialog) Calibration(uk.ac.sussex.gdsc.smlm.data.config.CalibrationProtos.Calibration) CalibrationReader(uk.ac.sussex.gdsc.smlm.data.config.CalibrationReader) CameraType(uk.ac.sussex.gdsc.smlm.data.config.CalibrationProtos.CameraType)

Example 15 with Calibration

use of uk.ac.sussex.gdsc.smlm.data.config.CalibrationProtos.Calibration in project GDSC-SMLM by aherbert.

the class SettingsManagerTest method canReadWriteConfiguration.

@SeededTest
void canReadWriteConfiguration(RandomSeed seed) throws IOException {
    final UniformRandomProvider rand = RngUtils.create(seed.getSeed());
    final Calibration.Builder builder = Calibration.newBuilder();
    builder.getCameraCalibrationBuilder().setBias(rand.nextDouble());
    final Calibration exp = builder.build();
    final String dir = SettingsManager.getSettingsDirectory();
    // System.out.println(dir);
    final File tmp = createTempDirectory(false);
    try {
        SettingsManager.setSettingsDirectory(tmp.getPath());
        Calibration obs = SettingsManager.readCalibration(SettingsManager.FLAG_SILENT | SettingsManager.FLAG_NO_DEFAULT);
        Assertions.assertTrue(obs == null, "Failed to read null");
        obs = SettingsManager.readCalibration(SettingsManager.FLAG_SILENT);
        Assertions.assertTrue(exp.getDefaultInstanceForType().equals(obs), "Failed to read default");
        Assertions.assertTrue(SettingsManager.writeSettings(exp), "Failed to write");
        obs = SettingsManager.readCalibration(0);
        Assertions.assertTrue(exp.equals(obs), "Not equal");
        SettingsManager.clearSettings(exp.getClass());
        obs = SettingsManager.readCalibration(SettingsManager.FLAG_SILENT | SettingsManager.FLAG_NO_DEFAULT);
        Assertions.assertTrue(obs == null, "Failed to clear");
    } finally {
        // Reset
        // Will work if the directory is empty
        tmp.delete();
        SettingsManager.setSettingsDirectory(dir);
    }
}
Also used : UniformRandomProvider(org.apache.commons.rng.UniformRandomProvider) Calibration(uk.ac.sussex.gdsc.smlm.data.config.CalibrationProtos.Calibration) File(java.io.File) SeededTest(uk.ac.sussex.gdsc.test.junit5.SeededTest)

Aggregations

Calibration (uk.ac.sussex.gdsc.smlm.data.config.CalibrationProtos.Calibration)17 CalibrationWriter (uk.ac.sussex.gdsc.smlm.data.config.CalibrationWriter)7 Rectangle (java.awt.Rectangle)5 ExtendedGenericDialog (uk.ac.sussex.gdsc.core.ij.gui.ExtendedGenericDialog)5 DistanceUnit (uk.ac.sussex.gdsc.smlm.data.config.UnitProtos.DistanceUnit)5 CalibrationReader (uk.ac.sussex.gdsc.smlm.data.config.CalibrationReader)4 PSF (uk.ac.sussex.gdsc.smlm.data.config.PSFProtos.PSF)4 MemoryPeakResults (uk.ac.sussex.gdsc.smlm.results.MemoryPeakResults)3 InvalidProtocolBufferException (com.google.protobuf.InvalidProtocolBufferException)2 Plot (ij.gui.Plot)2 Roi (ij.gui.Roi)2 LUT (ij.process.LUT)2 Color (java.awt.Color)2 EOFException (java.io.EOFException)2 IOException (java.io.IOException)2 NoSuchElementException (java.util.NoSuchElementException)2 Matcher (java.util.regex.Matcher)2 Pattern (java.util.regex.Pattern)2 Test (org.junit.jupiter.api.Test)2 WindowOrganiser (uk.ac.sussex.gdsc.core.ij.plugin.WindowOrganiser)2