Commit 6f561d8a authored by lombardo's avatar lombardo
Browse files

remove unneeded file

parent 7427cddd
#@ImagePlus (label="Sequence to analyze") imp0
#@int (label="channel to analyze") channel_analyze
#@File (label="Illumination region") roi_illumination_file
#@int (label="frame to start analysis") tStart
#@int (label="frame to measure signal increase") tIncreased
#@int (label="time averaging window") dtWindow
#@boolean (label="auto detect illumination frame") autoDetectIlluminationTime
#@int (label="manual illumination frame", required="false") t_illum
# author: Benoit Lombardot @ Scientific Computing Facility, MPI-CBG, Dresden Germany
# License: BSD, see the LICENSE.txt file coming with the script
# 2016-12-13 : version 1.2, correct measures of the background intensity and add measure of the bacground standard deviation, info are added to the script info
# 2017-02-28 : version 1.2, shift the window to create the mask 1 timestep earlier
# 2017-04-10 : version 1.2.1, remove nChannel parameter, (if necessary use bioformat to open the image correctly before starting the analysis)
from ij import IJ
from ij.io import RoiDecoder
from ij.plugin import ZProjector
from ij.plugin import Duplicator
from ij.plugin import ImageCalculator
from ij.plugin.frame import RoiManager
from ij.gui import Overlay
from ij.plugin import HyperStackConverter
from ij.plugin.filter import ParticleAnalyzer as PA
from ij.gui import Wand
from ij.process import FloodFiller
from ij.gui import PolygonRoi
from ij.gui import ShapeRoi
from ij.gui import Roi
from ij import ImagePlus
from ij.plugin.filter import ThresholdToSelection
from ij.process import ImageProcessor
from java.awt import Color
# load illumination region roi: roi_illum (blue)
# create image that will be segmented: impToSeg
# create a cell mask: roi_Cells (cyan)
# create positive mask in cells: roi_postive (red)
# and negative mask in cells: can be detected from the other region
###############################################
# reshape the image sequence if necessary #####
imp0.hide()
dims = imp0.getDimensions()
if dims[4]==1: #if only one time step exchange t and z
#reshape exchange time and z
imp0 = HyperStackConverter.toHyperStack(imp0, dims[2], dims[4], dims[3], "default", "Color");
IJ.log("Warning: only on time frame, time and z dimension were exchanged ");
imp0.show()
dims = imp0.getDimensions()
if tIncreased==-1 :
tIncreased=dims[4]
#############################################
# detect illumination time ##################
roi_illum = RoiDecoder.open(roi_illumination_file.getPath())
tRange = range(dims[4])
if autoDetectIlluminationTime :
t_illum= tRange[0]
I_illum_max = 0
z=1
IIllumSeq = []
for t in tRange:
idx1 = imp0.getStackIndex(channel_analyze,z,t+1)
ip = imp0.getStack().getProcessor(idx1)
ip.setRoi(roi_illum)
IIllumSeq.append( ip.getStatistics().mean )
aux = list(IIllumSeq) # duplicate the list
aux.sort()
IIllumMedian = aux[int(dims[4]/2)];
IIllumMax = max(IIllumSeq)
IThresh = IIllumMedian + (IIllumMax-IIllumMedian)/2
tOutlier = []
for t in tRange:
I = IIllumSeq[t]
if I>IThresh:
tOutlier.append( t )
t_illum = min(tOutlier)+1;
roi_illum.setProperty("tIllum", str(t_illum) )
dtWindow = min(dtWindow , t_illum-tStart, tIncreased-t_illum)
##################################################
# create an image of signal before segmentation ##
def doProjection(imp0, ch, tStart, tEnd):
imp = Duplicator().run(imp0,ch,ch,1,1,tStart, tEnd)
zproj = ZProjector(imp);
zproj.setMethod(ZProjector.AVG_METHOD);
zproj.setStartSlice(tStart);
zproj.setStopSlice(tEnd);
zproj.doProjection();
return zproj.getProjection();
imp_init0 = doProjection(imp0, channel_analyze, t_illum-dtWindow-1, t_illum-2)
impToSeg = Duplicator().run(imp_init0)
imp_end0 = doProjection(imp0, channel_analyze, tIncreased-dtWindow+1, tIncreased)
imp_dif0 = ImageCalculator().run('SUBTRACT create 32-bit', imp_end0, imp_init0)
stat_dif0 = imp_dif0.getStatistics();
#################################################
# create a mask of cell #########################
cellMask = Duplicator().run(impToSeg)
IJ.run(cellMask,"Log","")
IJ.setAutoThreshold(cellMask, "Triangle dark");
IJ.run(cellMask, "Convert to Mask", "");
IJ.run(cellMask, "Analyze Particles...", "size=200-Infinity show=Masks in_situ");
# clean the mask
'''
# 4 connectedness does not work one need to create a function with the wand and the foold filler
options = PA.FOUR_CONNECTED + PA.SHOW_MASKS + PA.IN_SITU_SHOW + PA.ADD_TO_MANAGER
measurements = 0
rt= None
minSize=20000
maxSize = imp0.getWidth()*imp0.getHeight()
analyzer = PA( options, measurements, rt, minSize, maxSize)
analyzer.analyze(cellMask)
cellMask = analyzer.getOutputImage()
cellMask.show()
'''
def flood_fill_labelling(imp, minSize, maxSize):
width = imp.getWidth();
ip = imp.getProcessor().convertToFloatProcessor();
wand = Wand(ip);
ff = FloodFiller(ip)
roiList = []
count=255;
pix_data = ip.getPixels();
for i in range(0, len(pix_data)):
#print pix_data[i]
if (pix_data[i]==255):
y = i/width;
x = i%width;
count = count+1;
ip.setValue(count);
#ip.fill(roi);
ff.fill(x,y)
wand.autoOutline(x,y,count,count, Wand.FOUR_CONNECTED)
if not wand.npoints == 0 :
roi = PolygonRoi(wand.xpoints, wand.ypoints, wand.npoints, PolygonRoi.POLYGON)
rect0 = roi.getBounds().width * roi.getBounds().height
rois = ShapeRoi(roi).getRois() # could be used to check which roi is the largest
roi_stats = roi.getStatistics()
if (roi_stats.area > minSize) & (roi_stats.area < maxSize) :
roiList.append(roi);
else :
ip.setValue(0)
ip.fill(roi)
for i in range(0, len(pix_data)):
if (pix_data[i]>255):
pix_data[i] = pix_data[i]-255;
return ImagePlus('labels', ip), roiList;
# remove small 4 connected regions in the foreground
# (otherwise applying a split later on with the roi manager will create tones of micro regions)
minSize = 500 # would be better to have the fraction of a cell size
maxSize = imp0.getWidth() * imp0.getHeight()
cellMask, roiList = flood_fill_labelling(cellMask, minSize, maxSize)
IJ.setThreshold(cellMask,1,cellMask.getStatistics().max)
IJ.run(cellMask, "Convert to Mask", "");
# remove small 4 connected regions in the background
IJ.run(cellMask, "Invert", "");
minSize = 20
maxSize = imp0.getWidth() * imp0.getHeight()
cellMask, roiList = flood_fill_labelling(cellMask, minSize, maxSize)
IJ.setThreshold(cellMask,1,cellMask.getStatistics().max)
IJ.run(cellMask, "Convert to Mask", "");
IJ.run(cellMask, "Invert", "");
ip = cellMask.getProcessor()
ip.setThreshold(128,255,ImageProcessor.NO_LUT_UPDATE)
roi_cells = ThresholdToSelection().convert(ip)
#######################################################
# create a mask of cell region with increased signal ##
imp_maskPositive = Duplicator().run(imp_dif0)
IJ.setThreshold(imp_maskPositive, 0, stat_dif0.max)
IJ.run(imp_maskPositive, "Convert to Mask", "")
imp_maskPositive = ImageCalculator().run('AND create', imp_maskPositive, cellMask)
ip = imp_maskPositive.getProcessor()
ip.setThreshold(128,255,ImageProcessor.NO_LUT_UPDATE)
roi_positive = ThresholdToSelection().convert(ip)
################################################
# display / output results of the preparation ##
roi_illum.setStrokeColor( Color.BLUE )
roi_illum.setName('illuminationRegion')
roi_cells.setStrokeColor( Color.CYAN )
roi_cells.setName('cellMask')
roi_positive.setFillColor( Color(1.0,0.0,0.0,0.3))
roi_positive.setName('increasedSignalRegion')
roiManager = RoiManager.getInstance()
if roiManager==None :
roiManager = RoiManager()
roiManager.reset()
roiManager.runCommand(imp0,"Show None");
roiManager.addRoi(roi_illum)
roiManager.addRoi(roi_positive)
roiManager.addRoi(roi_cells)
IJ.run(impToSeg, "Log", "")
IJ.run(impToSeg, "Enhance Contrast", "saturated=0.35")
impToSeg.setTitle("Cell segmentation")
impToSeg.show()
roiManager.runCommand(impToSeg,"Show All without labels")
# get Background intensity and set it in imp0 properties
# invert the roi
s1 = ShapeRoi(roi_cells);
s2 = ShapeRoi( Roi(0,0, imp0.getWidth(), imp0.getHeight()));
inv_roi_cell = s1.xor(s2)
def getMeanAndStdDev(imp, roi, frames, ch, z):
if len(frames) == 1 :
frame = frames[0]
idx = imp.getStackIndex(ch,z,frame)
ip = imp.getStack().getProcessor(idx)
ip.setRoi(roi)
stats = ip.getStatistics()
ip.resetRoi()
means = stats.mean
stDevs = stats.stdDev
elif len(frames)>1 :
stdDevs = []
means = []
for fr in frames:
idx = imp.getStackIndex(ch,z,fr)
ip = imp.getStack().getProcessor(idx)
ip.setRoi(roi)
stats = ip.getStatistics()
ip.resetRoi()
stdDevs.append(stats.stdDev)
means.append(stats.mean)
return means, stdDevs
def getMedian(data):
return data[len(data)/2]
z=1
ch = channel_analyze
means, stdDevs = getMeanAndStdDev(imp0, inv_roi_cell, range(t_illum-dtWindow-1, t_illum-2), ch, z)
bgMean = getMedian(means)
bgStdDev = getMedian(stdDevs)
print "bg (mean):"+str(bgMean)
print "bg (stdDev):"+str(bgStdDev)
print ("tillum: "+str(t_illum))
scriptInfo = imp0.getProperty("scriptInfo")
scriptInfo['background mean'] = bgMean
scriptInfo['background stdDev'] = bgStdDev
scriptInfo['frameIllumination'] = t_illum
imp0.setProperty('scriptInfo', scriptInfo)
#roiManager.addRoi(inv_roi_cell)
#print bgValue
Copyright 2017 Max Planck Institute of Molecular Cell Biology and Genetics, Dresden, Germany
Redistribution and use in source and binary forms, with or without modification, are permitted provided that the following conditions are met:
1 -
Redistributions of source code must retain the above copyright notice, this list of conditions and the following disclaimer.
2 - Redistributions in binary form must reproduce the above copyright notice, this list of conditions and the following disclaimer in the documentation and/or other materials provided with the distribution.
3 - Neither the name of the copyright holder nor the names of its contributors may be used to endorse or promote products derived from this software without specific prior written permission.
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
\ No newline at end of file
# @ImagePlus (label="Sequence to analyze") impSeq0
# @int (label="channel to analyze") channel_analyze
# author: Benoit Lombardot @ Scientific Computing Facility, MPI-CBG, Dresden Germany
# License: BSD, see the LICENSE.txt file coming with the script
# 2016-11-14 : update to version 1.2, 3 empty columns are added to simplify further procesing
# 2016-12-13 : version 1.3, do measure from frame=1,
# add median0 (for each region median of the mean in the 5 frame before t_illum ) (robust version of mean0 )
# make the measure compatible with multiple channel image
# add background stddev to the measures
# add the measure of mean intensity in the parent cell of the region of measure
# 2016-12-19: version 1.4, force integer valued measure to show with one decimal in the results table
# it should facilitate data interpretation in R
# Also Rois of the illumination region, the increased inscreased signal mask, and the rois of the cell
# to analyze should be in the the roiManager.
# the illumination time is given by the illumination roi in the roi manager
# 2 results table are created:
# * cell information: 1 row per cell indicating
# * origin image name
# * analysis time
# * roi name
# * status (in/out/both)
# * (parameter used)
# * sequence length
# * t illumination
# * intensity 1st slice (for normalisation)
# * cell sequence: 2 column per cell (for the region of increased/decreased signal)
resultPrecision = 1 # results table number of decimal
from ij import IJ
from ij.plugin import Duplicator
from ij.plugin import ImageCalculator
from ij.plugin.frame import RoiManager
from ij.process import ByteProcessor
from ij import ImagePlus
#from ij.gui import ShapeRoi
from ij.plugin.filter import ThresholdToSelection
from ij.process import ImageProcessor
from ij.measure import ResultsTable
from ij import WindowManager;
from ij.text import TextWindow;
from ij.gui import Roi
import time, datetime
ch = channel_analyze
z=1
Dt = 5
scriptInfo = impSeq0.getProperty('scriptInfo')
#######################################################
# collect rois from the rois manager ##################
roiManager = RoiManager.getInstance()
if roiManager==None :
roiManager = RoiManager()
rois = roiManager.getRoisAsArray()
cellRois = []
for roi in rois:
rName = roi.getName()
if rName=="illuminationRegion":
roiIllum = roi
tIllum = int( roiIllum.getProperty("tIllum") )
#print tIllum
elif rName=="increasedSignalRegion":
roiIncrease = roi
elif not (rName=="cellMask"):
cellRois.append(roi)
##########################################################
# determine roi position versus the region of interest ###
# create a mask however issue can arise with create selection that select the black part of the image
# which is not necessarily 255
def createMask(roi,w,h,title):
ip = ByteProcessor( w, h )
ip.setValue(255)
ip.fill(roi)
return ImagePlus(title, ip)
def getMean(roi, imp, frame, ch, z):
idx = imp.getStackIndex(ch,z,frame)
ip = imp.getStack().getProcessor(idx)
ip.setRoi(roi)
mean = ip.getStatistics().mean
ip.resetRoi()
return mean
def getMedianOfTheMeans(roi, imp, frames, ch, z):
if len(frames)==0:
return -1
means = []
for frame in frames:
means.append( getMean(roi, imp, frame, ch, z) )
means = sorted(means)
idxMedian = len(means)/2
median = means[idxMedian]
return median
w = impSeq0.getWidth()
h = impSeq0.getHeight()
maskIllum = createMask(roiIllum, w, h, "mask Illum")
cellStatus = []
for roi in cellRois:
maskIllum.setRoi(roi)
stats = maskIllum.getStatistics()
if (stats.max==0):
cellStatus.append(0)
elif (stats.min==255):
cellStatus.append(1)
else:
cellStatus.append(2)
# reject the cell that are both inside and outside from the rest of the analysis
cellRois = [roi for roi,status in zip(cellRois, cellStatus) if status<2]
cellStatus = [status for status in cellStatus if status<2]
cellList = []
for i in range(len(cellRois)) :
cellRoi = cellRois[i]
cellData = {}
cellData["name"] = cellRoi.getName()
cellData["status"] = cellStatus[i]
cellData["roi"] = cellRoi
cellList.append(cellData)
##########################################################
# create roi of increased/decreased signal in each cell ##
def isEmpty(roi):
if roi==None:
return True
if roi.getType()==Roi.COMPOSITE:
if len(roi.getRois())>0 :
return False
else:
return True
else :
if roi.getPolygon().npoints >0 :
return False
else:
return True
measureRegions = []
maskIncrease = createMask(roiIncrease, w, h, "region of increased signal")
for cellData in cellList :
# create roiInc, roiDec
roi = cellData["roi"]
mean0Cell = getMean(roi, impSeq0, tIllum-1, ch, z)
x= roi.getBounds().x
y= roi.getBounds().y
maskIncrease.setRoi(roi)
auxIncMask = Duplicator().run(maskIncrease)
roi.setLocation(0,0)
auxCellMask = createMask(roi,roi.getBounds().width,roi.getBounds().height,"")
roi.setLocation(x,y)
cellIncMask = ImageCalculator().run('AND create ', auxCellMask, auxIncMask)
ip = cellIncMask.getProcessor()
ip.setThreshold(128,255,ImageProcessor.NO_LUT_UPDATE)
roiInc = ThresholdToSelection().convert(ip)
if not isEmpty(roiInc):
roiInc.setLocation(x+roiInc.getBounds().x, y+roiInc.getBounds().y)
roiName = roi.getName()+"_inc"
roiInc.setName( roiName )
mean0 = getMean(roiInc, impSeq0, tIllum-1, ch, z)
median0 = getMedianOfTheMeans(roiInc, impSeq0, range(tIllum-1-Dt, tIllum) , ch, z)
measureRegionData = { "roi": roiInc , \
"roiCell": cellData["roi"] , \
"isIncreased": 'increased' , \
"cellName": cellData["name"] , \
"cellStatus": cellData["status"] , \
"regionName": roiName , \
"mean0": mean0 , \
"mean0Cell": mean0Cell , \
"median0": median0 , \
}
measureRegions.append( measureRegionData )
IJ.run(auxIncMask, "Invert", "")
cellDecMask = ImageCalculator().run('AND create ', auxCellMask, auxIncMask)
ip = cellDecMask.getProcessor()
ip.setThreshold(128,255,ImageProcessor.NO_LUT_UPDATE)
roiDec = ThresholdToSelection().convert(ip)
if not isEmpty(roiDec):
roiDec.setLocation(x+roiDec.getBounds().x, y+roiDec.getBounds().y)
roiName = roi.getName()+"_dec"
roiDec.setName(roiName)
mean0 = getMean(roiDec, impSeq0, tIllum-1, ch, z)
median0 = getMedianOfTheMeans(roiDec, impSeq0, range(tIllum-1-Dt, tIllum) , ch, z)
measureRegionData = { "roi": roiDec , \
"roiCell": cellData["roi"] , \
"isIncreased": 'notIncreased' , \
"cellName": cellData["name"] , \
"cellStatus": cellData["status"] , \
"regionName": roiName , \
"mean0": mean0 , \
"mean0Cell": mean0Cell , \
"median0": median0 , \
}
measureRegions.append( measureRegionData )
############################################################
# measure reference signal in cell, band, measured region ##
##########################################################
# measure signal sequences for each cells ################
dims = impSeq0.getDimensions()
timeStamp = datetime.datetime.fromtimestamp(time.time()).strftime('%Y_%m_%d-%H_%M_%S') # if the file is analyzed multiple times
measures = []
timeValues = scriptInfo["timeValues"]
for fr in range(1,dims[4]+1):
idx = impSeq0.getStackIndex(ch,z,fr)
ip = impSeq0.getStack().getProcessor(idx)
for region in measureRegions :
roi = region["roi"]
ip.setRoi(roi)
stats = ip.getStatistics()
roi = region["roiCell"]