Examples with Spot gdsc.smlm.tsf.TaggedSpotFile.Spot used on opensource projects

Example 1 with Spot

use of gdsc.smlm.tsf.TaggedSpotFile.Spot in project GDSC-SMLM by aherbert.

the class ResultsManagerTest method writeTSFMatchesRead.

private void writeTSFMatchesRead(int channels, int slices, int positions, int types) {
    String filename = createFile();
    FileOutputStream out = null;
    try {
        out = new FileOutputStream(filename);
    } catch (Exception e) {
        closeOutput(out);
        e.printStackTrace();
        Assert.fail(e.getMessage());
    }
    // Write the offsets used in the TSF format
    try {
        DataOutputStream dos = new DataOutputStream(out);
        dos.writeInt(0);
        dos.writeLong(0);
    } catch (IOException e) {
        closeOutput(out);
        e.printStackTrace();
        Assert.fail(e.getMessage());
    }
    // Generate random spots
    int size = 1000;
    Spot[] spots = new Spot[size];
    for (int i = 1; i <= size; i++) {
        Spot.Builder builder = Spot.newBuilder();
        builder.setChannel(1 + rand.nextInt(channels));
        builder.setSlice(1 + rand.nextInt(slices));
        builder.setPos(1 + rand.nextInt(positions));
        builder.setFluorophoreType(rand.nextInt(1, types));
        // This is a required field but is ignored when reading
        builder.setMolecule(i);
        builder.setCluster(rand.nextInt(10));
        builder.setFrame(rand.nextInt(1, 100));
        builder.setXPosition(rand.nextInt(50));
        builder.setYPosition(rand.nextInt(50));
        builder.setBackground(rand.next());
        builder.setIntensity(rand.next());
        builder.setX(rand.next());
        builder.setY(rand.next());
        builder.setWidth(TSFPeakResultsWriter.SD_TO_FWHM_FACTOR * rand.next());
        Spot spot = builder.build();
        spots[i - 1] = spot;
        try {
            spot.writeDelimitedTo(out);
        } catch (IOException e) {
            closeOutput(out);
            e.printStackTrace();
            Assert.fail(e.getMessage());
        }
    }
    // Write the header
    // Get the offset to the SpotList message
    long offset = 0;
    try {
        // The offset is the amount to skip forward after reading the int 
        // magic number (4 bytes) and long offset (8 bytes)
        //out.flush();
        offset = out.getChannel().position() - 12;
    } catch (IOException e) {
        closeOutput(out);
        e.printStackTrace();
        Assert.fail(e.getMessage());
    }
    // Record the SpotList message
    SpotList.Builder builder = SpotList.newBuilder();
    builder.setApplicationId(1);
    builder.setNrSpots(size);
    builder.setLocationUnits(LocationUnits.PIXELS);
    builder.setIntensityUnits(IntensityUnits.COUNTS);
    builder.setFitMode(FitMode.ONEAXIS);
    builder.setNrChannels(channels);
    builder.setNrSlices(slices);
    builder.setNrPos(positions);
    for (int type = 1; type <= types; type++) {
        FluorophoreType.Builder typeBuilder = FluorophoreType.newBuilder();
        typeBuilder.setId(type);
        typeBuilder.setDescription("Type " + type);
        typeBuilder.setIsFiducial(rand.next() < 0.5f);
        builder.addFluorophoreTypes(typeBuilder.build());
    }
    SpotList spotList = builder.build();
    try {
        spotList.writeDelimitedTo(out);
    } catch (IOException e) {
        e.printStackTrace();
        Assert.fail(e.getMessage());
    } finally {
        closeOutput(out);
    }
    // Write the offset to the SpotList message into the offset position
    RandomAccessFile f = null;
    try {
        f = new RandomAccessFile(new File(filename), "rw");
        f.seek(4);
        f.writeLong(offset);
    } catch (Exception e) {
        e.printStackTrace();
        Assert.fail(e.getMessage());
    } finally {
        if (f != null) {
            try {
                f.close();
            } catch (IOException e) {
            }
        }
    }
    // Read each combination 
    for (int channel = 1; channel <= channels; channel++) for (int slice = 1; slice <= slices; slice++) for (int position = 1; position <= positions; position++) for (int type = 1; type <= types; type++) {
        StringBuilder sb = new StringBuilder();
        sb.append(" channel=").append(channel);
        sb.append(" slice=").append(slice);
        sb.append(" position=").append(position);
        sb.append(" fluorophore_type=[").append(type).append(":Type ").append(type).append(":fiducial=").append(builder.getFluorophoreTypes(type - 1).getIsFiducial()).append(']');
        // This is needed to trick the Macro class into returning the options 
        // for the thread to the GenericDialog used in the ResultsManager
        Thread.currentThread().setName("Run$_");
        Macro.setOptions(sb.toString());
        ResultsManager.setInputFilename(filename);
        MemoryPeakResults in = ResultsManager.loadInputResults(ResultsManager.INPUT_FILE, false);
        checkEqual(spots, channel, slice, position, type, in);
    }
}

Example 2 with Spot

use of gdsc.smlm.tsf.TaggedSpotFile.Spot in project GDSC-SMLM by aherbert.

the class TSFPeakResultsWriter method addAll.

public void addAll(Collection<PeakResult> results) {
    if (out == null)
        return;
    Spot[] spots = new Spot[20];
    int count = 0;
    Spot.Builder builder = Spot.newBuilder();
    for (PeakResult result : results) {
        final float[] params = result.params;
        builder.setMolecule(id.incrementAndGet());
        builder.setChannel(1);
        builder.setFluorophoreType(1);
        builder.setFrame(result.getFrame());
        builder.setXPosition(result.origX);
        builder.setYPosition(result.origY);
        setBackground(builder, params[Gaussian2DFunction.BACKGROUND]);
        builder.setIntensity(params[Gaussian2DFunction.SIGNAL]);
        builder.setX(params[Gaussian2DFunction.X_POSITION]);
        builder.setY(params[Gaussian2DFunction.Y_POSITION]);
        setWidth(params, builder);
        if (result.hasPrecision()) {
            // Use the actual precision
            float precision = (float) result.getPrecision();
            builder.setXPrecision(precision);
            builder.setYPrecision(precision);
        } else if (canComputePrecision) {
            // Compute precision
            double s = (params[Gaussian2DFunction.X_SD] + params[Gaussian2DFunction.Y_SD]) * 0.5 * nmPerPixel;
            float precision = (float) PeakResult.getPrecision(nmPerPixel, s, params[Gaussian2DFunction.SIGNAL] / gain, result.noise / gain, isEmCCD);
            builder.setXPrecision(precision);
            builder.setYPrecision(precision);
        }
        if (result.hasId())
            builder.setCluster(result.getId());
        else
            builder.clearCluster();
        builder.setError(result.error);
        builder.setNoise(result.noise);
        if (result.hasEndFrame())
            builder.setEndFrame(result.getEndFrame());
        else
            builder.clearEndFrame();
        builder.setOriginalValue(result.origValue);
        addParamsStdDev(builder, result.paramsStdDev);
        spots[count++] = builder.build();
        // Flush the output to allow for very large input lists
        if (count >= spots.length) {
            writeResult(count, spots);
            if (!isActive())
                return;
            count = 0;
        }
    }
    writeResult(count, spots);
}

Example 3 with Spot

use of gdsc.smlm.tsf.TaggedSpotFile.Spot in project GDSC-SMLM by aherbert.

the class TSFPeakResultsReader method read.

/**
	 * Read the results from the TSF file into memory
	 * 
	 * @return The results set (or null if an error occurred)
	 */
public MemoryPeakResults read() {
    readHeader();
    if (spotList == null)
        return null;
    MemoryPeakResults results = createResults();
    // Read the messages that contain the spot data
    FileInputStream fi = null;
    try {
        fi = new FileInputStream(filename);
        // size of int + size of long
        fi.skip(12);
    } catch (Exception e) {
        System.err.println("Failed to read open TSF file: " + filename);
        e.printStackTrace();
        if (fi != null) {
            try {
                fi.close();
            } catch (IOException ex) {
            }
        }
        return null;
    }
    LocationUnits locationUnits = spotList.getLocationUnits();
    boolean locationUnitsWarning = false;
    float locationConversion = 1;
    switch(locationUnits) {
        case PIXELS:
            break;
        case NM:
            locationConversion = 1 / getNmPerPixel(results.getCalibration());
            break;
        case UM:
            float nmPerPixel = getNmPerPixel(results.getCalibration());
            if (nmPerPixel != 1) {
                float umPerPixel = nmPerPixel / 1000;
                locationConversion = 1 / umPerPixel;
            }
            break;
        default:
            System.err.println("Unsupported location units conversion: " + locationUnits);
    }
    IntensityUnits intensityUnits = spotList.getIntensityUnits();
    boolean intensityUnitsWarning = false;
    float intensityConversion = 1;
    switch(intensityUnits) {
        case COUNTS:
            break;
        case PHOTONS:
            intensityConversion = getGain(results.getCalibration());
            break;
        default:
            System.err.println("Unsupported intensity units conversion: " + intensityUnits);
    }
    final float bias = (results.getCalibration().hasBias()) ? (float) results.getCalibration().getBias() : 0f;
    ThetaUnits thetaUnits = spotList.getThetaUnits();
    float thetaConversion = 1;
    switch(thetaUnits) {
        case DEGREES:
            break;
        case RADIANS:
            thetaConversion = (float) (180.0 / Math.PI);
            break;
        default:
            System.err.println("Unsupported theta units conversion: " + thetaUnits);
    }
    FitMode fitMode = FitMode.ONEAXIS;
    if (spotList.hasFitMode())
        fitMode = spotList.getFitMode();
    final boolean filterPosition = position > 0;
    final boolean filterSlice = slice > 0;
    long expectedSpots = getExpectedSpots();
    try {
        long totalSpots = 0;
        while (fi.available() > 0 && (totalSpots < expectedSpots || expectedSpots == 0)) {
            totalSpots++;
            Spot spot = Spot.parseDelimitedFrom(fi);
            // Only read the specified channel, position, slice and fluorophore type
            if (spot.getChannel() != channel)
                continue;
            if (filterPosition && spot.hasPos() && spot.getPos() != position)
                continue;
            if (filterSlice && spot.hasSlice() && spot.getSlice() != slice)
                continue;
            if (spot.hasFluorophoreType() && spot.getFluorophoreType() != fluorophoreType)
                continue;
            if (spot.hasLocationUnits() && !locationUnitsWarning && spot.getLocationUnits() != locationUnits) {
                System.err.println("Spot message has different location units, the units will be ignored: " + spot.getLocationUnits());
                locationUnitsWarning = true;
            }
            if (spot.hasIntensityUnits() && !intensityUnitsWarning && spot.getIntensityUnits() != intensityUnits) {
                System.err.println("Spot message has different intensity units, the units will be ignored: " + spot.getIntensityUnits());
                intensityUnitsWarning = true;
            }
            // Required fields
            int frame = spot.getFrame();
            float[] params = new float[7];
            params[Gaussian2DFunction.X_POSITION] = spot.getX() * locationConversion;
            params[Gaussian2DFunction.Y_POSITION] = spot.getY() * locationConversion;
            params[Gaussian2DFunction.SIGNAL] = spot.getIntensity() * intensityConversion;
            if (spot.hasBackground()) {
                // The writer of the TSF file should have removed the bias (as documented in the format).
                // We can add it back for GSDC results
                params[Gaussian2DFunction.BACKGROUND] = spot.getBackground() * intensityConversion + bias;
            }
            // Support different Gaussian shapes
            if (fitMode == FitMode.ONEAXIS) {
                params[Gaussian2DFunction.X_SD] = params[Gaussian2DFunction.Y_SD] = spot.getWidth() / TSFPeakResultsWriter.SD_TO_FWHM_FACTOR;
            } else {
                if (!spot.hasA()) {
                    params[Gaussian2DFunction.X_SD] = params[Gaussian2DFunction.Y_SD] = spot.getWidth() / TSFPeakResultsWriter.SD_TO_FWHM_FACTOR;
                } else {
                    double a = Math.sqrt(spot.getA());
                    double sd = spot.getWidth() / TSFPeakResultsWriter.SD_TO_FWHM_FACTOR;
                    params[Gaussian2DFunction.X_SD] = (float) (sd * a);
                    params[Gaussian2DFunction.Y_SD] = (float) (sd / a);
                }
                if (fitMode == FitMode.TWOAXISANDTHETA && spot.hasTheta()) {
                    params[Gaussian2DFunction.SHAPE] = spot.getTheta() * thetaConversion;
                }
            }
            // We can use the original position in pixels used for fitting
            int origX = (spot.hasXPosition()) ? spot.getXPosition() : 0;
            int origY = (spot.hasYPosition()) ? spot.getYPosition() : 0;
            // Q. Should we use the required field 'molecule'?
            float origValue = 0;
            double error = 0;
            float noise = 0;
            float[] paramsStdDev = null;
            int id = 0;
            int endFrame = frame;
            if (isGDSC) {
                error = spot.getError();
                noise = spot.getNoise();
                id = spot.getCluster();
                origValue = spot.getOriginalValue();
                endFrame = spot.getEndFrame();
                if (spot.getParamsStdDevCount() != 0) {
                    paramsStdDev = new float[spot.getParamsStdDevCount()];
                    for (int i = 0; i < paramsStdDev.length; i++) paramsStdDev[i] = spot.getParamsStdDev(i);
                }
            } else // Use the standard cluster field for the ID
            if (spot.hasCluster()) {
                id = spot.getCluster();
            }
            // Allow storing any of the optional attributes
            AttributePeakResult peakResult = new AttributePeakResult(frame, origX, origY, origValue, error, noise, params, paramsStdDev);
            peakResult.setEndFrame(endFrame);
            peakResult.setId(id);
            if (spot.hasXPrecision() || spot.hasYPrecision()) {
                // Use the average. Note this is not the Euclidean distance since we divide by n!
                double p = 0;
                int n = 0;
                if (spot.hasXPrecision()) {
                    p = spot.getXPrecision() * spot.getXPrecision();
                    n++;
                }
                if (spot.hasYPrecision()) {
                    p += spot.getYPrecision() * spot.getYPrecision();
                    n++;
                }
                peakResult.setPrecision(Math.sqrt(p / n));
            }
            results.add(peakResult);
        }
    } catch (IOException e) {
        System.err.println("Failed to read Spot message");
        e.printStackTrace();
        // Only fail if there is a number of expected spots. 
        if (expectedSpots != 0) {
            System.err.println("Unexpected error in reading Spot messages, no results will be returned");
            return null;
        }
    } finally {
        if (fi != null) {
            try {
                fi.close();
            } catch (IOException e) {
            }
        }
    }
    return results;
}

Example 4 with Spot

use of gdsc.smlm.tsf.TaggedSpotFile.Spot in project GDSC-SMLM by aherbert.

the class TSFPeakResultsWriter method add.

/*
	 * (non-Javadoc)
	 * 
	 * @see gdsc.utils.fitting.results.PeakResults#add(int, int, int, float, double, float, float[], float[])
	 */
public void add(int peak, int origX, int origY, float origValue, double error, float noise, float[] params, float[] paramsStdDev) {
    if (out == null)
        return;
    Spot.Builder builder = Spot.newBuilder();
    builder.setMolecule(id.incrementAndGet());
    builder.setChannel(1);
    builder.setFluorophoreType(1);
    builder.setFrame(peak);
    builder.setXPosition(origX);
    builder.setYPosition(origY);
    setBackground(builder, params[Gaussian2DFunction.BACKGROUND]);
    builder.setIntensity(params[Gaussian2DFunction.SIGNAL]);
    builder.setX(params[Gaussian2DFunction.X_POSITION]);
    builder.setY(params[Gaussian2DFunction.Y_POSITION]);
    setWidth(params, builder);
    if (canComputePrecision) {
        double s = (params[Gaussian2DFunction.X_SD] + params[Gaussian2DFunction.Y_SD]) * 0.5 * nmPerPixel;
        float precision = (float) PeakResult.getPrecision(nmPerPixel, s, params[Gaussian2DFunction.SIGNAL] / gain, noise / gain, isEmCCD);
        builder.setXPrecision(precision);
        builder.setYPrecision(precision);
    }
    builder.setError(error);
    builder.setNoise(noise);
    builder.setOriginalValue(origValue);
    if (paramsStdDev != null)
        addNewParamsStdDev(builder, paramsStdDev);
    Spot spot = builder.build();
    writeResult(1, spot);
}

Example 5 with Spot

use of gdsc.smlm.tsf.TaggedSpotFile.Spot in project GDSC-SMLM by aherbert.

the class TSFPeakResultsWriter method add.

@Override
public void add(PeakResult result) {
    final float[] params = result.params;
    Spot.Builder builder = Spot.newBuilder();
    builder.setMolecule(id.incrementAndGet());
    builder.setChannel(1);
    builder.setFluorophoreType(1);
    builder.setFrame(result.getFrame());
    builder.setXPosition(result.origX);
    builder.setYPosition(result.origY);
    setBackground(builder, params[Gaussian2DFunction.BACKGROUND]);
    builder.setIntensity(params[Gaussian2DFunction.SIGNAL]);
    builder.setX(params[Gaussian2DFunction.X_POSITION]);
    builder.setY(params[Gaussian2DFunction.Y_POSITION]);
    setWidth(params, builder);
    if (result.hasPrecision()) {
        // Use the actual precision
        float precision = (float) result.getPrecision();
        builder.setXPrecision(precision);
        builder.setYPrecision(precision);
    } else if (canComputePrecision) {
        // Compute precision
        double s = (params[Gaussian2DFunction.X_SD] + params[Gaussian2DFunction.Y_SD]) * 0.5 * nmPerPixel;
        float precision = (float) PeakResult.getPrecision(nmPerPixel, s, params[Gaussian2DFunction.SIGNAL] / gain, result.noise / gain, isEmCCD);
        builder.setXPrecision(precision);
        builder.setYPrecision(precision);
    }
    if (result.hasId())
        builder.setCluster(result.getId());
    builder.setError(result.error);
    builder.setNoise(result.noise);
    if (result.hasEndFrame())
        builder.setEndFrame(result.getEndFrame());
    builder.setOriginalValue(result.origValue);
    if (result.paramsStdDev != null)
        addNewParamsStdDev(builder, result.paramsStdDev);
    Spot spot = builder.build();
    writeResult(1, spot);
}