Search in sources :

Example 6 with Random

use of gdsc.core.utils.Random in project GDSC-SMLM by aherbert.

the class PSFEstimator method calculateStatistics.

private boolean calculateStatistics(PeakFit fitter, double[] params, double[] params_dev) {
    debug("  Fitting PSF");
    swapStatistics();
    // Create the fit engine using the PeakFit plugin
    FitConfiguration fitConfig = config.getFitConfiguration();
    fitConfig.setInitialAngle((float) params[0]);
    fitConfig.setInitialPeakStdDev0((float) params[1]);
    fitConfig.setInitialPeakStdDev1((float) params[2]);
    ImageStack stack = imp.getImageStack();
    Rectangle roi = stack.getProcessor(1).getRoi();
    ImageSource source = new IJImageSource(imp);
    // Allow interlaced data by wrapping the image source
    if (interlacedData) {
        source = new InterlacedImageSource(source, dataStart, dataBlock, dataSkip);
    }
    // Allow frame aggregation by wrapping the image source
    if (integrateFrames > 1) {
        source = new AggregatedImageSource(source, integrateFrames);
    }
    fitter.initialiseImage(source, roi, true);
    fitter.addPeakResults(this);
    fitter.initialiseFitting();
    FitEngine engine = fitter.createFitEngine();
    // Use random slices
    int[] slices = new int[stack.getSize()];
    for (int i = 0; i < slices.length; i++) slices[i] = i + 1;
    Random rand = new Random();
    rand.shuffle(slices);
    IJ.showStatus("Fitting ...");
    // Use multi-threaded code for speed
    int i;
    for (i = 0; i < slices.length; i++) {
        int slice = slices[i];
        //debug("  Processing slice = %d\n", slice);
        IJ.showProgress(size(), settings.numberOfPeaks);
        ImageProcessor ip = stack.getProcessor(slice);
        // stack processor does not set the bounds required by ImageConverter
        ip.setRoi(roi);
        FitJob job = new FitJob(slice, ImageConverter.getData(ip), roi);
        engine.run(job);
        if (sampleSizeReached() || Utils.isInterrupted()) {
            break;
        }
    }
    if (Utils.isInterrupted()) {
        IJ.showProgress(1);
        engine.end(true);
        return false;
    }
    // Wait until we have enough results
    while (!sampleSizeReached() && !engine.isQueueEmpty()) {
        IJ.showProgress(size(), settings.numberOfPeaks);
        try {
            Thread.sleep(50);
        } catch (InterruptedException e) {
            break;
        }
    }
    // End now if we have enough samples
    engine.end(sampleSizeReached());
    IJ.showStatus("");
    IJ.showProgress(1);
    // This count will be an over-estimate given that the provider is ahead of the consumer
    // in this multi-threaded system
    debug("  Processed %d/%d slices (%d peaks)", i, slices.length, size());
    setParams(ANGLE, params, params_dev, sampleNew[ANGLE]);
    setParams(X, params, params_dev, sampleNew[X]);
    setParams(Y, params, params_dev, sampleNew[Y]);
    if (settings.showHistograms) {
        int[] idList = new int[NAMES.length];
        int count = 0;
        boolean requireRetile = false;
        for (int ii = 0; ii < 3; ii++) {
            if (sampleNew[ii].getN() == 0)
                continue;
            StoredDataStatistics stats = new StoredDataStatistics(sampleNew[ii].getValues());
            idList[count++] = Utils.showHistogram(TITLE, stats, NAMES[ii], 0, 0, settings.histogramBins, "Mean = " + Utils.rounded(stats.getMean()) + ". Median = " + Utils.rounded(sampleNew[ii].getPercentile(50)));
            requireRetile = requireRetile || Utils.isNewWindow();
        }
        if (requireRetile && count > 0) {
            new WindowOrganiser().tileWindows(Arrays.copyOf(idList, count));
        }
    }
    if (size() < 2) {
        log("ERROR: Insufficient number of fitted peaks, terminating ...");
        return false;
    }
    return true;
}
Also used : InterlacedImageSource(gdsc.smlm.results.InterlacedImageSource) AggregatedImageSource(gdsc.smlm.results.AggregatedImageSource) ImageStack(ij.ImageStack) Rectangle(java.awt.Rectangle) StoredDataStatistics(gdsc.core.utils.StoredDataStatistics) WindowOrganiser(ij.plugin.WindowOrganiser) IJImageSource(gdsc.smlm.ij.IJImageSource) ImageProcessor(ij.process.ImageProcessor) FitEngine(gdsc.smlm.engine.FitEngine) Random(gdsc.core.utils.Random) FitConfiguration(gdsc.smlm.fitting.FitConfiguration) InterlacedImageSource(gdsc.smlm.results.InterlacedImageSource) ImageSource(gdsc.smlm.results.ImageSource) AggregatedImageSource(gdsc.smlm.results.AggregatedImageSource) IJImageSource(gdsc.smlm.ij.IJImageSource) FitJob(gdsc.smlm.engine.FitJob)

Example 7 with Random

use of gdsc.core.utils.Random in project GDSC-SMLM by aherbert.

the class CreateData method run.

/*
	 * (non-Javadoc)
	 * 
	 * @see ij.plugin.PlugIn#run(java.lang.String)
	 */
public void run(String arg) {
    SMLMUsageTracker.recordPlugin(this.getClass(), arg);
    extraOptions = Utils.isExtraOptions();
    simpleMode = (arg != null && arg.contains("simple"));
    benchmarkMode = (arg != null && arg.contains("benchmark"));
    spotMode = (arg != null && arg.contains("spot"));
    trackMode = (arg != null && arg.contains("track"));
    if ("load".equals(arg)) {
        loadBenchmarkData();
        return;
    }
    // Each localisation is a simulated emission of light from a point in space and time
    List<LocalisationModel> localisations = null;
    // Each localisation set is a collection of localisations that represent all localisations
    // with the same ID that are on in the same image time frame (Note: the simulation
    // can create many localisations per fluorophore per time frame which is useful when 
    // modelling moving particles)
    List<LocalisationModelSet> localisationSets = null;
    // Each fluorophore contains the on and off times when light was emitted 
    List<? extends FluorophoreSequenceModel> fluorophores = null;
    if (simpleMode || benchmarkMode || spotMode) {
        if (!showSimpleDialog())
            return;
        resetMemory();
        // 1 second frames
        settings.exposureTime = 1000;
        areaInUm = settings.size * settings.pixelPitch * settings.size * settings.pixelPitch / 1e6;
        // Number of spots per frame
        int n = 0;
        int[] nextN = null;
        SpatialDistribution dist;
        if (benchmarkMode) {
            // --------------------
            // BENCHMARK SIMULATION
            // --------------------
            // Draw the same point on the image repeatedly
            n = 1;
            dist = createFixedDistribution();
            reportAndSaveFittingLimits(dist);
        } else if (spotMode) {
            // ---------------
            // SPOT SIMULATION
            // ---------------
            // The spot simulation draws 0 or 1 random point per frame. 
            // Ensure we have 50% of the frames with a spot.
            nextN = new int[settings.particles * 2];
            Arrays.fill(nextN, 0, settings.particles, 1);
            Random rand = new Random();
            rand.shuffle(nextN);
            // Only put spots in the central part of the image
            double border = settings.size / 4.0;
            dist = createUniformDistribution(border);
        } else {
            // -----------------
            // SIMPLE SIMULATION
            // -----------------
            // The simple simulation draws n random points per frame to achieve a specified density.
            // No points will appear in multiple frames.
            // Each point has a random number of photons sampled from a range.
            // We can optionally use a mask. Create his first as it updates the areaInUm
            dist = createDistribution();
            // Randomly sample (i.e. not uniform density in all frames)
            if (settings.samplePerFrame) {
                final double mean = areaInUm * settings.density;
                System.out.printf("Mean samples = %f\n", mean);
                if (mean < 0.5) {
                    GenericDialog gd = new GenericDialog(TITLE);
                    gd.addMessage("The mean samples per frame is low: " + Utils.rounded(mean) + "\n \nContinue?");
                    gd.enableYesNoCancel();
                    gd.hideCancelButton();
                    gd.showDialog();
                    if (!gd.wasOKed())
                        return;
                }
                PoissonDistribution poisson = new PoissonDistribution(createRandomGenerator(), mean, PoissonDistribution.DEFAULT_EPSILON, PoissonDistribution.DEFAULT_MAX_ITERATIONS);
                StoredDataStatistics samples = new StoredDataStatistics(settings.particles);
                while (samples.getSum() < settings.particles) {
                    samples.add(poisson.sample());
                }
                nextN = new int[samples.getN()];
                for (int i = 0; i < nextN.length; i++) nextN[i] = (int) samples.getValue(i);
            } else {
                // Use the density to get the number per frame
                n = (int) FastMath.max(1, Math.round(areaInUm * settings.density));
            }
        }
        RandomGenerator random = null;
        localisations = new ArrayList<LocalisationModel>(settings.particles);
        localisationSets = new ArrayList<LocalisationModelSet>(settings.particles);
        final int minPhotons = (int) settings.photonsPerSecond;
        final int range = (int) settings.photonsPerSecondMaximum - minPhotons + 1;
        if (range > 1)
            random = createRandomGenerator();
        // Add frames at the specified density until the number of particles has been reached
        int id = 0;
        int t = 0;
        while (id < settings.particles) {
            // Allow the number per frame to be specified
            if (nextN != null) {
                if (t >= nextN.length)
                    break;
                n = nextN[t];
            }
            // Simulate random positions in the frame for the specified density
            t++;
            for (int j = 0; j < n; j++) {
                final double[] xyz = dist.next();
                // Ignore within border. We do not want to draw things we cannot fit.
                //if (!distBorder.isWithinXY(xyz))
                //	continue;
                // Simulate random photons
                final int intensity = minPhotons + ((random != null) ? random.nextInt(range) : 0);
                LocalisationModel m = new LocalisationModel(id, t, xyz, intensity, LocalisationModel.CONTINUOUS);
                localisations.add(m);
                // Each localisation can be a separate localisation set
                LocalisationModelSet set = new LocalisationModelSet(id, t);
                set.add(m);
                localisationSets.add(set);
                id++;
            }
        }
    } else {
        if (!showDialog())
            return;
        resetMemory();
        areaInUm = settings.size * settings.pixelPitch * settings.size * settings.pixelPitch / 1e6;
        int totalSteps;
        double correlation = 0;
        ImageModel imageModel;
        if (trackMode) {
            // ----------------
            // TRACK SIMULATION
            // ----------------
            // In track mode we create fixed lifetime fluorophores that do not overlap in time.
            // This is the simplest simulation to test moving molecules.
            settings.seconds = (int) Math.ceil(settings.particles * (settings.exposureTime + settings.tOn) / 1000);
            totalSteps = 0;
            final double simulationStepsPerFrame = (settings.stepsPerSecond * settings.exposureTime) / 1000.0;
            imageModel = new FixedLifetimeImageModel(settings.stepsPerSecond * settings.tOn / 1000.0, simulationStepsPerFrame);
        } else {
            // ---------------
            // FULL SIMULATION
            // ---------------
            // The full simulation draws n random points in space.
            // The same molecule may appear in multiple frames, move and blink.
            //
            // Points are modelled as fluorophores that must be activated and then will 
            // blink and photo-bleach. The molecules may diffuse and this can be simulated 
            // with many steps per image frame. All steps from a frame are collected
            // into a localisation set which can be drawn on the output image.
            SpatialIllumination activationIllumination = createIllumination(settings.pulseRatio, settings.pulseInterval);
            // Generate additional frames so that each frame has the set number of simulation steps
            totalSteps = (int) Math.ceil(settings.seconds * settings.stepsPerSecond);
            // Since we have an exponential decay of activations
            // ensure half of the particles have activated by 30% of the frames.
            double eAct = totalSteps * 0.3 * activationIllumination.getAveragePhotons();
            // Q. Does tOn/tOff change depending on the illumination strength?
            imageModel = new ActivationEnergyImageModel(eAct, activationIllumination, settings.stepsPerSecond * settings.tOn / 1000.0, settings.stepsPerSecond * settings.tOffShort / 1000.0, settings.stepsPerSecond * settings.tOffLong / 1000.0, settings.nBlinksShort, settings.nBlinksLong);
            imageModel.setUseGeometricDistribution(settings.nBlinksGeometricDistribution);
            // Only use the correlation if selected for the distribution
            if (PHOTON_DISTRIBUTION[PHOTON_CORRELATED].equals(settings.photonDistribution))
                correlation = settings.correlation;
        }
        imageModel.setRandomGenerator(createRandomGenerator());
        imageModel.setPhotonBudgetPerFrame(true);
        imageModel.setDiffusion2D(settings.diffuse2D);
        imageModel.setRotation2D(settings.rotate2D);
        IJ.showStatus("Creating molecules ...");
        SpatialDistribution distribution = createDistribution();
        List<CompoundMoleculeModel> compounds = createCompoundMolecules();
        if (compounds == null)
            return;
        List<CompoundMoleculeModel> molecules = imageModel.createMolecules(compounds, settings.particles, distribution, settings.rotateInitialOrientation);
        // Activate fluorophores
        IJ.showStatus("Creating fluorophores ...");
        // Note: molecules list will be converted to compounds containing fluorophores
        fluorophores = imageModel.createFluorophores(molecules, totalSteps);
        if (fluorophores.isEmpty()) {
            IJ.error(TITLE, "No fluorophores created");
            return;
        }
        IJ.showStatus("Creating localisations ...");
        // TODO - Output a molecule Id for each fluorophore if using compound molecules. This allows analysis
        // of the ratio of trimers, dimers, monomers, etc that could be detected.
        totalSteps = checkTotalSteps(totalSteps, fluorophores);
        if (totalSteps == 0)
            return;
        imageModel.setPhotonDistribution(createPhotonDistribution());
        imageModel.setConfinementDistribution(createConfinementDistribution());
        // This should be optimised
        imageModel.setConfinementAttempts(10);
        localisations = imageModel.createImage(molecules, settings.fixedFraction, totalSteps, (double) settings.photonsPerSecond / settings.stepsPerSecond, correlation, settings.rotateDuringSimulation);
        // Re-adjust the fluorophores to the correct time
        if (settings.stepsPerSecond != 1) {
            final double scale = 1.0 / settings.stepsPerSecond;
            for (FluorophoreSequenceModel f : fluorophores) f.adjustTime(scale);
        }
        // Integrate the frames
        localisationSets = combineSimulationSteps(localisations);
        localisationSets = filterToImageBounds(localisationSets);
    }
    datasetNumber++;
    localisations = drawImage(localisationSets);
    if (localisations == null || localisations.isEmpty()) {
        IJ.error(TITLE, "No localisations created");
        return;
    }
    fluorophores = removeFilteredFluorophores(fluorophores, localisations);
    double signalPerFrame = showSummary(fluorophores, localisations);
    if (!benchmarkMode) {
        boolean fullSimulation = (!(simpleMode || spotMode));
        saveSimulationParameters(localisations.size(), fullSimulation, signalPerFrame);
    }
    IJ.showStatus("Saving data ...");
    //convertRelativeToAbsolute(molecules);
    saveFluorophores(fluorophores);
    saveImageResults(results);
    saveLocalisations(localisations);
    // The settings for the filenames may have changed
    SettingsManager.saveSettings(globalSettings);
    IJ.showStatus("Done");
}
Also used : PoissonDistribution(org.apache.commons.math3.distribution.PoissonDistribution) ActivationEnergyImageModel(gdsc.smlm.model.ActivationEnergyImageModel) CompoundMoleculeModel(gdsc.smlm.model.CompoundMoleculeModel) RandomGenerator(org.apache.commons.math3.random.RandomGenerator) Random(gdsc.core.utils.Random) GenericDialog(ij.gui.GenericDialog) SpatialIllumination(gdsc.smlm.model.SpatialIllumination) SpatialDistribution(gdsc.smlm.model.SpatialDistribution) StoredDataStatistics(gdsc.core.utils.StoredDataStatistics) LocalisationModel(gdsc.smlm.model.LocalisationModel) FluorophoreSequenceModel(gdsc.smlm.model.FluorophoreSequenceModel) FixedLifetimeImageModel(gdsc.smlm.model.FixedLifetimeImageModel) LocalisationModelSet(gdsc.smlm.model.LocalisationModelSet) ImageModel(gdsc.smlm.model.ImageModel) ActivationEnergyImageModel(gdsc.smlm.model.ActivationEnergyImageModel) FixedLifetimeImageModel(gdsc.smlm.model.FixedLifetimeImageModel)

Example 8 with Random

use of gdsc.core.utils.Random in project GDSC-SMLM by aherbert.

the class ImageSourceTest method memoryImageSourceCanReturnDataWithGet.

@Test
public void memoryImageSourceCanReturnDataWithGet() {
    int w = 5;
    int h = 3;
    float[][] data = createData(w, h, 15);
    MemoryImageSource source = new MemoryImageSource(w, h, data);
    int[] frames = new int[data.length];
    for (int i = 0; i < data.length; i++) frames[i] = i + 1;
    Random rand = new Random();
    rand.shuffle(frames);
    Assert.assertTrue(source.open());
    for (int i = 0; i < data.length; i++) {
        int frame = frames[i];
        Assert.assertTrue(source.isValid(frame));
        Assert.assertArrayEquals(data[frame - 1], source.get(frame), 0);
    }
    Assert.assertFalse(source.isValid(0));
    Assert.assertFalse(source.isValid(data.length + 1));
}
Also used : Random(gdsc.core.utils.Random) Test(org.junit.Test)

Example 9 with Random

use of gdsc.core.utils.Random in project GDSC-SMLM by aherbert.

the class ImageSourceTest method aggregatedImageSourceCanReturnDataWithGet.

@Test
public void aggregatedImageSourceCanReturnDataWithGet() {
    int w = 5;
    int h = 3;
    int aggregate = 3;
    float[][] data = createData(w, h, 15);
    ImageSource source = new AggregatedImageSource(new MemoryImageSource(w, h, data), aggregate);
    int[] frames = new int[data.length / 3];
    for (int i = 0, frame = 1; i < frames.length; i++, frame += 3) frames[i] = frame;
    Random rand = new Random();
    rand.shuffle(frames);
    Assert.assertTrue(source.open());
    for (int i = 0; i < frames.length; i++) {
        int frame = frames[i];
        Assert.assertTrue("Invalid frame " + frame, source.isValid(frame));
        float[] d = source.get(frame);
        Assert.assertEquals(frame, source.getStartFrameNumber());
        Assert.assertEquals(frame + 2, source.getEndFrameNumber());
        float[] all = combine(data[frame - 1], data[frame], data[frame + 1]);
        Assert.assertArrayEquals("Invalid frame data " + frame, all, d, 0);
    }
    Assert.assertFalse(source.isValid(0));
    Assert.assertFalse(source.isValid(data.length + 1));
}
Also used : Random(gdsc.core.utils.Random) Test(org.junit.Test)

Example 10 with Random

use of gdsc.core.utils.Random in project GDSC-SMLM by aherbert.

the class ImageSourceTest method aggregatedImageSourceCanReturnCroppedDataWithGet.

@Test
public void aggregatedImageSourceCanReturnCroppedDataWithGet() {
    int w = 5;
    int h = 3;
    int aggregate = 3;
    float[][] data = createData(w, h, 15);
    Rectangle bounds = new Rectangle(2, 1, 3, 1);
    ImageSource source = new AggregatedImageSource(new MemoryImageSource(w, h, data), aggregate);
    int[] frames = new int[data.length / 3];
    for (int i = 0, frame = 1; i < frames.length; i++, frame += 3) frames[i] = frame;
    Random rand = new Random();
    rand.shuffle(frames);
    Assert.assertTrue(source.open());
    for (int i = 0; i < frames.length; i++) {
        int frame = frames[i];
        Assert.assertTrue(source.isValid(frame));
        float[] d = source.get(frame, bounds);
        Assert.assertEquals(frame, source.getStartFrameNumber());
        Assert.assertEquals(frame + 2, source.getEndFrameNumber());
        float[] all = combine(crop(data[frame - 1], w, bounds), crop(data[frame], w, bounds), crop(data[frame + 1], w, bounds));
        Assert.assertArrayEquals("Invalid frame data " + frame, all, d, 0);
    }
    Assert.assertFalse(source.isValid(0));
    Assert.assertFalse(source.isValid(data.length + 1));
}
Also used : Random(gdsc.core.utils.Random) Rectangle(java.awt.Rectangle) Test(org.junit.Test)

Aggregations

Random (gdsc.core.utils.Random)13 Test (org.junit.Test)11 Rectangle (java.awt.Rectangle)4 StoredDataStatistics (gdsc.core.utils.StoredDataStatistics)2 FitEngine (gdsc.smlm.engine.FitEngine)1 FitJob (gdsc.smlm.engine.FitJob)1 FitConfiguration (gdsc.smlm.fitting.FitConfiguration)1 IJImageSource (gdsc.smlm.ij.IJImageSource)1 ActivationEnergyImageModel (gdsc.smlm.model.ActivationEnergyImageModel)1 CompoundMoleculeModel (gdsc.smlm.model.CompoundMoleculeModel)1 FixedLifetimeImageModel (gdsc.smlm.model.FixedLifetimeImageModel)1 FluorophoreSequenceModel (gdsc.smlm.model.FluorophoreSequenceModel)1 ImageModel (gdsc.smlm.model.ImageModel)1 LocalisationModel (gdsc.smlm.model.LocalisationModel)1 LocalisationModelSet (gdsc.smlm.model.LocalisationModelSet)1 SpatialDistribution (gdsc.smlm.model.SpatialDistribution)1 SpatialIllumination (gdsc.smlm.model.SpatialIllumination)1 AggregatedImageSource (gdsc.smlm.results.AggregatedImageSource)1 ImageSource (gdsc.smlm.results.ImageSource)1 InterlacedImageSource (gdsc.smlm.results.InterlacedImageSource)1