Search in sources :

Example 6 with LocalisationModel

use of uk.ac.sussex.gdsc.smlm.model.LocalisationModel in project GDSC-SMLM by aherbert.

the class BlinkEstimatorTest method estimateBlinking.

private TIntHashSet estimateBlinking(UniformRandomProvider rg, double blinkingRate, double ton, double toff, int particles, double fixedFraction, boolean timeAtLowerBound, boolean doAssert) {
    Assumptions.assumeTrue(TestSettings.allow(TestComplexity.MAXIMUM));
    final SpatialIllumination activationIllumination = new UniformIllumination(100);
    int totalSteps = 100;
    final double eAct = totalSteps * 0.3 * activationIllumination.getAveragePhotons();
    final ImageModel imageModel = new ActivationEnergyImageModel(eAct, activationIllumination, ton, 0, toff, 0, blinkingRate, rg);
    final double[] max = new double[] { 256, 256, 32 };
    final double[] min = new double[3];
    final SpatialDistribution distribution = new UniformDistribution(min, max, rg.nextInt());
    final List<CompoundMoleculeModel> compounds = new ArrayList<>(1);
    final CompoundMoleculeModel c = new CompoundMoleculeModel(1, 0, 0, 0, Arrays.asList(new MoleculeModel(0, 0, 0, 0)));
    c.setDiffusionRate(diffusionRate);
    c.setDiffusionType(DiffusionType.RANDOM_WALK);
    compounds.add(c);
    final List<CompoundMoleculeModel> molecules = imageModel.createMolecules(compounds, particles, distribution, false);
    // Activate fluorophores
    final List<? extends FluorophoreSequenceModel> fluorophores = imageModel.createFluorophores(molecules, totalSteps);
    totalSteps = checkTotalSteps(totalSteps, fluorophores);
    final List<LocalisationModel> localisations = imageModel.createImage(molecules, fixedFraction, totalSteps, photons, 0.5, false);
    // // Remove localisations to simulate missed counts.
    // List<LocalisationModel> newLocalisations = new
    // ArrayList<LocalisationModel>(localisations.size());
    // boolean[] id = new boolean[fluorophores.size() + 1];
    // Statistics photonStats = new Statistics();
    // for (LocalisationModel l : localisations)
    // {
    // photonStats.add(l.getIntensity());
    // // Remove by intensity threshold and optionally at random.
    // if (l.getIntensity() < minPhotons || rand.nextDouble() < pDelete)
    // continue;
    // newLocalisations.add(l);
    // id[l.getId()] = true;
    // }
    // localisations = newLocalisations;
    // logger.info("Photons = %f", photonStats.getMean());
    // 
    // List<FluorophoreSequenceModel> newFluorophores = new
    // ArrayList<FluorophoreSequenceModel>(fluorophores.size());
    // for (FluorophoreSequenceModel f : fluorophores)
    // {
    // if (id[f.getId()])
    // newFluorophores.add(f);
    // }
    // fluorophores = newFluorophores;
    final MemoryPeakResults results = new MemoryPeakResults();
    final CalibrationWriter calibration = new CalibrationWriter();
    calibration.setNmPerPixel(pixelPitch);
    calibration.setExposureTime(msPerFrame);
    calibration.setCountPerPhoton(1);
    results.setCalibration(calibration.getCalibration());
    results.setPsf(PsfHelper.create(PSFType.ONE_AXIS_GAUSSIAN_2D));
    final float b = 0;
    float intensity;
    final float z = 0;
    for (final LocalisationModel l : localisations) {
        // Remove by intensity threshold and optionally at random.
        if (l.getIntensity() < minPhotons || rg.nextDouble() < probabilityDelete) {
            continue;
        }
        final int frame = l.getTime();
        intensity = (float) l.getIntensity();
        final float x = (float) l.getX();
        final float y = (float) l.getY();
        final float[] params = Gaussian2DPeakResultHelper.createParams(b, intensity, x, y, z, psfWidth);
        results.add(frame, 0, 0, 0, 0, 0, 0, params, null);
    }
    // Add random localisations
    // Intensity doesn't matter at the moment for tracing
    intensity = (float) photons;
    for (int i = (int) (localisations.size() * probabilityAdd); i-- > 0; ) {
        final int frame = 1 + rg.nextInt(totalSteps);
        final float x = (float) (rg.nextDouble() * max[0]);
        final float y = (float) (rg.nextDouble() * max[1]);
        final float[] params = Gaussian2DPeakResultHelper.createParams(b, intensity, x, y, z, psfWidth);
        results.add(frame, 0, 0, 0, 0, 0, 0, params, null);
    }
    // Get actual simulated stats ...
    final Statistics statsNBlinks = new Statistics();
    final Statistics statsTOn = new Statistics();
    final Statistics statsTOff = new Statistics();
    final Statistics statsSampledNBlinks = new Statistics();
    final Statistics statsSampledTOn = new Statistics();
    final StoredDataStatistics statsSampledTOff = new StoredDataStatistics();
    for (final FluorophoreSequenceModel f : fluorophores) {
        statsNBlinks.add(f.getNumberOfBlinks());
        statsTOn.add(f.getOnTimes());
        statsTOff.add(f.getOffTimes());
        final int[] on = f.getSampledOnTimes();
        statsSampledNBlinks.add(on.length);
        statsSampledTOn.add(on);
        statsSampledTOff.add(f.getSampledOffTimes());
    }
    logger.info(FunctionUtils.getSupplier("N = %d (%d), N-blinks = %f, tOn = %f, tOff = %f, Fixed = %f", fluorophores.size(), localisations.size(), blinkingRate, ton, toff, fixedFraction));
    logger.info(FunctionUtils.getSupplier("Actual N-blinks = %f (%f), tOn = %f (%f), tOff = %f (%f), 95%% = %f, max = %f", statsNBlinks.getMean(), statsSampledNBlinks.getMean(), statsTOn.getMean(), statsSampledTOn.getMean(), statsTOff.getMean(), statsSampledTOff.getMean(), statsSampledTOff.getStatistics().getPercentile(95), statsSampledTOff.getStatistics().getMax()));
    logger.info("-=-=--=-");
    final BlinkEstimator be = new BlinkEstimator();
    be.setMaxDarkTime((int) (toff * 10));
    be.setMsPerFrame(msPerFrame);
    be.setRelativeDistance(false);
    final double d = ImageModel.getRandomMoveDistance(diffusionRate);
    be.setSearchDistance((fixedFraction < 1) ? Math.sqrt(2 * d * d) * 3 : 0);
    be.setTimeAtLowerBound(timeAtLowerBound);
    // Assertions.assertTrue("Max dark time must exceed the dark time of the data (otherwise no
    // plateau)",
    // be.maxDarkTime > statsSampledTOff.getStatistics().getMax());
    final int nMolecules = fluorophores.size();
    if (usePopulationStatistics) {
        blinkingRate = statsNBlinks.getMean();
        toff = statsTOff.getMean();
    } else {
        blinkingRate = statsSampledNBlinks.getMean();
        toff = statsSampledTOff.getMean();
    }
    // See if any fitting regime gets a correct answer
    final TIntHashSet ok = new TIntHashSet();
    for (int numberOfFittedPoints = MIN_FITTED_POINTS; numberOfFittedPoints <= MAX_FITTED_POINTS; numberOfFittedPoints++) {
        be.setNumberOfFittedPoints(numberOfFittedPoints);
        be.computeBlinkingRate(results, true);
        final double moleculesError = DoubleEquality.relativeError(nMolecules, be.getNMolecules());
        final double blinksError = DoubleEquality.relativeError(blinkingRate, be.getNBlinks());
        final double offError = DoubleEquality.relativeError(toff * msPerFrame, be.getTOff());
        logger.info(FunctionUtils.getSupplier("Error %d: N = %f, blinks = %f, tOff = %f : %f", numberOfFittedPoints, moleculesError, blinksError, offError, (moleculesError + blinksError + offError) / 3));
        if (moleculesError < relativeError && blinksError < relativeError && offError < relativeError) {
            ok.add(numberOfFittedPoints);
            logger.info("-=-=--=-");
            logger.info(FunctionUtils.getSupplier("*** Correct at %d fitted points ***", numberOfFittedPoints));
            if (doAssert) {
                break;
            }
        }
    // if (!be.isIncreaseNFittedPoints())
    // break;
    }
    logger.info("-=-=--=-");
    if (doAssert) {
        Assertions.assertFalse(ok.isEmpty());
    }
    // relativeError);
    return ok;
}
Also used : ActivationEnergyImageModel(uk.ac.sussex.gdsc.smlm.model.ActivationEnergyImageModel) CompoundMoleculeModel(uk.ac.sussex.gdsc.smlm.model.CompoundMoleculeModel) ArrayList(java.util.ArrayList) TIntHashSet(gnu.trove.set.hash.TIntHashSet) MoleculeModel(uk.ac.sussex.gdsc.smlm.model.MoleculeModel) CompoundMoleculeModel(uk.ac.sussex.gdsc.smlm.model.CompoundMoleculeModel) SpatialIllumination(uk.ac.sussex.gdsc.smlm.model.SpatialIllumination) CalibrationWriter(uk.ac.sussex.gdsc.smlm.data.config.CalibrationWriter) MemoryPeakResults(uk.ac.sussex.gdsc.smlm.results.MemoryPeakResults) SpatialDistribution(uk.ac.sussex.gdsc.smlm.model.SpatialDistribution) UniformDistribution(uk.ac.sussex.gdsc.smlm.model.UniformDistribution) StoredDataStatistics(uk.ac.sussex.gdsc.core.utils.StoredDataStatistics) StoredDataStatistics(uk.ac.sussex.gdsc.core.utils.StoredDataStatistics) Statistics(uk.ac.sussex.gdsc.core.utils.Statistics) UniformIllumination(uk.ac.sussex.gdsc.smlm.model.UniformIllumination) LocalisationModel(uk.ac.sussex.gdsc.smlm.model.LocalisationModel) FluorophoreSequenceModel(uk.ac.sussex.gdsc.smlm.model.FluorophoreSequenceModel) ActivationEnergyImageModel(uk.ac.sussex.gdsc.smlm.model.ActivationEnergyImageModel) ImageModel(uk.ac.sussex.gdsc.smlm.model.ImageModel)

Example 7 with LocalisationModel

use of uk.ac.sussex.gdsc.smlm.model.LocalisationModel in project GDSC-SMLM by aherbert.

the class CreateData method drawImage.

// StoredDataStatistics rawPhotons = new StoredDataStatistics();
// StoredDataStatistics drawPhotons = new StoredDataStatistics();
// private synchronized void addRaw(double d)
// {
// //rawPhotons.add(d);
// }
// 
// private synchronized void addDraw(double d)
// {
// //drawPhotons.add(d);
// }
/**
 * Create an image from the localisations using the configured PSF width. Draws a new stack image.
 *
 * <p>Note that the localisations are filtered using the signal. The input list of localisations
 * will be updated.
 *
 * @param localisationSets the localisation sets
 * @return The localisations
 */
private List<LocalisationModel> drawImage(final List<LocalisationModelSet> localisationSets) {
    if (localisationSets.isEmpty()) {
        return null;
    }
    // Create a new list for all localisation that are drawn (i.e. pass the signal filters)
    List<LocalisationModelSet> newLocalisations = Collections.synchronizedList(new ArrayList<LocalisationModelSet>(localisationSets.size()));
    photonsRemoved = new AtomicInteger();
    removedT1 = new AtomicInteger();
    removedTn = new AtomicInteger();
    photonStats = new SummaryStatistics();
    // Add drawn spots to memory
    results = new MemoryPeakResults();
    final CalibrationWriter c = new CalibrationWriter();
    c.setDistanceUnit(DistanceUnit.PIXEL);
    c.setIntensityUnit(IntensityUnit.PHOTON);
    c.setNmPerPixel(settings.getPixelPitch());
    c.setExposureTime(settings.getExposureTime());
    c.setCameraType(settings.getCameraType());
    if (settings.getCameraType() == CameraType.SCMOS) {
        c.setCameraModelName(settings.getCameraModelName());
    } else {
        final CreateDataSettingsHelper helper = new CreateDataSettingsHelper(settings);
        c.setCountPerPhoton(helper.getTotalGainSafe());
        c.setBias(settings.getBias());
        c.setReadNoise(helper.getReadNoiseInCounts());
        c.setQuantumEfficiency(helper.getQuantumEfficiency());
    }
    results.setCalibration(c.getCalibration());
    results.setSortAfterEnd(true);
    results.begin();
    maxT = localisationSets.get(localisationSets.size() - 1).getTime();
    // Display image
    final ImageStack stack = new ImageStack(settings.getSize(), settings.getSize(), maxT);
    final double psfSd = getPsfSd();
    if (psfSd <= 0) {
        return null;
    }
    final PsfModel psfModel = createPsfModel(localisationSets);
    if (psfModel == null) {
        return null;
    }
    // Create the camera noise model
    createPerPixelCameraModelData(cameraModel);
    createBackgroundPixels();
    final int threadCount = Prefs.getThreads();
    IJ.showStatus("Drawing image ...");
    // Multi-thread for speed
    final PeakResults syncResults = SynchronizedPeakResults.create(results, threadCount);
    final ExecutorService threadPool = Executors.newFixedThreadPool(threadCount);
    final List<Future<?>> futures = new LinkedList<>();
    // Count all the frames to process
    final Ticker ticker = ImageJUtils.createTicker(maxT, threadCount);
    // Collect statistics on the number of photons actually simulated
    // Process all frames
    int index = 0;
    int lastT = -1;
    for (final LocalisationModelSet l : localisationSets) {
        if (ImageJUtils.isInterrupted()) {
            break;
        }
        if (l.getTime() != lastT) {
            lastT = l.getTime();
            futures.add(threadPool.submit(new ImageGenerator(localisationSets, newLocalisations, index, lastT, psfModel.copy(), syncResults, stack, poissonNoise, createRandomGenerator(), ticker)));
        }
        index++;
    }
    // Finish processing data
    ConcurrencyUtils.waitForCompletionUnchecked(futures);
    futures.clear();
    if (ImageJUtils.isInterrupted()) {
        IJ.showProgress(1);
        return null;
    }
    // Do all the frames that had no localisations
    float[] limits = null;
    for (int t = 1; t <= maxT; t++) {
        if (ImageJUtils.isInterrupted()) {
            break;
        }
        final Object pixels = stack.getPixels(t);
        if (pixels == null) {
            futures.add(threadPool.submit(new ImageGenerator(localisationSets, newLocalisations, maxT, t, null, syncResults, stack, poissonNoise, createRandomGenerator(), ticker)));
        } else if (limits == null) {
            limits = MathUtils.limits((float[]) pixels);
        }
    }
    // Finish
    ConcurrencyUtils.waitForCompletionUnchecked(futures);
    futures.clear();
    threadPool.shutdown();
    IJ.showProgress(1);
    if (ImageJUtils.isInterrupted() || limits == null) {
        return null;
    }
    results.end();
    if (photonsRemoved.get() > 0) {
        ImageJUtils.log("Removed %d localisations with less than %.1f rendered photons", photonsRemoved.get(), settings.getMinPhotons());
    }
    if (removedT1.get() > 0) {
        ImageJUtils.log("Removed %d localisations with no neighbours @ SNR %.2f", removedT1.get(), settings.getMinSnrT1());
    }
    if (removedTn.get() > 0) {
        ImageJUtils.log("Removed %d localisations with valid neighbours @ SNR %.2f", removedTn.get(), settings.getMinSnrTN());
    }
    if (photonStats.getN() > 0) {
        ImageJUtils.log("Average photons rendered = %s +/- %s", MathUtils.rounded(photonStats.getMean()), MathUtils.rounded(photonStats.getStandardDeviation()));
    }
    // System.out.printf("rawPhotons = %f\n", rawPhotons.getMean());
    // System.out.printf("drawPhotons = %f\n", drawPhotons.getMean());
    // new HistogramPlotBuilder("draw photons", drawPhotons, "photons", true, 0, 1000);
    // Update with all those localisation that have been drawn
    localisationSets.clear();
    localisationSets.addAll(newLocalisations);
    newLocalisations = null;
    IJ.showStatus("Displaying image ...");
    ImageStack newStack = stack;
    if (!settings.getRawImage()) {
        // Get the global limits and ensure all values can be represented
        final Object[] imageArray = stack.getImageArray();
        limits = MathUtils.limits((float[]) imageArray[0]);
        for (int j = 1; j < imageArray.length; j++) {
            limits = MathUtils.limits(limits, (float[]) imageArray[j]);
        }
        // float limits0 = limits[0];
        // Leave bias in place
        final float limits0 = 0;
        // Check if the image will fit in a 16-bit range
        if ((limits[1] - limits0) < 65535) {
            // Convert to 16-bit
            newStack = new ImageStack(stack.getWidth(), stack.getHeight(), stack.getSize());
            // Account for rounding
            final float min = (float) (limits0 - 0.5);
            for (int j = 0; j < imageArray.length; j++) {
                final float[] image = (float[]) imageArray[j];
                final short[] pixels = new short[image.length];
                for (int k = 0; k < pixels.length; k++) {
                    pixels[k] = (short) (image[k] - min);
                }
                newStack.setPixels(pixels, j + 1);
                // Free memory
                imageArray[j] = null;
                // Attempt to stay within memory (check vs 32MB)
                if (MemoryUtils.getFreeMemory() < 33554432L) {
                    MemoryUtils.runGarbageCollectorOnce();
                }
            }
            for (int k = 2; k-- > 0; ) {
                limits[k] = (float) Math.floor(limits[k] - min);
            }
        } else {
            // Keep as 32-bit but round to whole numbers
            for (int j = 0; j < imageArray.length; j++) {
                final float[] pixels = (float[]) imageArray[j];
                for (int k = 0; k < pixels.length; k++) {
                    pixels[k] = Math.round(pixels[k]);
                }
            }
            for (int k = 2; k-- > 0; ) {
                limits[k] = Math.round(limits[k]);
            }
        }
    }
    // Show image
    final ImagePlus imp = ImageJUtils.display(CREATE_DATA_IMAGE_TITLE, newStack);
    final ij.measure.Calibration cal = new ij.measure.Calibration();
    String unit = "nm";
    double unitPerPixel = settings.getPixelPitch();
    if (unitPerPixel > 100) {
        unit = "um";
        unitPerPixel /= 1000.0;
    }
    cal.setUnit(unit);
    cal.pixelHeight = cal.pixelWidth = unitPerPixel;
    final Rectangle bounds = cameraModel.getBounds();
    if (bounds != null) {
        cal.xOrigin = -bounds.x;
        cal.yOrigin = -bounds.y;
    }
    imp.setCalibration(cal);
    imp.setDimensions(1, 1, newStack.getSize());
    imp.setDisplayRange(limits[0], limits[1]);
    imp.updateAndDraw();
    saveImage(imp);
    final IJImageSource imageSource = new IJImageSource(imp);
    // Shift simulation image source to correct location
    results.setSource(imageSource);
    results.setName(CREATE_DATA_IMAGE_TITLE + " (" + TITLE + ")");
    // Bounds are relative to the image source
    results.setBounds(new Rectangle(settings.getSize(), settings.getSize()));
    results.setPsf(createPsf(psfSd));
    MemoryPeakResults.addResults(results);
    setBenchmarkResults(imp, results);
    if (benchmarkMode && benchmarkParameters != null) {
        benchmarkParameters.setPhotons(results);
    }
    final List<LocalisationModel> localisations = toLocalisations(localisationSets);
    savePulses(localisations, results);
    // Saved the fixed and moving localisations into different datasets
    saveFixedAndMoving(results);
    saveCompoundMolecules(results);
    return localisations;
}
Also used : Rectangle(java.awt.Rectangle) AiryPsfModel(uk.ac.sussex.gdsc.smlm.model.AiryPsfModel) GaussianPsfModel(uk.ac.sussex.gdsc.smlm.model.GaussianPsfModel) ImagePsfModel(uk.ac.sussex.gdsc.smlm.model.ImagePsfModel) PsfModel(uk.ac.sussex.gdsc.smlm.model.PsfModel) MemoryPeakResults(uk.ac.sussex.gdsc.smlm.results.MemoryPeakResults) TextFilePeakResults(uk.ac.sussex.gdsc.smlm.results.TextFilePeakResults) ImmutableMemoryPeakResults(uk.ac.sussex.gdsc.smlm.results.ImmutableMemoryPeakResults) PeakResults(uk.ac.sussex.gdsc.smlm.results.PeakResults) SynchronizedPeakResults(uk.ac.sussex.gdsc.smlm.results.SynchronizedPeakResults) CalibrationWriter(uk.ac.sussex.gdsc.smlm.data.config.CalibrationWriter) MemoryPeakResults(uk.ac.sussex.gdsc.smlm.results.MemoryPeakResults) ImmutableMemoryPeakResults(uk.ac.sussex.gdsc.smlm.results.ImmutableMemoryPeakResults) ImageStack(ij.ImageStack) Ticker(uk.ac.sussex.gdsc.core.logging.Ticker) SummaryStatistics(org.apache.commons.math3.stat.descriptive.SummaryStatistics) ImagePlus(ij.ImagePlus) ReadHint(uk.ac.sussex.gdsc.smlm.results.ImageSource.ReadHint) LinkedList(java.util.LinkedList) IJImageSource(uk.ac.sussex.gdsc.smlm.ij.IJImageSource) LocalisationModel(uk.ac.sussex.gdsc.smlm.model.LocalisationModel) AtomicInteger(java.util.concurrent.atomic.AtomicInteger) CreateDataSettingsHelper(uk.ac.sussex.gdsc.smlm.ij.settings.CreateDataSettingsHelper) ExecutorService(java.util.concurrent.ExecutorService) Future(java.util.concurrent.Future) LocalisationModelSet(uk.ac.sussex.gdsc.smlm.model.LocalisationModelSet)

Example 8 with LocalisationModel

use of uk.ac.sussex.gdsc.smlm.model.LocalisationModel in project GDSC-SMLM by aherbert.

the class CreateData method removeFilteredFluorophores.

/**
 * Remove all fluorophores which were not drawn.
 *
 * @param fluorophores the fluorophores
 * @param localisations the localisations
 * @return the filtered list
 */
private List<? extends FluorophoreSequenceModel> removeFilteredFluorophores(List<? extends FluorophoreSequenceModel> fluorophores, List<LocalisationModel> localisations) {
    if (fluorophores == null) {
        return null;
    }
    // movingMolecules will be created with an initial capacity to hold all the unique IDs
    final TIntHashSet idSet = new TIntHashSet((movingMolecules != null) ? movingMolecules.capacity() : 0);
    for (final LocalisationModel l : localisations) {
        idSet.add(l.getId());
    }
    final List<FluorophoreSequenceModel> newFluorophores = new ArrayList<>(idSet.size());
    for (final FluorophoreSequenceModel f : fluorophores) {
        if (idSet.contains(f.getId())) {
            newFluorophores.add(f);
        }
    }
    return newFluorophores;
}
Also used : LocalisationModel(uk.ac.sussex.gdsc.smlm.model.LocalisationModel) FluorophoreSequenceModel(uk.ac.sussex.gdsc.smlm.model.FluorophoreSequenceModel) ArrayList(java.util.ArrayList) TIntArrayList(gnu.trove.list.array.TIntArrayList) TFloatArrayList(gnu.trove.list.array.TFloatArrayList) TIntHashSet(gnu.trove.set.hash.TIntHashSet)

Example 9 with LocalisationModel

use of uk.ac.sussex.gdsc.smlm.model.LocalisationModel in project GDSC-SMLM by aherbert.

the class CreateData method saveLocalisations.

/**
 * Save the localisations to a text file.
 *
 * @param localisations the localisations
 */
private void saveLocalisations(List<LocalisationModel> localisations) {
    if (!settings.getSaveLocalisations()) {
        return;
    }
    sortLocalisationsByTime(localisations);
    // Collections.sort(localisations, new Comparator<LocalisationModel>(){
    // 
    // public int compare(LocalisationModel o1, LocalisationModel o2)
    // {
    // int cellx1 = (int)(o1.getX() / settings.cellSize);
    // int cellx2 = (int)(o2.getX() / settings.cellSize);
    // int result = cellx2 - cellx1;
    // if (result != 0)
    // return result;
    // int celly1 = (int)(o1.getY() / settings.cellSize);
    // int celly2 = (int)(o2.getY() / settings.cellSize);
    // result = celly2 - celly1;
    // if (result != 0)
    // return result;
    // return (o1.getZ() == o2.getZ()) ? 0 : (o1.getZ() == 0) ? -1 : 1;
    // }});
    final LoadLocalisationsSettings.Builder settings = SettingsManager.readLoadLocalisationsSettings(0).toBuilder();
    final String[] path = ImageJUtils.decodePath(settings.getLocalisationsFilename());
    final OpenDialog chooser = new OpenDialog("Localisations_File", path[0], path[1]);
    if (chooser.getFileName() != null) {
        settings.setLocalisationsFilename(FileUtils.replaceExtension(chooser.getDirectory() + chooser.getFileName(), "xls"));
        SettingsManager.writeSettings(settings.build());
        try (BufferedWriter output = Files.newBufferedWriter(Paths.get(settings.getLocalisationsFilename()))) {
            output.write(createResultsFileHeader());
            output.write("#T\tId\tX\tY\tZ\tIntensity");
            output.newLine();
            for (final LocalisationModel l : localisations) {
                final StringBuilder sb = new StringBuilder();
                sb.append(l.getTime()).append('\t');
                sb.append(l.getId()).append('\t');
                sb.append(l.getX()).append('\t');
                sb.append(l.getY()).append('\t');
                sb.append(l.getZ()).append('\t');
                sb.append(l.getIntensity());
                output.write(sb.toString());
                output.newLine();
            }
        } catch (final Exception ex) {
            // Q. Add better handling of errors?
            ex.printStackTrace();
            IJ.log("Failed to save localisations to file: " + settings.getLocalisationsFilename());
        }
    }
}
Also used : LoadLocalisationsSettings(uk.ac.sussex.gdsc.smlm.ij.settings.GUIProtos.LoadLocalisationsSettings) LocalisationModel(uk.ac.sussex.gdsc.smlm.model.LocalisationModel) ConfigurationException(uk.ac.sussex.gdsc.smlm.data.config.ConfigurationException) IOException(java.io.IOException) DataException(uk.ac.sussex.gdsc.core.data.DataException) ConversionException(uk.ac.sussex.gdsc.core.data.utils.ConversionException) NullArgumentException(org.apache.commons.math3.exception.NullArgumentException) OpenDialog(ij.io.OpenDialog) BufferedWriter(java.io.BufferedWriter)

Example 10 with LocalisationModel

use of uk.ac.sussex.gdsc.smlm.model.LocalisationModel in project GDSC-SMLM by aherbert.

the class CreateData method combineSimulationSteps.

private List<LocalisationModelSet> combineSimulationSteps(List<LocalisationModel> localisations) {
    // Allow fractional integration steps
    final double simulationStepsPerFrame = (settings.getStepsPerSecond() * settings.getExposureTime()) / 1000.0;
    final List<LocalisationModelSet> newLocalisations = new ArrayList<>((int) (localisations.size() / simulationStepsPerFrame));
    // System.out.printf("combineSimulationSteps @ %f\n", simulationStepsPerFrame);
    // final double gain = new CreateDataSettingsHelper(settings).getTotalGainSafe();
    sortLocalisationsByIdThenTime(localisations);
    final int[] idList = getIds(localisations);
    movingMolecules = new TIntHashSet(idList.length);
    int index = 0;
    for (final int id : idList) {
        final int fromIndex = findIndexById(localisations, index, id);
        if (fromIndex > -1) {
            final int toIndex = findLastIndexById(localisations, fromIndex, id);
            final List<LocalisationModel> subset = localisations.subList(fromIndex, toIndex + 1);
            index = toIndex;
            // Store the IDs of any moving molecules
            if (isMoving(subset)) {
                movingMolecules.add(id);
            }
            // The frames may be longer or shorter than the simulation steps. Allocate the step
            // proportionately to each frame it overlaps:
            // 
            // Steps: |-- 0 --|-- 1 --|-- 2 --|--
            // Frames: |--- 0 ---|--- 1 ---|--- 2 ---|
            // 
            // ^ ^
            // | |
            // | End frame
            // |
            // Start frame
            final double firstFrame = getStartFrame(subset.get(0), simulationStepsPerFrame);
            final double lastFrame = getEndFrame(subset.get(subset.size() - 1), simulationStepsPerFrame);
            // Get the first frame offset and allocate space to store all potential frames
            final int intFirstFrame = (int) firstFrame;
            final int intLastFrame = (int) Math.ceil(lastFrame);
            final LocalisationModelSet[] sets = new LocalisationModelSet[intLastFrame - intFirstFrame + 1];
            // Process each step
            for (final LocalisationModel l : subset) {
                // Get the fractional start and end frames
                final double startFrame = getStartFrame(l, simulationStepsPerFrame);
                final double endFrame = getEndFrame(l, simulationStepsPerFrame);
                // Round down to get the actual frames that are overlapped
                final int start = (int) startFrame;
                int end = (int) endFrame;
                // that frame
                if (end > start && endFrame == end) {
                    // E.g. convert
                    // Steps: |-- 0 --|
                    // Frames: |- 0 -|- 1 -|- 2 -|
                    // to
                    // Steps: |-- 0 --|
                    // Frames: |- 0 -|- 1 -|
                    end--;
                }
                if (start == end) {
                    // If the step falls within one frame then add it to the set
                    final int tIndex = start - intFirstFrame;
                    if (sets[tIndex] == null) {
                        sets[tIndex] = new LocalisationModelSet(id, start);
                    }
                    sets[tIndex].add(l);
                } else {
                    // Add the localisation to all the frames that the step spans
                    final double total = endFrame - startFrame;
                    // double t = 0;
                    for (int frame = start; frame <= end; frame++) {
                        // Get the fraction to allocate to this frame
                        double fraction;
                        int state = (l.isContinuous() ? LocalisationModel.CONTINUOUS : 0);
                        if (frame == start) {
                            state |= LocalisationModel.NEXT | (l.hasPrevious() ? LocalisationModel.PREVIOUS : 0);
                            // start
                            if (startFrame == start) {
                                fraction = 1;
                            } else {
                                fraction = (Math.ceil(startFrame) - startFrame);
                            }
                        } else if (frame == end) {
                            state |= LocalisationModel.PREVIOUS | (l.hasNext() ? LocalisationModel.NEXT : 0);
                            // |=====|----|
                            // | |
                            // | endFrame
                            // end
                            fraction = (endFrame - end);
                        } else {
                            state |= LocalisationModel.CONTINUOUS;
                            fraction = 1;
                        }
                        // t += fraction;
                        // Add to the set
                        final int tIndex = frame - intFirstFrame;
                        if (sets[tIndex] == null) {
                            sets[tIndex] = new LocalisationModelSet(id, frame);
                        }
                        sets[tIndex].add(getFraction(l, fraction / total, state));
                    }
                // if (t < total * 0.98)
                // {
                // System.out.printf("Total error %g < %g : %f (%d) -> %f (%d)\n", t, total, startFrame,
                // start, endFrame, end);
                // }
                }
            }
            LocalisationModelSet previous = null;
            for (int i = 0; i < sets.length; i++) {
                if (sets[i] != null) {
                    sets[i].setPrevious(previous);
                    // Create a data array and store the current intensity.
                    // This is used later to filter based on SNR
                    sets[i].setData(new double[] { // * gain
                    0, // * gain
                    0, // * gain
                    0, // * gain
                    0, // * gain
                    sets[i].getIntensity() });
                    newLocalisations.add(sets[i]);
                }
                previous = sets[i];
            }
        }
    }
    // Sort by time
    Collections.sort(newLocalisations, (r1, r2) -> Integer.compare(r1.getTime(), r2.getTime()));
    return newLocalisations;
}
Also used : LocalisationModel(uk.ac.sussex.gdsc.smlm.model.LocalisationModel) ArrayList(java.util.ArrayList) TIntArrayList(gnu.trove.list.array.TIntArrayList) TFloatArrayList(gnu.trove.list.array.TFloatArrayList) LocalisationModelSet(uk.ac.sussex.gdsc.smlm.model.LocalisationModelSet) TIntHashSet(gnu.trove.set.hash.TIntHashSet) ReadHint(uk.ac.sussex.gdsc.smlm.results.ImageSource.ReadHint)

Aggregations

LocalisationModel (uk.ac.sussex.gdsc.smlm.model.LocalisationModel)12 ReadHint (uk.ac.sussex.gdsc.smlm.results.ImageSource.ReadHint)8 LocalisationModelSet (uk.ac.sussex.gdsc.smlm.model.LocalisationModelSet)6 TIntArrayList (gnu.trove.list.array.TIntArrayList)4 TIntHashSet (gnu.trove.set.hash.TIntHashSet)4 IOException (java.io.IOException)4 ArrayList (java.util.ArrayList)4 NullArgumentException (org.apache.commons.math3.exception.NullArgumentException)4 DataException (uk.ac.sussex.gdsc.core.data.DataException)4 ConversionException (uk.ac.sussex.gdsc.core.data.utils.ConversionException)4 ConfigurationException (uk.ac.sussex.gdsc.smlm.data.config.ConfigurationException)4 FluorophoreSequenceModel (uk.ac.sussex.gdsc.smlm.model.FluorophoreSequenceModel)4 MemoryPeakResults (uk.ac.sussex.gdsc.smlm.results.MemoryPeakResults)4 TFloatArrayList (gnu.trove.list.array.TFloatArrayList)3 ImagePlus (ij.ImagePlus)3 StoredDataStatistics (uk.ac.sussex.gdsc.core.utils.StoredDataStatistics)3 CalibrationWriter (uk.ac.sussex.gdsc.smlm.data.config.CalibrationWriter)3 ActivationEnergyImageModel (uk.ac.sussex.gdsc.smlm.model.ActivationEnergyImageModel)3 CompoundMoleculeModel (uk.ac.sussex.gdsc.smlm.model.CompoundMoleculeModel)3 ImageModel (uk.ac.sussex.gdsc.smlm.model.ImageModel)3