Skip to content
Snippets Groups Projects
Commit 57cdbba6 authored by Niko Ehrenfeuchter's avatar Niko Ehrenfeuchter :keyboard:
Browse files

Initial commit.

parents
No related branches found
No related tags found
No related merge requests found
target/
pom.xml 0 → 100644
<?xml version="1.0" encoding="UTF-8"?>
<project xmlns="http://maven.apache.org/POM/4.0.0"
xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance"
xsi:schemaLocation="http://maven.apache.org/POM/4.0.0
http://maven.apache.org/xsd/maven-4.0.0.xsd">
<modelVersion>4.0.0</modelVersion>
<parent>
<groupId>org.scijava</groupId>
<artifactId>pom-scijava</artifactId>
<version>17.1.1</version>
<relativePath />
</parent>
<groupId>ch.unibas.biozentrum.imcf</groupId>
<artifactId>NoiSee</artifactId>
<version>0.1.0-SNAPSHOT</version>
<name>NoiSee</name>
<description>An easy-to-use ImageJ macro for measuring SNR.</description>
<url>http://www.biozentrum.unibas.ch/noisee</url>
<inceptionYear>2018</inceptionYear>
<organization>
<name>University of Basel</name>
<url>http://www.unibas.ch/</url>
</organization>
<licenses>
<license>
<name>GNU General Public License (GPL) v3</name>
<distribution>repo</distribution>
</license>
</licenses>
<developers>
<!-- See https://imagej.net/Team -->
<developer>
<id>ehrenfeu</id>
<name>Niko Ehrenfeuchter</name>
<url>http://imagej.net/User:Ehrenfeu</url>
<roles>
<role>founder</role>
<role>lead</role>
<role>developer</role>
<role>debugger</role>
<role>reviewer</role>
<role>support</role>
<role>maintainer</role>
</roles>
</developer>
</developers>
<contributors>
<contributor>
<name>N/A</name>
</contributor>
</contributors>
<mailingLists>
<mailingList>
<name>ImageJ Forum</name>
<archive>http://forum.imagej.net/</archive>
</mailingList>
</mailingLists>
<scm>
<connection>scm:git:git://github.com/imcf/noisee</connection>
<developerConnection>scm:git:git@github.com:imcf/noisee</developerConnection>
<tag>HEAD</tag>
<url>https://github.com/imcf/noisee</url>
</scm>
<issueManagement>
<system>GitHub Issues</system>
<url>https://github.com/imcf/example-script-collection/issues</url>
</issueManagement>
<ciManagement>
<system>None</system>
</ciManagement>
<properties>
<package-name>ch.unibas.biozentrum.imcf.noisee</package-name>
<main-class>ch.unibas.biozentrum.imcf.noisee</main-class>
<license.licenseName>gpl_3</license.licenseName>
<license.copyrightOwners>University of Basel, Switzerland</license.copyrightOwners>
</properties>
<repositories>
<repository>
<id>imagej.public</id>
<url>http://maven.imagej.net/content/groups/public</url>
</repository>
</repositories>
<dependencies>
<dependency>
<groupId>net.imagej</groupId>
<artifactId>imagej</artifactId>
</dependency>
</dependencies>
</project>
// @Boolean(label="NoiSee - run SNR bead analysis",description="all open images will be closed",value="true") run_macro
// @File(label="Beads time-series image",description="2D time-lapse acquisition of fluorescent beads") beadsimage
// @Integer(label="Beads diameter (in pixels)",description="approximate bead diameter (in pixels)",value=15) beads_diameter
// @Integer(label="Find Maxima noise tolerance",description="typical values: [PMT=50] [HyD (photon counting)=10] [Camera=500]",value=50) beads_noisetolerance
// @Boolean(label="Create Kymographs",description="create kymograph plots from selected beads",value="true") make_kymographs
// @Boolean(label="Save additional measurements",description="store 'StdDev', 'SNR', 'Mean' and 'bleaching' measurements",value="false") save_measurements
// @Boolean(label="Save results as PDF",description="generate a PDF with images and plots",value="true") save_pdf
// @Boolean(label="Keep ROI images open",description="ROI visualizations will only be added to PDF otherwise",value="false") keep_roiimages
// @Integer(label="Log level",description="higher number means more messages",min=0,max=2,style="scroll bar") LOGLEVEL
// @Boolean(label="Save log messages",description="save contents of the 'Log' window in a text file",value="false") save_log
//////////// NoiSee ///////////////////////////////////////////////////////
// SNR evaluation macro, written by Kai Schleicher, Niko Ehrenfeuchter, IMCF Basel
// licence: GPLv3
/// TODOs:
// - test images with other dynamic ranges (16 bit)
// - consider using lacan's LibInstaller: https://github.com/lacan/LibInstaller
/// Naming conventions
// - "imgs_" - image IDs of stacks
// - "img_" - image IDs of 2D images
// - "rgb_" - image IDs of 2D RGB images
// - "roi_" - index number of a ROI manager entry
if (run_macro == false)
exit("Please select the 'run analysis' option to execute the macro!");
////////////////// function definitions ///////////////////////////////
function duplicateImage(id_orig, new_name, single){
// duplicate an image (given by ID), give it a new name and return the new ID
selectImage(id_orig);
run("Select None"); // clear selection to avoid messing up the duplication
if (single == true) {
run("Duplicate...", "use"); // duplicate single slice only
} else {
run("Duplicate...", "duplicate"); // duplicate the entire stack
}
rename(new_name);
return getImageID();
}
function duplicateAndClose(id_orig) {
// helper function for re-ordering images when creating a PDF
// the PDF exporter places images one per page in the order when they were created, so the only
// way to change that order is to duplicate an existing image into a new one, closing the old...
selectImage(id_orig);
new_id = duplicateImage(id_orig, getTitle(), false);
closeImage(id_orig);
return new_id;
}
function plotResult(Title, XLabel, YLabel, col) {
// create a high-res plot of the values in column "col" of the results table
// using the given title and labels
// returns the ID of the newly created plot image
if (nResults == 0) {
print("Unable to create plot [" + Title + "], results table empty!");
return;
}
print("Creating plot [" + Title + "] for " + nResults + " results.");
Yvalues=newArray(nResults);
Ymin = getResult(col, 0);
Ymax = Ymin;
// store Y-values in an array, get min and max:
for (i=0; i < nResults; i++) {
Yvalues[i] = getResult(col, i);
Ymax = maxOf(Ymax, Yvalues[i]);
Ymin = minOf(Ymin, Yvalues[i]);
}
// make sure not to over-emphasize data that has only little variance by
// setting the plot limits to min/max +/- 10:
Ymin = round(Ymin - 10.5);
Ymax = round(Ymax + 10.5);
print("Plot Y-axis range: " + Ymin + " - " + Ymax);
temp_plot_name = "NoiSee temporary plot window";
Plot.create(temp_plot_name, XLabel, YLabel);
Plot.setLimits(0, nResults, Ymin, Ymax);
Plot.setLineWidth(2);
Plot.setColor("#6688dd");
Plot.add("line",Yvalues);
Plot.setColor("red");
Plot.add("circles", Yvalues);
Plot.show();
Plot.makeHighResolution("plot of " + col, 2);
closeImage(temp_plot_name);
selectWindow("plot of " + col);
plot_id = getImageID();
addTextToImage(plot_id, true, "Center", 48, Title);
addTextToImage(plot_id, true, "Center", 24, "");
return plot_id;
}
function erode(image_id, num_iterations) {
// run the "erode" operation multiple times on an image / stack
selectImage(image_id);
logd("Eroding current image by " + num_iterations + " iterations...");
for (i=0; i < num_iterations; i++) {
run("Erode", "stack");
}
}
function binaryToSelection(image_id, roi_name, savepath) {
// create a selection from a binary image, add it to the ROI Manager with the given name and
// return the ID of the newly created ROI
selectImage(image_id);
if (is("binary") == false)
exit("Image " + getTitle() + " is not binary, not creating a selection!");
// "Create Selection" requires a threshold to be set, see the source code in
// https://imagej.nih.gov/ij/developer/source/ij/plugin/filter/ThresholdToSelection.java.html
// particularly for method selected(x, y) as the documentation is not very clear about this!
setThreshold(255, 255);
run("Create Selection");
return selectionToROI(roi_name, savepath);
}
function selectionToROI(roi_name, savepath) {
// add an existing selection as a new ROI with the given name and return the ID
roiManager("Add");
roi_id = roiManager("count") - 1; // ROIs are appended, hence "count" == newest ID + 1
roiManager("Select", roi_id);
roiManager("Rename", roi_name);
roiManager("Show None"); // clean up display behavior of the ROI manager
if (savepath != "")
roiManager("save selected", savepath + "/roi-" + roi_name + ".zip");
return roi_id;
}
function pickResults(count) {
// generate a list of index numbers creating a (more or less) equally spaced subset of the
// current results
if (nResults < count) {
setBatchMode("exit and display"); // exit batch mode and show images
exit("less than " + count + " results found, check your image!");
}
indexes = newArray(count);
halfstep = nResults / (count * 2);
for (i=0; i<count; i++) {
row = round((i*2+1) * halfstep);
indexes[i] = row;
}
return indexes;
}
function createCenterLines(image_id, res, len, namepfx) {
// create line selections of length "2*len" for the center of mass of a given subset of the
// current results ("res" is an array with the index numbers of results to be used)
roiManager("Set Color", "magenta");
run("Select None");
roiManager("reset"); // make sure the ROI Manager only contains our newly created line ROIs
selectImage(image_id);
for (i=0; i<res.length; i++) {
index = res[i];
cx = getResult("XM", index - 1);
cy = getResult("YM", index - 1);
logd("draw line for result " + index + " (c-o-m: " + cx + " / " + cy + ")");
makeLine(cx - len, cy, cx + len, cy);
roiManager("add");
roiManager("select", roiManager("count") - 1);
roiManager("rename", namepfx + index);
}
roiManager("Show all without labels");
roiManager("Set Line Width", 3);
}
function lineKymograph(image_id) {
// expects a number of "line" ROIs in the ROI Manager, creates a Kymograph (stack or orthogonal
// profile) for each of them, returning the newly created image IDs as an array
num_rois = roiManager("count");
ortho_imgs = newArray(num_rois);
for (i=0; i<num_rois; i++) {
selectImage(image_id);
roiManager("select", i);
roi_name = call("ij.plugin.frame.RoiManager.getName", i);
run("Reslice [/]...", "output=1.000 slice_count=1");
ortho_imgs[i] = getImageID();
rename("kymograph-" + roi_name);
}
return ortho_imgs;
}
function makeKymographMontage(imgs_beadmask, imgs_beads, diameter, basepath) {
// calculate the expected bead area to set the size parameter for "Analyze Particles" below:
beads_area = pow(diameter / 2, 2) * 3.14;
lo = beads_area * 0.6; // lower limit
hi = beads_area * 1.4; // upper limit
logi("range of expected beads areas: " + lo + " - " + hi);
img_tmp = duplicateImage(imgs_beadmask, "beadmask-T0", true);
run("Set Measurements...", "area center"); // center-of-mass gives us the bead location
run("Analyze Particles...", "size=" + lo + "-" + hi+ " display clear add");
resetThreshold();
closeImage(img_tmp);
saveAs("Results", basepath + "/individual-beads.txt");
logi("found " + nResults + " individual beads (excluding clustered / clumped ones)");
indexes = pickResults(4); // select a subset of four beads from the results
imgs_tmp = duplicateImage(imgs_beads, "tmp-stack", false);
createCenterLines(imgs_tmp, indexes, diameter, "bead-");
lineKymograph(imgs_tmp);
closeImage(imgs_tmp);
run("Images to Stack", "name=[kymograph-stack] title=[kymograph-] use");
imgs_tmp = getImageID();
run("Scale...", "x=16 y=16 z=1.0 depth=4 interpolation=None process create");
imgs_tmp2 = getImageID();
closeImage(imgs_tmp);
run("Make Montage...", "columns=2 rows=2 scale=1 border=40 font=16 label use");
closeImage(imgs_tmp2);
rename("Kymographs Montage");
run("cool");
return getImageID();
}
function findROI(roi_name) {
// scan for an ROI with the given name and return its ID
// WARNING: will terminate the macro if no ROI with the name exists!
for (roi_id=0; roi_id<roiManager("count"); roi_id++) {
cur_name = call("ij.plugin.frame.RoiManager.getName", roi_id);
if (cur_name == roi_name)
return roi_id;
}
exit("ERROR: ROI with name '" + roi_name + "' doesn't exist!");
}
function makeFilledROIImage(image_id, roi_name, inverse) {
roi_id = findROI(roi_name);
new_id = duplicateImage(image_id, "ROI - " + roi_name, true);
logd("creating filled ROI image for [" + getTitle() + "] and ROI [" + roi_name + "]");
run("RGB Color");
setForegroundColor(255, 57, 0);
roiManager("Select", roi_id);
if (inverse == true) {
run("Make Inverse");
}
// roiManager("Fill"); // WARNING: this command doesn't respect the "make inverse" from above!
run("Fill", "slice");
return new_id;
}
function drawROIsOnImage(image_id, namepfx, color) {
// use a given image (or stack) to create an RGB image where all ROIs starting with the given
// prefix are rendered in the specified color
// returns the ID of the new RGB image
tmp_id = duplicateImage(image_id, "tmp_image", true);
run("RGB Color");
// scan the ROI manager for the given prefix and assemble an array with corresponding IDs:
num_rois = roiManager("count");
roi_ids = newArray();
for (i=0; i<num_rois; i++) {
roiManager("select", i);
roi_name = call("ij.plugin.frame.RoiManager.getName", i);
if (startsWith(roi_name, namepfx)) {
roi_ids = Array.concat(roi_ids, i);
}
}
// Array.print(roi_ids);
// now select the IDs, combine them and flatten them into the new image:
roiManager("Select", roi_ids);
roiManager("Combine");
roiManager("Add", "00ff00", 0);
// the new ROI-ID equals "num_rois":
roiManager("Select", num_rois);
roiManager("Set Fill Color", color);
run("Flatten");
new_id = getImageID();
rename("Kymographs Selections");
roiManager("Delete");
closeImage(tmp_id);
return new_id;
}
function measureSelection(image_id, title, basepath, multi, savetxt) {
// measure (using the existing ROIs) on a given image and save the results
// in a text file at the given path if requested
selectImage(image_id);
run("Restore Selection");
run("Clear Results");
if (multi) {
roiManager("Multi Measure"); // measure on all slices of a stack
} else {
// run("Measure"); // is there any difference to roiManager("Measure") ??
roiManager("Measure");
}
if (savetxt == false)
return;
File.makeDirectory(basepath);
if (File.isDirectory(basepath) == false)
exit("ERROR creating directory [" + basepath + "], stopping.");
saveAs("Results", basepath + "/" + title + ".txt");
}
function flattenOverlay(image_id) {
selectImage(image_id);
title = getTitle();
run("Flatten");
new_id = getImageID();
closeImage(image_id);
selectImage(new_id);
rename(title);
return new_id;
}
function closeImage(image_id) {
if (isOpen(image_id) == false)
return;
selectImage(image_id);
close();
}
function dressImage(image_id, lut_name, enhance, desc) {
// applies a given LUT to an image, applies the "Enhance contrast" command to it if the
// "enhance" parameter is non-negative, adds a calibration bar (non-overlay) and optionally
// adds a text description below the image
// returns the ID of the new image, closing the old one
selectImage(image_id);
logi("dressImage(" + getTitle() + ", " + lut_name + ", " + enhance + ")");
run("Select None");
logd("applying LUT " + lut_name);
run(lut_name); // WARNING: applying a LUT changes the selected image!!
selectImage(image_id);
if (enhance >= 0) {
logd("enhancing contrast (saturated=" + enhance + ")");
run("Enhance Contrast", "saturated=" + enhance);
}
// sometimes calibration bar doesn't work depending on when / where the image window appears on
// the screen, so make sure the right window is selected / active:
selectImage(image_id);
logd("adding calibration bar");
title = getTitle();
run("Calibration Bar...","location=[Upper Right] fill=None " +
"label=White number=5 decimal=0 font=12 zoom=3.3 bold");
new_id = getImageID();
rename(title);
if (desc != "")
addTextToImage(new_id, false, "Center", 22, desc);
closeImage(image_id);
return new_id;
}
function createTable(title, pairs, fname) {
// create a new 2-column table with the given title, using the key-value tuples from the array
// given in "pairs" to fill the cells, saving it as to the text file "fname"
titleb = "[" + title + "]";
if (isOpen(title))
print(titleb, "\\Clear");
else
run("Table...", "name=" + titleb + " width=320 height=340");
print(titleb, "\\Headings:Key\tValue");
for (i = 0; i < pairs.length; i = i + 2) {
print(titleb, pairs[i] + "\t" + pairs[i+1]);
}
selectWindow(title);
if (fname != "")
saveAs("Results", fname);
}
function createTableImage(title, pairs) {
// create an image with a 2-column table from the given array of pairs
// returns the ID of the new image
newImage("tmp_table_col1", "RGB white", 1, 1, 1);
img_col1 = getImageID();
newImage("tmp_table_col2", "RGB white", 1, 1, 1);
img_col2 = getImageID();
for (i = 0; i < pairs.length; i = i + 2) {
addTextToImage(img_col1, false, "Left", 28, pairs[i]);
addTextToImage(img_col2, false, "Left", 28, pairs[i+1]);
}
selectImage(img_col2);
col2h = getHeight();
col2w = getWidth() + 28;
run("Canvas Size...", "width=" + col2w + " height=" + col2h + " position=Top-Right");
run("Combine...", "stack1=tmp_table_col1 stack2=tmp_table_col2");
img_table = getImageID();
rename(title);
addTextToImage(img_table, true, "Center", 18, "");
addTextToImage(img_table, true, "Center", 36, title);
return img_table;
}
function addTextToImage(image_id, above, justify, font_size, text) {
selectImage(image_id);
logd("addTextToImage(" + getTitle() + ", ''" + text + "'')");
run("Select None");
curw = getWidth();
neww = curw;
curh = getHeight();
newh = curh + font_size + 2;
setFont("SansSerif", font_size);
textw = getStringWidth(text);
if (textw > curw)
neww = textw;
if (above == true) {
canvpos = "Bottom-" + justify;
textpos = font_size + 2;
} else {
canvpos = "Top-" + justify;
textpos = newh;
}
setForegroundColor(30, 30, 30);
setBackgroundColor(255, 255, 255);
run("Canvas Size...", "width=" + neww + " height=" + newh + " position=" + canvpos);
drawString(text, 0, textpos);
}
function arrangeImages(ids, y_coord) {
// try to place the image windows given in the "ids" array horizontally on the screen,
// arranging them in equal distance and adjusting the zoom level if necessary
selectImage(ids[0]);
getLocationAndSize(_, _, window_width, window_height);
px_width = getWidth();
px_height = getHeight();
border_x = window_width - px_width;
border_y = window_height - px_height;
// try varius zoom levels until the windows fit on the screen and calculate the spacing:
for (s=0; s<4; s++) {
zoom = 1 - 0.25 * s;
width = px_width * zoom + border_x;
height = px_height * zoom + border_y;
total_width = ids.length * width;
logd("width " + width + " (zoom " + zoom + ", total width " + total_width + ")");
remaining = screenWidth - total_width;
if (remaining > 0) {
spacing = remaining / (ids.length + 1);
logd("zoom / width / spacing: " + zoom + " / " + width + " / " + spacing);
break;
} else {
continue;
}
// if we reach this point we can't seem to calculate a zoom factor, so we give up without
// placing / relocating the image windows (doesn't make sense):
return;
}
// now place the images:
for (i=0; i<ids.length; i++) {
selectImage(ids[i]);
x_coord = width * i + spacing * (i+1);
logd("placing image at " + x_coord + " / " + y_coord);
setLocation(x_coord, y_coord, width, height);
getLocationAndSize(_, _, window_width, window_height);
if (window_width != width) {
delta = (width - window_width) / 2;
setLocation(x_coord + delta, y_coord, width, height);
}
}
}
function stripOmeSuffix(orig) {
if (endsWith(orig, ".ome")) {
index = lastIndexOf(orig, ".ome");
orig = substring(orig, 0, index);
}
return orig;
}
function exit_show() {
setBatchMode("exit and display");
exit();
}
function log_formatter(message) {
print('image "' + getTitle() + "'' [type=" + bitDepth() + "] [id=" + getImageID() + "]: " +
message);
}
function logi(message) {
if (LOGLEVEL < 1)
return;
log_formatter(message);
}
function logd(message) {
if (LOGLEVEL < 2)
return;
log_formatter(message);
}
//////////////// set the user environment /////////////////////////////////////////////////////
// results will potentially be screwed up if other images are open, so close all:
if (nImages > 0) {
run("Close All");
}
print("\\Clear"); // clear the Log window
print("============================================");
print("NoiSee is published in Ferrand, Schleicher & Biehlmaier et al. 2018");
print("============================================");
print("running on ImageJ version " + getVersion);
run("Options...", "pad"); // pad edges when eroding
setOption("black background", true); // set "Black Background" in "Binary Options"
roiManager("reset"); // results are only correct if no previous ROI exists
run("Clear Results"); // empty the results table
run("Appearance...", " menu=0 16-bit=Automatic"); // disable inverting LUT
run("Colors...", "foreground=white background=black selection=red");
run("Input/Output...", "file=.txt"); // set saving format to .txt files
process_beads();
//////////////// bead method //////////////////////////////////////////////////////////////////
function process_beads() {
setBatchMode(true);
////////// open the image stack and prepare processing ////////// ////////// //////////
bfopts = " color_mode=Default";
bfopts += " view=Hyperstack";
bfopts += " stack_order=XYCZT";
run("Bio-Formats Importer", "open=[" + beadsimage + "]" + bfopts);
imgs_orig = getImageID();
fpath = File.getParent(beadsimage); // path only
fname = File.getName(beadsimage); // filename only
fname_nosuffix = stripOmeSuffix(File.nameWithoutExtension); // filename without extension
respath = fpath + "/" + fname_nosuffix + "_NoiSee-results"; // path for additional results
File.makeDirectory(respath);
// TODO: check if image dimensions meet our expectations (z=1, t>1)
print("processing image: " + fname + " (location: [" + fpath + "])");
if (bitDepth() > 8) {
print("image type " + bitDepth() + " bit detected, converting to 8 bit...");
resetMinAndMax();
run("8-bit");
}
// remove the scaling so all units (measurements, coordinates, ...) are pixel-based:
run("Set Scale...", "pixel=1 unit=pixel");
// make sure to reset any color LUT:
run("Grays");
getDimensions(width, height, channels, slices, frames);
area_full = width * height;
logi("dimensions: " + width + " x " + height + " (" + area_full + ")");
////////// ////////// ////////// ////////// ////////// ////////// ////////// //////////
////////// find beads and create a mask (binary image) ////////// ////////// //////////
imgs_beadmask = duplicateImage(imgs_orig, "beadmask", false);
run("Convert to Mask", "method=IsoData background=Dark calculate black");
roi_beads = binaryToSelection(imgs_beadmask, "beads", respath);
// now the background is 0 and beads are 255
////////// ////////// ////////// ////////// ////////// ////////// ////////// //////////
////////// create kymographs for a number of beads ////////// ////////// //////////
if (make_kymographs) {
rgb_kymo = makeKymographMontage(imgs_beadmask, imgs_orig, beads_diameter, respath);
addTextToImage(rgb_kymo, false, "Center", 22,
"Kymographs of selected beads. Each of them should be well-aligned from top to bottom.");
addTextToImage(rgb_kymo, false, "Center", 14, "Bleaching or variations in excitation " +
"power can be identified from changes in the intensities,");
addTextToImage(rgb_kymo, false, "Center", 14, "drifts in X result in a // or \\\\ shape " +
"of the Kymograph, drifts in Y will show as /\\ or \\/ shape.");
}
////////// ////////// ////////// ////////// ////////// ////////// ////////// //////////
////////// identify background and create a background-subtracted image ////////// //////////
imgs_bg = duplicateImage(imgs_beadmask, "background", false);
run("Invert", "stack"); // we want the background, not the beads...
erode(imgs_bg, 3); // erode the bead mask by 3 pixels to segment only real background
/// transfer selection to original image, measure mean and standard-deviation:
roi_bg = binaryToSelection(imgs_bg, "background", respath);
selectImage(imgs_orig);
run("Restore Selection");
getStatistics(area_bg, mean_bg, _, _, std_bg);
areapct_full_bg = area_bg / area_full * 100;
print("background stats: mean=" + mean_bg + " std=" + std_bg + " (" +
areapct_full_bg + "% of full area)");
/// create the background-subtracted image
imgs_bgsub = duplicateImage(imgs_orig, "background-subtracted", false);
run("Subtract...", "value=" + mean_bg);
////////// ////////// ////////// ////////// ////////// ////////// ////////// //////////
////////// calculate StdDev, Mean and SNR images ////////// ////////// ////////// //////////
selectImage(imgs_bgsub);
run("Z Project...", "projection=[Standard Deviation]");
img_std = getImageID();
rename("Standard Deviation (StdDev)");
selectImage(imgs_bgsub);
run("Z Project...", "projection=[Average Intensity]");
img_avg = getImageID();
rename("Average Intensity (Mean)");
// calculate the SNR as an image:
imageCalculator("Divide create 32-bit", img_avg, img_std);
img_snr = getImageID();
rename("Mean divided by StdDev (SNR)");
////////// ////////// ////////// ////////// ////////// ////////// ////////// //////////
////////// find saturated beads and exclude them ////////// ////////// ////////// //////////
/// TODO: discuss whether looking for saturation at this point makes sense,
/// or if this should be done on the original image (before the
/// background is subtracted and the time-series is averaged)
img_satmask = duplicateImage(img_avg, "saturation-mask", false);
run("8-bit"); // make sure the approach also works if input is not 8-bit
setThreshold(0, 253);
run("Convert to Mask", "method=Default background=Default");
// check if any close-to-saturation pixels exist at all:
getStatistics(_, _, min_satmask);
avg_has_saturated_pixels = true;
if (min_satmask == 255) {
print("No saturated pixels found in Average Intensity image, no need to exclude beads.");
avg_has_saturated_pixels = false;
} else {
print("Average Intensity image has saturated pixels, using a mask to exclude those beads!");
// erode by half of the given bead diameter, adding 30% safety margin:
erode(img_satmask, round((beads_diameter/2)*1.3));
roi_satmask = binaryToSelection(img_satmask, "saturation_mask", respath);
run("Divide...", "value=255 stack"); // create a binary mask with 1 and 0
setMinAndMax(0, 1);
// now apply the mask to the mean/average image:
imageCalculator("Multiply", img_avg, img_satmask);
}
closeImage(img_satmask); // the saturation mask is not required any more
////////// ////////// ////////// ////////// ////////// ////////// ////////// //////////
////////// find local maxima to determine locations for SNR evaluation ////////// //////////
// prepare measurements:
run("Set Measurements...", "min redirect=None decimal=2");
/// Find Maxima - noise parameter: Maxima are ignored if they do not stand out from the
/// surroundings by more than this value (calibrated units for calibrated images). In other
/// words, a threshold is set at the maximum value minus noise tolerance and the contiguous
/// area around the maximum above the threshold is analyzed. For accepting a maximum, this
/// area must not contain any point with a value higher at than the maximum. Only one maximum
/// within this area is accepted.
selectImage(img_avg);
run("Find Maxima...", "noise=" + beads_noisetolerance + " output=[Point Selection]");
roi_peaks = selectionToROI("peaks_from_mean", respath);
////////// ////////// ////////// ////////// ////////// ////////// ////////// //////////
////////// measure (evaulate) using the above determined points ////////// ////////// //////////
// DISCUSS: should this be normalized (bit depth) to allow comparing different sensors?
measureSelection(img_std, "StdDev", respath, false, save_measurements);
measureSelection(img_snr, "SNR", respath, false, save_measurements);
getStatistics(_, mean_snr, _, _, std_snr); // mean and std from the SNR
print("SNR stats: mean=" + mean_snr + " std=" + std_snr);
measureSelection(img_avg, "mean", respath, false, save_measurements);
getStatistics(_, mean_sig, _, _, std_sig); // mean and std from the signal
print("signal stats: mean=" + mean_sig + " std=" + std_sig);
mean_sbr = mean_sig / mean_bg; // calculate signal-to-background (SBR)
print("SBR (signal-to-background): " + mean_sbr);
// NOTE about error-propagation: all components contributing to the SBR contain an error and
// thus the SBR standard-deviation should be taken into consideration when further analyzing /
// using the mean SBR value!
nBeads = nResults;
print("number of beads detected / measured: " + nBeads);
////////// ////////// ////////// ////////// ////////// ////////// ////////// //////////
////////// measure bleaching ////////// ////////// ////////// ////////// ////////// //////////
// prepare measurements: area/mean/stddev to allow calculating the SEM = StdDev / sqrt(npixels)
run("Set Measurements...", "area mean standard redirect=None decimal=2");
imgs_bleach = duplicateImage(imgs_beadmask, "bleaching", false);
// erode the mask a bit to compensate for a possible xy-drift:
erode(imgs_bleach, 3);
// NOTE: if a scatter plot of stddev vs. time reveals an increasing trend this might indicate
// an xy-drift. try dilating the bead mask more in that case, if that does not help the
// time-series should be stabilized (e.g. using stackreg) before analyzing (make sure NOT to
// have sub-pixel correction as this will change the SNR!!)
roi_bleach = binaryToSelection(imgs_bleach, "bead_intensity_vs_time", respath);
closeImage(imgs_bleach); // we have the regions in the ROI manager, no need to keep the image
measureSelection(imgs_orig, "bead_intensity_vs_time", respath, true, save_measurements);
////////// ////////// ////////// ////////// ////////// ////////// ////////// //////////
/////////// create images showing detected ROIs //////////// ////////// ////////// //////////
rgb_roi_bg = makeFilledROIImage(imgs_orig, "background", false);
addTextToImage(rgb_roi_bg, false, "Center", 22,
"Red is representing the area detected as background.");
if (avg_has_saturated_pixels) {
rgb_roi_satmask = makeFilledROIImage(imgs_orig, "saturation_mask", true);
addTextToImage(rgb_roi_satmask, false, "Center", 22,
"Red areas contain saturation and are being ignored.");
} else {
rgb_roi_satmask = duplicateImage(imgs_orig, "ROI - saturation_mask", true);
addTextToImage(rgb_roi_satmask, false, "Center", 22,
"No saturated pixels found, no need to exclude any beads from analysis.");
run("RGB Color");
}
rgb_roi_intensity = makeFilledROIImage(imgs_orig, "bead_intensity_vs_time", false);
addTextToImage(rgb_roi_intensity, false, "Center", 22,
"Red pixels used for Mean / StdDev plots over time.");
if (make_kymographs) {
rgb_kymo_lines = drawROIsOnImage(imgs_orig, "bead-", "green");
addTextToImage(rgb_kymo_lines, false, "Center", 22,
"Green lines showing selections used to create Kymographs.");
}
////////// ////////// ////////// ////////// ////////// ////////// ////////// //////////
/////////// decorate images, add calibration bars //////////// ////////// ////////// //////////
rgb_t0 = dressImage(imgs_orig, "Glow", -1, "First timepoint (T0) of the original data.");
rgb_avg = dressImage(img_avg, "Fire", -1, "Average projection of the original data.");
rgb_std = dressImage(img_std, "Fire", -1, "StdDev projection of the original data.");
rgb_snr = dressImage(img_snr, "physics", 0.35, "Signal-To-Noise ratio for all pixels.");
// adjust name for original image (now showing the first time point):
addTextToImage(rgb_t0, false, "Center", 12, "Original file name:");
addTextToImage(rgb_t0, false, "Center", 12, fname);
rename("Original Data");
// clean up image windows:
closeImage(imgs_bgsub);
closeImage(imgs_bg);
closeImage(imgs_beadmask);
////////// ////////// ////////// ////////// ////////// ////////// ////////// //////////
////////// create plots from bleaching measurements ////////// ////////// ////////// //////////
// linear decay = z-drift, exponential decay = bleaching
col_mean = "Mean(bead_intensity_vs_time)";
plot_mean = plotResult("Z-drift or bleaching?", "Frame", "bead intensity mean", col_mean);
addTextToImage(plot_mean, false, "Center", 22,
"A linear decay would be an indicator for a z-drift, an exponential decay for bleaching.");
addTextToImage(plot_mean, false, "Center", 96, "");
// linear increase = xy-drift
// this reveals only drift that is relvant for the regions that measure bleaching or z-drift.
// could in general be used to reveal xy-drift as it also interferes with SNR
col_std = "StdDev(bead_intensity_vs_time)";
plot_std = plotResult("Stable Z-drift / bleaching measurement?", "Frame",
"bead intensity StdDev", col_std);
addTextToImage(plot_std, false, "Center", 22,
"Increase in the StdDev indicates a drift in X/Y.");
run("Combine...", "stack1=[plot of " + col_mean + "] stack2=[plot of " + col_std + "] combine");
rename("Plots of intensity vs. time");
rgb_plots = getImageID();
////////// ////////// ////////// ////////// ////////// ////////// ////////// //////////
/////////// write summary of results and settings into a new table //////////// //////////
title = "NoiSee Results Summary";
pairs = newArray();
pairs = Array.concat(pairs, "SNR mean", mean_snr);
pairs = Array.concat(pairs, "SNR StdDev", std_snr);
pairs = Array.concat(pairs, "SBR mean", mean_sbr);
pairs = Array.concat(pairs, "Signal mean", mean_sig);
pairs = Array.concat(pairs, "Signal StdDev (noise)", std_sig);
pairs = Array.concat(pairs, "Background mean", mean_bg);
pairs = Array.concat(pairs, "Background StdDev", std_bg);
pairs = Array.concat(pairs, "Background area", area_bg);
pairs = Array.concat(pairs, "Background area %", areapct_full_bg);
pairs = Array.concat(pairs, "Number of beads", nBeads);
pairs = Array.concat(pairs, "", "");
pairs = Array.concat(pairs, "-- Settings --", "");
pairs = Array.concat(pairs, "Beads diameter", beads_diameter);
pairs = Array.concat(pairs, "Noise tolerance", beads_noisetolerance);
pairs = Array.concat(pairs, "", "");
pairs = Array.concat(pairs, "-- Input --", "");
pairs = Array.concat(pairs, "Filename", fname);
createTable(title, pairs, respath + "/" + "NoiSee-summary.txt");
img_summary = createTableImage(title, pairs);
////////// ////////// ////////// ////////// ////////// ////////// ////////// //////////
////////// adjust order of the images and create a PDF ////////// ////////// //////////
// the PDF exporter of ImageJ simply concatenates all open images in ascending order,
// starting with the oldest one first, therefore we need to duplicate all images into
// new ones (closing the old ones) in the order which is desired in the resulting PDF:
rgb_t0 = duplicateAndClose(rgb_t0);
rgb_avg = duplicateAndClose(rgb_avg);
rgb_std = duplicateAndClose(rgb_std);
rgb_snr = duplicateAndClose(rgb_snr);
if (make_kymographs){
rgb_kymo = duplicateAndClose(rgb_kymo);
rgb_kymo_lines = duplicateAndClose(rgb_kymo_lines);
}
rgb_roi_bg = duplicateAndClose(rgb_roi_bg);
rgb_roi_satmask = duplicateAndClose(rgb_roi_satmask);
rgb_roi_intensity = duplicateAndClose(rgb_roi_intensity);
rgb_plots = duplicateAndClose(rgb_plots);
img_summary = duplicateAndClose(img_summary);
////////// ////////// ////////// ////////// ////////// ////////// ////////// //////////
////////// save the log and arrange the windows on the screen ////////// ////////// //////////
if (save_log) {
selectWindow("Log");
saveAs("Text", respath + "/Log.txt");
}
setBatchMode("exit and display"); // exit batch mode and show images
// wait(100); // give the OS some time to display all image windows
if (save_pdf) {
fname_pdf = respath + "/" + "NoiSee-report.pdf";
print("Creating PDF from images: " + fname_pdf);
run("PDF ... ", "show show scale save one save=[" + fname_pdf + "]");
}
closeImage(img_summary);
row_upper = newArray(rgb_t0, rgb_avg, rgb_plots);
row_lower = newArray(rgb_std, rgb_snr);
if (make_kymographs) {
row_lower = Array.concat(row_lower, rgb_kymo);
}
if (keep_roiimages) {
row_upper = Array.concat(row_upper, rgb_roi_intensity, rgb_roi_bg);
if (make_kymographs) {
row_lower = Array.concat(row_lower, rgb_kymo_lines);
}
row_lower = Array.concat(row_lower, rgb_roi_satmask);
} else {
closeImage(rgb_kymo_lines);
closeImage(rgb_roi_bg);
closeImage(rgb_roi_satmask);
closeImage(rgb_roi_intensity);
}
arrangeImages(row_upper, screenHeight * 0.1);
arrangeImages(row_lower, screenHeight * 0.5);
////////// ////////// ////////// ////////// ////////// ////////// ////////// //////////
}
// @Boolean(label="NoiSee - run Fluorescein SNR analysis",value="true") run_macro
// @File(label="Dark image",description="dark field image") darkimage
// @File(label="Fluorescein image",description="fluorescein image") fluoimage
// @Boolean(label="Save results as PDF",description="generate a PDF with images and plots",value="true") save_pdf
// @Integer(label="Log level",description="higher number means more messages",min=0,max=2,style="scroll bar") LOGLEVEL
//////////// NoiSee ///////////////////////////////////////////////////////
// SNR evaluation macro, written by Kai Schleicher, Niko Ehrenfeuchter, IMCF Basel
// licence: GPLv3
/// TODOs:
// - create a summary image with results to be included in PDF
/// Naming conventions
// - "img_" - image IDs of 2D images
// - "rgb_" - image IDs of 2D RGB images
if (run_macro == false)
exit("Please select the 'run analysis' option to execute the macro!");
////////////////// function definitions ///////////////////////////////
function duplicateImage(id_orig, new_name, single){
// duplicate an image (given by ID), give it a new name and return the new ID
selectImage(id_orig);
run("Select None"); // clear selection to avoid messing up the duplication
if (single == true) {
run("Duplicate...", "use"); // duplicate single slice only
} else {
run("Duplicate...", "duplicate"); // duplicate the entire stack
}
rename(new_name);
return getImageID();
}
function duplicateAndClose(id_orig) {
// helper function for re-ordering images when creating a PDF
// the PDF exporter places images one per page in the order when they were created, so the only
// way to change that order is to duplicate an existing image into a new one, closing the old...
selectImage(id_orig);
new_id = duplicateImage(id_orig, getTitle(), false);
closeImage(id_orig);
return new_id;
}
function closeImage(image_id) {
if (isOpen(image_id) == false)
return;
selectImage(image_id);
close();
}
function makeHistogram(image_id, lut_name, logarithmic) {
// applies a given LUT to the image specified with image_id, and creates a
// new image containing the histogram of the given one, optionally enabling
// the "logarithmic" view
selectImage(image_id);
orig_title = getTitle();
logi("makeHistogram(" + orig_title + ", " + lut_name + ")");
run("Select None");
logd("applying LUT " + lut_name);
run(lut_name); // WARNING: applying a LUT changes the selected image!!
selectImage(image_id);
if (logarithmic == true)
setKeyDown("shift");
run("Histogram");
setKeyDown("none");
// we have to rename the histogram as by default it will discard everything
// after the first space from the original title, so:
rename("Histogram of " + orig_title);
return getImageID();
}
function dressImage(image_id, lut_name, enhance) {
// applies a given LUT to an image, applies the "Enhance contrast" command to it if the
// "enhance" parameter is non-negative, adds a calibration bar (non-overlay)
// returns the ID of the new image, closing the old one
selectImage(image_id);
logi("dressImage(" + getTitle() + ", " + lut_name + ", " + enhance + ")");
run("Select None");
logd("applying LUT " + lut_name);
run(lut_name); // WARNING: applying a LUT changes the selected image!!
selectImage(image_id);
if (enhance >= 0) {
logd("enhancing contrast (saturated=" + enhance + ")");
run("Enhance Contrast", "saturated=" + enhance);
}
// sometimes calibration bar doesn't work depending on when / where the image window appears on
// the screen, so make sure the right window is selected / active:
selectImage(image_id);
logd("adding calibration bar");
title = getTitle();
run("Calibration Bar...","location=[Upper Right] fill=None " +
"label=White number=5 decimal=0 font=12 zoom=3.3 bold");
new_id = getImageID();
rename(title);
closeImage(image_id);
return new_id;
}
function arrangeImages(ids, y_coord) {
// try to place the image windows given in the "ids" array horizontally on the screen,
// arranging them in equal distance and adjusting the zoom level if necessary
selectImage(ids[0]);
getLocationAndSize(_, _, window_width, window_height);
px_width = getWidth();
px_height = getHeight();
border_x = window_width - px_width;
border_y = window_height - px_height;
// try varius zoom levels until the windows fit on the screen and calculate the spacing:
for (s=0; s<4; s++) {
zoom = 1 - 0.25 * s;
width = px_width * zoom + border_x;
height = px_height * zoom + border_y;
total_width = ids.length * width;
logd("width " + width + " (zoom " + zoom + ", total width " + total_width + ")");
remaining = screenWidth - total_width;
if (remaining > 0) {
spacing = remaining / (ids.length + 1);
logd("zoom / width / spacing: " + zoom + " / " + width + " / " + spacing);
break;
} else {
continue;
}
// if we reach this point we can't seem to calculate a zoom factor, so we give up without
// placing / relocating the image windows (doesn't make sense):
return;
}
// now place the images:
for (i=0; i<ids.length; i++) {
selectImage(ids[i]);
x_coord = width * i + spacing * (i+1);
logd("placing image at " + x_coord + " / " + y_coord);
setLocation(x_coord, y_coord, width, height);
}
}
function stripOmeSuffix(orig) {
if (endsWith(orig, ".ome")) {
index = lastIndexOf(orig, ".ome");
orig = substring(orig, 0, index);
}
return orig;
}
function log_formatter(message) {
print('image "' + getTitle() + "'' [type=" + bitDepth() + "] [id=" + getImageID() + "]: " +
message);
}
function logi(message) {
if (LOGLEVEL < 1)
return;
log_formatter(message);
}
function logd(message) {
if (LOGLEVEL < 2)
return;
log_formatter(message);
}
//////////////// set the user environment /////////////////////////////////////////////////////
// results will potentially be screwed up if other images are open, so close all:
if (nImages > 0) {
run("Close All");
}
print("\\Clear"); // clear the Log window
print("============================================");
print("NoiSee is published in Ferrand, Schleicher & Biehlmaier et al. 2017");
print("============================================");
print("running on ImageJ version " + getVersion);
run("Options...", "pad"); // pad edges when eroding
setOption("black background", true); // set "Black Background" in "Binary Options"
// roiManager("reset"); // results are only correct if no previous ROI exists
run("Clear Results"); // empty the results table
run("Appearance...", " menu=0 16-bit=Automatic"); // disable inverting LUT
run("Colors...", "foreground=white background=black selection=red");
run("Input/Output...", "file=.txt"); // set saving format to .txt files
process_fluo();
//////////////// fluorescein method ////////////////////////////////////////////////////////////////
function process_fluo() {
setBatchMode(true);
////////// open the image stack and prepare processing ////////// ////////// //////////
bfopts = " color_mode=Default";
bfopts += " view=Hyperstack";
bfopts += " stack_order=XYCZT";
run("Bio-Formats Importer", "open=[" + darkimage + "]" + bfopts);
img_bg = getImageID();
bg_fname = File.getName(darkimage); // filename only
bg_fname_nosuffix = stripOmeSuffix(File.nameWithoutExtension); // filename without extension
run("Bio-Formats Importer", "open=[" + fluoimage + "]" + bfopts);
img_fluo = getImageID();
fluo_fname = File.getName(fluoimage); // filename only
fluo_fname_nosuffix = stripOmeSuffix(File.nameWithoutExtension); // filename without extension
fpath = File.getParent(fluoimage); // path only
respath = fpath + "/" + fluo_fname_nosuffix + "_NoiSee-results"; // path for additional results
File.makeDirectory(respath);
print("processing images in location: [" + fpath + "]");
print(" - dark image: [" + bg_fname + "]");
print(" - fluorescein image: [" + fluo_fname + "]");
logd("\n\nfull paths of images (dark + fluorescein):\n" +
darkimage + "\n \n" +
fluoimage + "\n ");
print("path for saving results: [" + respath + "]");
////////// ////////// ////////// ////////// ////////// ////////// ////////// //////////
////////// get stats, calculate normalized image, calculate snr + sbr ////////// //////////
selectImage(img_bg);
getStatistics(_, mean_bg, _, _, std_bg);
print("darkimage stats: mean=" + mean_bg + " std=" + std_bg);
img_sig = duplicateImage(img_fluo, "absolute signal (fluorescein - background)", true);
run("Subtract...", "value=" + mean_bg);
getStatistics(_, mean_sig, _, _, std_sig);
print("absolute signal image stats: mean=" + mean_sig + " std=" + std_sig);
snr = mean_sig / std_sig;
sbr = mean_sig / mean_bg;
print("snr (signal-to-noise): " + snr);
print("sbr (signal-to-background): " + sbr);
////////// ////////// ////////// ////////// ////////// ////////// ////////// //////////
/////////// write summary of results and settings into a new table //////////// //////////
title = "NoiSee Results Summary";
titleb = "[" + title + "]";
if (isOpen(title))
print(titleb, "\\Clear");
else
run("Table...", "name=" + titleb + " width=320 height=340");
print(titleb, "\\Headings:Label\tValue");
print(titleb, "Signal mean\t" + mean_sig);
print(titleb, "Signal StdDev (noise)\t" + std_sig);
print(titleb, "SNR\t" + snr);
print(titleb, "Background mean\t" + mean_bg);
print(titleb, "SBR\t" + sbr);
print(titleb, "Background StdDev\t" + std_bg);
print(titleb, "\t");
print(titleb, "dark image\t" + bg_fname);
print(titleb, "fluorescein image\t" + fluo_fname);
selectWindow(title);
saveAs("Results", respath + "/" + "NoiSee-summary.txt");
////////// ////////// ////////// ////////// ////////// ////////// ////////// //////////
/////////// decorate images, add calibration bars //////////// ////////// ////////// //////////
rgb_hist = makeHistogram(img_sig, "physics", true);
rgb_sig = dressImage(img_sig, "physics", 0.35);
rgb_fluo = dressImage(img_fluo, "Glow", -1);
rgb_dark = dressImage(img_bg, "Glow", -1);
////////// ////////// ////////// ////////// ////////// ////////// ////////// //////////
////////// arrange the windows and create a PDF of images and plots ////////// //////////
rgb_sig = duplicateAndClose(rgb_sig);
rgb_hist = duplicateAndClose(rgb_hist);
setBatchMode("exit and display"); // exit batch mode and show images
// wait(100); // give the OS some time to display all image windows
if (save_pdf) {
fname_pdf = respath + "/" + "NoiSee-report.pdf";
print("Creating PDF from images: " + fname_pdf);
run("PDF ... ", "show show scale save one save=[" + fname_pdf + "]");
}
image_ids = newArray(rgb_sig, rgb_fluo);
arrangeImages(image_ids, screenHeight * 0.1);
image_ids = newArray(rgb_dark, rgb_hist);
arrangeImages(image_ids, screenHeight * 0.55);
////////// ////////// ////////// ////////// ////////// ////////// ////////// //////////
}
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment