use of uk.ac.sussex.gdsc.smlm.tsf.TSFProtos.Spot in project GDSC-SMLM by aherbert.
the class TsfPeakResultsWriter method add.
@Override
public void add(int peak, int origX, int origY, float origValue, double error, float noise, float meanIntensity, float[] params, float[] paramsStdDev) {
if (out == null) {
return;
}
final Spot.Builder builder = Spot.newBuilder();
builder.setMolecule(id.incrementAndGet());
builder.setChannel(1);
builder.setFluorophoreType(1);
builder.setFrame(peak);
builder.setXPosition(origX);
builder.setYPosition(origY);
setParam(params, builder);
builder.setError(error);
builder.setNoise(noise);
builder.setMeanIntensity(meanIntensity);
builder.setOriginalValue(origValue);
if (paramsStdDev != null) {
addNewParamStdDevs(builder, paramsStdDev);
}
final Spot spot = builder.build();
writeResult(1, spot);
}
use of uk.ac.sussex.gdsc.smlm.tsf.TSFProtos.Spot in project GDSC-SMLM by aherbert.
the class ResultsManagerTest method checkWriteTsfMatchesRead.
@SuppressWarnings("null")
private static void checkWriteTsfMatchesRead(RandomSeed seed, int channels, int slices, int positions, int types) {
Assumptions.assumeFalse(java.awt.GraphicsEnvironment.isHeadless());
final String filename = createFile();
FileOutputStream out = null;
try {
out = new FileOutputStream(filename);
} catch (final Exception ex) {
closeOutput(out);
ex.printStackTrace();
Assertions.fail(ex.getMessage());
}
// Write the offsets used in the Tsf format
try {
@SuppressWarnings("resource") final DataOutputStream dos = new DataOutputStream(out);
dos.writeInt(0);
dos.writeLong(0);
} catch (final IOException ex) {
closeOutput(out);
ex.printStackTrace();
Assertions.fail(ex.getMessage());
}
// Generate random spots
final UniformRandomProvider rand = RngUtils.create(seed.getSeed());
final int size = 100;
final Spot[] spots = new Spot[size];
for (int i = 1; i <= size; i++) {
final Spot.Builder builder = Spot.newBuilder();
builder.setChannel(nextInt(rand, channels));
builder.setSlice(nextInt(rand, slices));
builder.setPos(nextInt(rand, positions));
builder.setFluorophoreType(nextInt(rand, types));
// This is a required field but is ignored when reading
builder.setMolecule(i);
builder.setCluster(rand.nextInt(10));
builder.setFrame(nextInt(rand, 100));
builder.setXPosition(rand.nextInt(50));
builder.setYPosition(rand.nextInt(50));
builder.setBackground(rand.nextFloat());
builder.setIntensity(rand.nextFloat());
builder.setX(rand.nextFloat());
builder.setY(rand.nextFloat());
builder.setZ(rand.nextFloat());
builder.setWidth((float) (Gaussian2DFunction.SD_TO_FWHM_FACTOR * rand.nextDouble()));
final Spot spot = builder.build();
spots[i - 1] = spot;
try {
spot.writeDelimitedTo(out);
} catch (final IOException ex) {
closeOutput(out);
ex.printStackTrace();
Assertions.fail(ex.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 (final IOException ex) {
closeOutput(out);
ex.printStackTrace();
Assertions.fail(ex.getMessage());
}
// Record the SpotList message
final 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++) {
final FluorophoreType.Builder typeBuilder = FluorophoreType.newBuilder();
typeBuilder.setId(type);
typeBuilder.setDescription("Type " + type);
typeBuilder.setIsFiducial(rand.nextDouble() < 0.5);
builder.addFluorophoreTypes(typeBuilder.build());
}
final SpotList spotList = builder.build();
try {
spotList.writeDelimitedTo(out);
} catch (final IOException ex) {
ex.printStackTrace();
Assertions.fail(ex.getMessage());
} finally {
closeOutput(out);
}
// Write the offset to the SpotList message into the offset position
try (RandomAccessFile f = new RandomAccessFile(new File(filename), "rw")) {
f.seek(4);
f.writeLong(offset);
} catch (final Exception ex) {
ex.printStackTrace();
Assertions.fail(ex.getMessage());
}
// 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++) {
final 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());
final MemoryPeakResults in = ResultsManager.loadInputResults(ResultsManager.INPUT_FILE, false, null, null, new ResultsManager.FilenameLoadOption(filename));
checkEqual(spots, channel, slice, position, type, in);
}
}
}
}
}
use of uk.ac.sussex.gdsc.smlm.tsf.TSFProtos.Spot in project GDSC-SMLM by aherbert.
the class ResultsManagerTest method extract.
private static MemoryPeakResults extract(Spot[] spots, int channel, int slice, int position, int type) {
final MemoryPeakResults results = new MemoryPeakResults(PsfHelper.create(PSFType.ONE_AXIS_GAUSSIAN_2D));
for (final Spot spot : spots) {
if (spot.getChannel() == channel && spot.getSlice() == slice && spot.getPos() == position && spot.getFluorophoreType() == type) {
final int id = spot.getCluster();
final int startFrame = spot.getFrame();
final int origX = spot.getXPosition();
final int origY = spot.getYPosition();
final float origValue = 0;
final double error = 0;
final float noise = 0;
final float[] params = Gaussian2DPeakResultHelper.createOneAxisParams(spot.getBackground(), spot.getIntensity(), spot.getX(), spot.getY(), spot.getZ(), (float) (spot.getWidth() / Gaussian2DFunction.SD_TO_FWHM_FACTOR));
final float[] paramsStdDev = null;
final IdPeakResult peak = new IdPeakResult(startFrame, origX, origY, origValue, error, noise, 0, params, paramsStdDev, id);
results.add(peak);
}
}
return results;
}
use of uk.ac.sussex.gdsc.smlm.tsf.TSFProtos.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;
}
final MemoryPeakResults results = createResults();
// Used in the exception handler to check the correct number of spots were read
long expectedSpots = -1;
// Read the messages that contain the spot data
try (BufferedInputStream fi = new BufferedInputStream(new FileInputStream(filename))) {
// size of int + size of long
FileUtils.skip(fi, 12);
final LocationUnits locationUnits = spotList.getLocationUnits();
boolean locationUnitsWarning = false;
final IntensityUnits intensityUnits = spotList.getIntensityUnits();
boolean intensityUnitsWarning = false;
final FitMode fitMode = (spotList.hasFitMode()) ? spotList.getFitMode() : FitMode.ONEAXIS;
final boolean filterPosition = position > 0;
final boolean filterSlice = slice > 0;
// Set up to read two-axis and theta data
final PSF psf = PsfHelper.create(PSFType.TWO_AXIS_AND_THETA_GAUSSIAN_2D);
final int[] indices = PsfHelper.getGaussian2DWxWyIndices(psf);
final int isx = indices[0];
final int isy = indices[1];
final int ia = PsfHelper.getGaussian2DAngleIndex(psf);
int nparams = PeakResult.STANDARD_PARAMETERS;
switch(fitMode) {
case ONEAXIS:
nparams += 1;
break;
case TWOAXIS:
nparams += 2;
break;
case TWOAXISANDTHETA:
nparams += 3;
break;
default:
logger.log(Level.WARNING, () -> "Unknown fit mode: " + fitMode);
return null;
}
expectedSpots = getExpectedSpots();
long totalSpots = 0;
while (fi.available() > 0 && (totalSpots < expectedSpots || expectedSpots == 0)) {
totalSpots++;
final 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) {
logger.warning(() -> "Spot message has different location units, the units will be ignored: " + spot.getLocationUnits());
locationUnitsWarning = true;
}
if (spot.hasIntensityUnits() && !intensityUnitsWarning && spot.getIntensityUnits() != intensityUnits) {
logger.warning(() -> "Spot message has different intensity units, the units will be ignored: " + spot.getIntensityUnits());
intensityUnitsWarning = true;
}
// Required fields
final int frame = spot.getFrame();
final float[] params = new float[nparams];
params[PeakResult.X] = spot.getX();
params[PeakResult.Y] = spot.getY();
params[PeakResult.INTENSITY] = spot.getIntensity();
if (spot.hasBackground()) {
params[PeakResult.BACKGROUND] = spot.getBackground();
}
if (spot.hasZ()) {
params[PeakResult.Z] = spot.getZ();
}
// Support different Gaussian shapes
if (fitMode == FitMode.ONEAXIS) {
params[isx] = (float) (spot.getWidth() / Gaussian2DFunction.SD_TO_FWHM_FACTOR);
} else {
if (!spot.hasA()) {
params[isx] = params[isy] = (float) (spot.getWidth() / Gaussian2DFunction.SD_TO_FWHM_FACTOR);
} else {
final double a = Math.sqrt(spot.getA());
final double sd = spot.getWidth() / Gaussian2DFunction.SD_TO_FWHM_FACTOR;
params[isx] = (float) (sd * a);
params[isy] = (float) (sd / a);
}
if (fitMode == FitMode.TWOAXISANDTHETA && spot.hasTheta()) {
params[ia] = spot.getTheta();
}
}
// We can use the original position in pixels used for fitting
final int origX = (spot.hasXPosition()) ? spot.getXPosition() : 0;
final 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 meanIntensity = 0;
float[] paramsStdDev = null;
int id = Integer.MIN_VALUE;
int category = Integer.MIN_VALUE;
int endFrame = Integer.MIN_VALUE;
if (isGdsc) {
error = spot.getError();
noise = spot.getNoise();
meanIntensity = spot.getMeanIntensity();
id = spot.getCluster();
category = spot.getCategory();
origValue = spot.getOriginalValue();
endFrame = spot.getEndFrame();
if (spot.getParamStdDevsCount() != 0) {
paramsStdDev = new float[spot.getParamStdDevsCount()];
for (int i = 0; i < paramsStdDev.length; i++) {
paramsStdDev[i] = spot.getParamStdDevs(i);
}
}
} else if (spot.hasCluster()) {
// Use the standard cluster field for the ID.
// Note: The GDSC SMLM results do not distinguish molecule or cluster so here
// we pick cluster first.
id = spot.getCluster();
} else if (spot.hasMolecule()) {
// If no cluster then use the molecule.
// Note: This may just be a natural series with a unique number for each localisation.
id = spot.getMolecule();
}
// Allow storing any of the optional attributes
final AttributePeakResult peakResult = new AttributePeakResult(frame, origX, origY, origValue, error, noise, meanIntensity, params, paramsStdDev);
if (id != Integer.MIN_VALUE) {
peakResult.setId(id);
}
if (category != Integer.MIN_VALUE) {
peakResult.setCategory(category);
}
if (endFrame != Integer.MIN_VALUE) {
peakResult.setEndFrame(endFrame);
}
if (spot.hasXPrecision() || spot.hasYPrecision()) {
// Use the average. Note this is not the Euclidean distance since we divide by n
double sumSq = 0;
int count = 0;
if (spot.hasXPrecision()) {
sumSq = spot.getXPrecision() * spot.getXPrecision();
count++;
}
if (spot.hasYPrecision()) {
sumSq += spot.getYPrecision() * spot.getYPrecision();
count++;
}
peakResult.setPrecision(Math.sqrt(sumSq / count));
}
results.add(peakResult);
}
} catch (final IOException ex) {
logger.log(Level.WARNING, ex, () -> "Failed to read TSF file: " + filename);
if (expectedSpots == -1) {
// The exception was created during set-up.
return null;
}
// Only fail if there is a number of expected spots.
if (expectedSpots != 0) {
logger.warning(() -> "Unexpected error in reading Spot messages, no results will be returned");
return null;
}
}
return results;
}
use of uk.ac.sussex.gdsc.smlm.tsf.TSFProtos.Spot in project GDSC-SMLM by aherbert.
the class TsfPeakResultsWriter method addAll.
@Override
public void addAll(PeakResult[] results) {
if (out == null) {
return;
}
final Spot[] spots = new Spot[20];
int count = 0;
final Spot.Builder builder = Spot.newBuilder();
for (final PeakResult result : results) {
addStandardFields(builder, result);
addParamStdDevs(builder, result.getParameterDeviations());
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);
}
Aggregations