Examples with PSFOffset gdsc.smlm.ij.settings.PSFOffset used on opensource projects

Example 1 with PSFOffset

use of gdsc.smlm.ij.settings.PSFOffset in project GDSC-SMLM by aherbert.

the class PSFDrift method createImagePSF.

private ImagePSFModel createImagePSF(int lower, int upper) {
    int zCentre = psfSettings.zCentre;
    final double unitsPerPixel = 1.0 / scale;
    // So we can move from -depth to depth
    final double unitsPerSlice = 1;
    // Extract data uses index not slice number as arguments so subtract 1
    double noiseFraction = 1e-3;
    ImagePSFModel model = new ImagePSFModel(CreateData.extractImageStack(imp, lower - 1, upper - 1), zCentre - lower, unitsPerPixel, unitsPerSlice, psfSettings.fwhm, noiseFraction);
    // Add the calibrated centres
    if (psfSettings.offset != null && useOffset) {
        int sliceOffset = lower;
        for (PSFOffset offset : psfSettings.offset) {
            model.setRelativeCentre(offset.slice - sliceOffset, offset.cx, offset.cy);
        }
    }
    return model;
}

Example 2 with PSFOffset

use of gdsc.smlm.ij.settings.PSFOffset in project GDSC-SMLM by aherbert.

the class PSFDrift method computeDrift.

private void computeDrift() {
    // Create a grid of XY offset positions between 0-1 for PSF insert
    final double[] grid = new double[gridSize];
    for (int i = 0; i < grid.length; i++) grid[i] = (double) i / gridSize;
    // Configure fitting region
    final int w = 2 * regionSize + 1;
    centrePixel = w / 2;
    // Check region size using the image PSF
    double newPsfWidth = (double) imp.getWidth() / scale;
    if (Math.ceil(newPsfWidth) > w)
        Utils.log(TITLE + ": Fitted region size (%d) is smaller than the scaled PSF (%.1f)", w, newPsfWidth);
    // Create robust PSF fitting settings
    final double a = psfSettings.nmPerPixel * scale;
    final double sa = PSFCalculator.squarePixelAdjustment(psfSettings.nmPerPixel * (psfSettings.fwhm / Gaussian2DFunction.SD_TO_FWHM_FACTOR), a);
    fitConfig.setInitialPeakStdDev(sa / a);
    fitConfig.setBackgroundFitting(backgroundFitting);
    fitConfig.setNotSignalFitting(false);
    fitConfig.setComputeDeviations(false);
    fitConfig.setDisableSimpleFilter(true);
    // Create the PSF over the desired z-depth
    int depth = (int) Math.round(zDepth / psfSettings.nmPerSlice);
    int startSlice = psfSettings.zCentre - depth;
    int endSlice = psfSettings.zCentre + depth;
    int nSlices = imp.getStackSize();
    startSlice = (startSlice < 1) ? 1 : (startSlice > nSlices) ? nSlices : startSlice;
    endSlice = (endSlice < 1) ? 1 : (endSlice > nSlices) ? nSlices : endSlice;
    ImagePSFModel psf = createImagePSF(startSlice, endSlice);
    int minz = startSlice - psfSettings.zCentre;
    int maxz = endSlice - psfSettings.zCentre;
    final int nZ = maxz - minz + 1;
    final int gridSize2 = grid.length * grid.length;
    total = nZ * gridSize2;
    // Store all the fitting results
    int nStartPoints = getNumberOfStartPoints();
    results = new double[total * nStartPoints][];
    // TODO - Add ability to iterate this, adjusting the current offset in the PSF
    // each iteration
    // Create a pool of workers
    int nThreads = Prefs.getThreads();
    BlockingQueue<Job> jobs = new ArrayBlockingQueue<Job>(nThreads * 2);
    List<Worker> workers = new LinkedList<Worker>();
    List<Thread> threads = new LinkedList<Thread>();
    for (int i = 0; i < nThreads; i++) {
        Worker worker = new Worker(jobs, psf, w, fitConfig);
        Thread t = new Thread(worker);
        workers.add(worker);
        threads.add(t);
        t.start();
    }
    // Fit 
    Utils.showStatus("Fitting ...");
    final int step = Utils.getProgressInterval(total);
    outer: for (int z = minz, i = 0; z <= maxz; z++) {
        for (int x = 0; x < grid.length; x++) for (int y = 0; y < grid.length; y++, i++) {
            if (IJ.escapePressed()) {
                break outer;
            }
            put(jobs, new Job(z, grid[x], grid[y], i));
            if (i % step == 0) {
                IJ.showProgress(i, total);
            }
        }
    }
    // If escaped pressed then do not need to stop the workers, just return
    if (Utils.isInterrupted()) {
        IJ.showProgress(1);
        return;
    }
    // Finish all the worker threads by passing in a null job
    for (int i = 0; i < threads.size(); i++) {
        put(jobs, new Job());
    }
    // Wait for all to finish
    for (int i = 0; i < threads.size(); i++) {
        try {
            threads.get(i).join();
        } catch (InterruptedException e) {
            e.printStackTrace();
        }
    }
    threads.clear();
    IJ.showProgress(1);
    IJ.showStatus("");
    // Plot the average and SE for the drift curve
    // Plot the recall
    double[] zPosition = new double[nZ];
    double[] avX = new double[nZ];
    double[] seX = new double[nZ];
    double[] avY = new double[nZ];
    double[] seY = new double[nZ];
    double[] recall = new double[nZ];
    for (int z = minz, i = 0; z <= maxz; z++, i++) {
        Statistics statsX = new Statistics();
        Statistics statsY = new Statistics();
        for (int s = 0; s < nStartPoints; s++) {
            int resultPosition = i * gridSize2 + s * total;
            final int endResultPosition = resultPosition + gridSize2;
            while (resultPosition < endResultPosition) {
                if (results[resultPosition] != null) {
                    statsX.add(results[resultPosition][0]);
                    statsY.add(results[resultPosition][1]);
                }
                resultPosition++;
            }
        }
        zPosition[i] = z * psfSettings.nmPerSlice;
        avX[i] = statsX.getMean();
        seX[i] = statsX.getStandardError();
        avY[i] = statsY.getMean();
        seY[i] = statsY.getStandardError();
        recall[i] = (double) statsX.getN() / (nStartPoints * gridSize2);
    }
    // Find the range from the z-centre above the recall limit 
    int centre = 0;
    for (int slice = startSlice, i = 0; slice <= endSlice; slice++, i++) {
        if (slice == psfSettings.zCentre) {
            centre = i;
            break;
        }
    }
    if (recall[centre] < recallLimit)
        return;
    int start = centre, end = centre;
    for (int i = centre; i-- > 0; ) {
        if (recall[i] < recallLimit)
            break;
        start = i;
    }
    for (int i = centre; ++i < recall.length; ) {
        if (recall[i] < recallLimit)
            break;
        end = i;
    }
    int iterations = 1;
    LoessInterpolator loess = null;
    if (smoothing > 0)
        loess = new LoessInterpolator(smoothing, iterations);
    double[][] smoothx = displayPlot("Drift X", "X (nm)", zPosition, avX, seX, loess, start, end);
    double[][] smoothy = displayPlot("Drift Y", "Y (nm)", zPosition, avY, seY, loess, start, end);
    displayPlot("Recall", "Recall", zPosition, recall, null, null, start, end);
    WindowOrganiser wo = new WindowOrganiser();
    wo.tileWindows(idList);
    // Ask the user if they would like to store them in the image
    GenericDialog gd = new GenericDialog(TITLE);
    gd.enableYesNoCancel();
    gd.hideCancelButton();
    startSlice = psfSettings.zCentre - (centre - start);
    endSlice = psfSettings.zCentre + (end - centre);
    gd.addMessage(String.format("Save the drift to the PSF?\n \nSlices %d (%s nm) - %d (%s nm) above recall limit", startSlice, Utils.rounded(zPosition[start]), endSlice, Utils.rounded(zPosition[end])));
    gd.addMessage("Optionally average the end points to set drift outside the limits.\n(Select zero to ignore)");
    gd.addSlider("Number_of_points", 0, 10, positionsToAverage);
    gd.showDialog();
    if (gd.wasOKed()) {
        positionsToAverage = Math.abs((int) gd.getNextNumber());
        ArrayList<PSFOffset> offset = new ArrayList<PSFOffset>();
        final double pitch = psfSettings.nmPerPixel;
        int j = 0, jj = 0;
        for (int i = start, slice = startSlice; i <= end; slice++, i++) {
            j = findCentre(zPosition[i], smoothx, j);
            if (j == -1) {
                Utils.log("Failed to find the offset for depth %.2f", zPosition[i]);
                continue;
            }
            // The offset should store the difference to the centre in pixels so divide by the pixel pitch
            double cx = smoothx[1][j] / pitch;
            double cy = smoothy[1][j] / pitch;
            jj = findOffset(slice, jj);
            if (jj != -1) {
                cx += psfSettings.offset[jj].cx;
                cy += psfSettings.offset[jj].cy;
            }
            offset.add(new PSFOffset(slice, cx, cy));
        }
        addMissingOffsets(startSlice, endSlice, nSlices, offset);
        psfSettings.offset = offset.toArray(new PSFOffset[offset.size()]);
        psfSettings.addNote(TITLE, String.format("Solver=%s, Region=%d", PeakFit.getSolverName(fitConfig), regionSize));
        imp.setProperty("Info", XmlUtils.toXML(psfSettings));
    }
}

Example 3 with PSFOffset

use of gdsc.smlm.ij.settings.PSFOffset in project GDSC-SMLM by aherbert.

the class CreateData method createImagePSF.

/**
	 * Create a PSF model from the image that contains all the z-slices needed to draw the given localisations
	 * 
	 * @param localisationSets
	 * @return
	 */
private ImagePSFModel createImagePSF(List<LocalisationModelSet> localisationSets) {
    ImagePlus imp = WindowManager.getImage(settings.psfImageName);
    if (imp == null) {
        IJ.error(TITLE, "Unable to create the PSF model from image: " + settings.psfImageName);
        return null;
    }
    try {
        Object o = XmlUtils.fromXML(imp.getProperty("Info").toString());
        if (!(o != null && o instanceof PSFSettings))
            throw new RuntimeException("Unknown PSF settings for image: " + imp.getTitle());
        PSFSettings psfSettings = (PSFSettings) o;
        // Check all the settings have values
        if (psfSettings.nmPerPixel <= 0)
            throw new RuntimeException("Missing nmPerPixel calibration settings for image: " + imp.getTitle());
        if (psfSettings.nmPerSlice <= 0)
            throw new RuntimeException("Missing nmPerSlice calibration settings for image: " + imp.getTitle());
        if (psfSettings.zCentre <= 0)
            throw new RuntimeException("Missing zCentre calibration settings for image: " + imp.getTitle());
        if (psfSettings.fwhm <= 0)
            throw new RuntimeException("Missing FWHM calibration settings for image: " + imp.getTitle());
        // To save memory construct the Image PSF using only the slices that are within 
        // the depth of field of the simulation
        double minZ = Double.POSITIVE_INFINITY, maxZ = Double.NEGATIVE_INFINITY;
        for (LocalisationModelSet l : localisationSets) {
            for (LocalisationModel m : l.getLocalisations()) {
                final double z = m.getZ();
                if (minZ > z)
                    minZ = z;
                if (maxZ < z)
                    maxZ = z;
            }
        }
        int nSlices = imp.getStackSize();
        // z-centre should be an index and not the ImageJ slice number so subtract 1
        int zCentre = psfSettings.zCentre - 1;
        // Calculate the start/end slices to cover the depth of field
        // This logic must match the ImagePSFModel.
        final double unitsPerSlice = psfSettings.nmPerSlice / settings.pixelPitch;
        // We assume the PSF was imaged axially with increasing z-stage position (moving the stage 
        // closer to the objective). Thus higher z-coordinate are for higher slice numbers.
        int lower = (int) Math.round(minZ / unitsPerSlice) + zCentre;
        int upper = (int) Math.round(maxZ / unitsPerSlice) + zCentre;
        upper = (upper < 0) ? 0 : (upper >= nSlices) ? nSlices - 1 : upper;
        lower = (lower < 0) ? 0 : (lower >= nSlices) ? nSlices - 1 : lower;
        // Image PSF requires the z-centre for normalisation
        if (!(lower <= zCentre && upper >= zCentre)) {
            // Ensure we include the zCentre
            lower = Math.min(lower, zCentre);
            upper = Math.max(upper, zCentre);
        }
        final double noiseFraction = 1e-3;
        final ImagePSFModel model = new ImagePSFModel(extractImageStack(imp, lower, upper), zCentre - lower, psfSettings.nmPerPixel / settings.pixelPitch, unitsPerSlice, psfSettings.fwhm, noiseFraction);
        // Add the calibrated centres
        if (psfSettings.offset != null) {
            final int sliceOffset = lower + 1;
            for (PSFOffset offset : psfSettings.offset) {
                model.setRelativeCentre(offset.slice - sliceOffset, offset.cx, offset.cy);
            }
        }
        // Initialise the HWHM table so that it can be cloned
        model.initialiseHWHM();
        return model;
    } catch (Exception e) {
        IJ.error(TITLE, "Unable to create the image PSF model:\n" + e.getMessage());
        return null;
    }
}