Revert thing to a working version

ImageJ/script_mosaic_telemos.py

@@ -29,23 +29,97 @@ from ij import WindowManager
 import re
 import sys
 
-# New way to import based on IJ documentation
-# https://imagej.net/scripting/jython/
-from java.lang.System import getProperty
-path.append(getProperty('user.home') + '/Documents/ImageJ-Script')
-sys.path.append(basedir)
+try:
+   basedir = os.path.dirname(__file__)
+except:
+    basedir = os.path.expanduser("~/Documents/ImageJ-Script")
+    print(basedir)
 
+sys.path.append(basedir)
 from utils.ImagesLoaders import load_image_sequence
 
 # Import to do operation with ImageJ2 and ImgMath
 from net.imglib2.util import Intervals
 from net.imglib2.img.display.imagej import ImageJFunctions as IL
-from net.imglib2.algorithm.math.ImgMath import computeIntoFloat, sub, div, GT, LT, IF, THEN, ELSE, minimum
+from net.imglib2.algorithm.math.ImgMath import computeIntoFloats, sub, div, GT, LT, IF, THEN, ELSE, minimum
 from net.imglib2.type.numeric.real import FloatType
 from net.imglib2.view import Views 
 from net.imagej import Dataset
 from net.imagej.axis import Axes
 
+def getFirstMetadata(metadata):
+    """
+    Return the first keys starting with Metadata
+    """
+    ml = list(metadata)
+    cpt = 0
+    for l in ml:
+        cpt += 1
+        if l.startswith('Metadata'):
+            break
+    
+    return cpt
+    
+def getMMversion(metadata):
+    """
+    Get MM version from metadata 
+    """
+
+    if int(metadata['Summary']["MicroManagerVersion"].split('.')[0]) >= 2:
+        version = 2
+    else:
+        version = 1
+
+    return version
+
+
+def getRoi(metadata, mm_version=2):
+    """
+    Return the ROI rectangle from metadata file
+    """
+
+    if mm_version == 2:
+        # Use the first image to get the ROI
+        i = getFirstMetadata(metadata)
+        roi = metadata[metadata.keys()[i]]['ROI']
+    else:
+        roi = metadata['Summary']['ROI']
+
+    roi = [int(r) for r in roi.split('-')]
+    return roi
+
+
+def getScale(metadata):
+    """
+    Return the scale from metadata of MM
+    """
+
+    try:
+        # MM2 version
+        i = getFirstMetadata(metadata)
+        scale=metadata[metadata.keys()[i]]['PixelSizeUm']
+    except:
+        # MM1 version
+        scale = metadata['Summary']['PixelSize_um']
+        
+    return scale
+
+def getWidthHeight(metadata):
+    """
+    Get the width and height from metadata
+    """
+    try:
+        # MM2 version
+        i = getFirstMetadata(metadata)
+        height = metadata[metadata.keys()[i]]['Height']
+        width = metadata[metadata.keys()[i]]['Width']
+    except:
+        # MM1 version
+        height = metadata['Summary']['Height']
+        width = metadata['Summary']['Width']
+        
+    return width, height
+    
 def load_image_tiles_stack(images_dir):
     """
     Fonction pour charge le stack d'images depuis
@@ -116,12 +190,20 @@ def load_image_tiles_stack(images_dir):
         file_pattern = selected_roi+'-Pos_{rows:003d}_{cols:003d}'
         mmetapath = os.path.join(base_dir,file_pattern.format(rows=0, cols=0),'metadata.txt')
 
+    # Chargement des métadonnées MM (juste le resumée de la premiere tuile suifit)
+    metadata = json.load(open(mmetapath))
+    mmversion = getMMversion(metadata)
+    
+    # Test de version sur MicroManager
+    if mmversion >= 2:
+        # Add MM2.0 file structure
+        file_pattern = os.path.join(file_pattern, 'img_channel{filter}_position{tile-1:003d}_time000000000_z{z:003d}.tif')
+        hstack, good_names = load_image_sequence(base_dir, file_pattern, return_filenames=True)
+    else:
     	# Add MM1.4 file structure
     	file_pattern = os.path.join(file_pattern, 'img_000000000_{filter}_{z:003d}.tif')
     	hstack, good_names = load_image_sequence(base_dir, file_pattern, return_filenames=True)
 
-    # Chargement des métadonnées MM (juste le resumée de la premiere tuile suifit)
-    metadata = json.load(open(mmetapath))['Summary']
     
     hstack.setTitle("Data")
 
@@ -129,13 +211,18 @@ def load_image_tiles_stack(images_dir):
 
 def load_images(images_dir, roi=None, virtual=False):
     d = os.path.abspath(images_dir)
-    # MM1.4 pattern for white
+    # MM1.4 pattern for white has a pos0
+    if os.path.isdir(os.path.join(d,'Pos0')):
         white_pattern = os.path.join('Pos0','img_000000000_{filter}_{z:003d}.tif')
-    imp = load_image_sequence(d, white_pattern)
+    else:
+        # Use MM2 pattern img_channel000_position000_time000000000_z000.tif
+        white_pattern = os.path.join('Default','img_channel{filter}_position000_time000000000_z{z:003d}.tif')
         
+    imp = load_image_sequence(d, white_pattern)
     if roi is not None and len(roi) == 4 and imp.width > roi[2] and imp.height > roi[3]:
+        IJ.log('Roi %s on image size (%i, %i)' % (str(roi[2])+','+str(roi[3]), imp.width, imp.height))
         imp.setRoi(*roi)
-       #IJ.run(img, "Crop", "")
+        #IJ.run(imp, "Crop", "")
 
     return imp
 
@@ -207,7 +294,7 @@ def process_white(implus_white, implus_darkofwhite, zofwhite, roi=None):
 
     # Substract white and dark
     opsus = sub(img_whiteC, img_darkC)
-    whitemindark = computeIntoFloat(opsus)
+    whitemindark = computeIntoFloats(opsus)
     maxwhite = float(ops.run("stats.max", whitemindark).toString())
 
     # Compute the minimum value of the ratio
@@ -219,7 +306,7 @@ def process_white(implus_white, implus_darkofwhite, zofwhite, roi=None):
             THEN(div(opsus,maxwhite)), 
             ELSE(ratiomin))
 
-    result = computeIntoFloat(op)
+    result = computeIntoFloats(op)
     # Convert back to ImagePlus
     img_white_corr = IL.wrap(result, 
                              "white corrected from dark (z: %i)" % zofwhite)
@@ -252,7 +339,7 @@ def process_tiles(implus_tiles_stack, channel, z, implus_white_corr, implus_dark
         The stack containing the dark images for each channels
 
     - metadata: dict,
-        The metadata of MM 1.4 (only summary)
+        The metadata of MM 1.4 or MM2
 
     return
     ------
@@ -265,11 +352,31 @@ def process_tiles(implus_tiles_stack, channel, z, implus_white_corr, implus_dark
     # Extract the dark for the color position of the dataset:
     # TODO: Fix selecting channel when dark and white does not
     # have the same amount of filter
-    dark = IL.wrap(implus_dark_stack)
+    dark = convertService.convert(implus_dark_stack, Dataset)
+    if dark.numDimensions() > 2:
         dark = Views.hyperSlice(dark, 2, channel-1)
+    else:
+        # Need to apply crop on dark
+        roi = getRoi(metadata)
+        if roi is not None and len(roi) == 4 and implus_dark_stack.width > roi[2] and implus_dark_stack.height > roi[3]:
+            xmin, ymin, xmax, ymax = roi
+            xmax = xmin + xmax 
+            ymax = ymin + ymax 
+        else:
+            xmin = ymin = 0
+            xmax = implus_dark_stack.width
+            ymax = implus_dark_stack.height
+            
+        
+        dark = ops.transform().crop(dark, 
+                              Intervals.createMinMax(xmin, ymin,
+                                                           xmax-1, ymax-1))
+                                                           
+        IJ.log('dark is not an hyperslice, only one dark, keep this one')
+
     
     # Prepare the white
-    white = IL.wrap(implus_white_corr)
+    white = convertService.convert(implus_white_corr, Dataset)
 
     # Extract the zposition and channel of the data-to-procced
     data = convertService.convert(implus_tiles_stack, Dataset)
@@ -296,22 +403,24 @@ def process_tiles(implus_tiles_stack, channel, z, implus_white_corr, implus_dark
     if inter is not None:
         data = ops.transform().crop(data, inter)
 
+
     # Operation
     op = [None] * Ntiles
     cpt = 0
 
     for tile in xrange(Ntiles):
         IJ.showProgress(tile, Ntiles)
+        IJ.log('Compute %i / %i' % (tile, Ntiles))
         data_views = Views.hyperSlice(data, 2, tile)
         opsus = sub(data_views, dark)
-        op[cpt] = IF(GT(opsus, 0), THEN(div(opsus, white)), ELSE(0))
+        op[cpt] = IF(GT(opsus, 0.0), THEN(div(opsus, white)), ELSE(0.0))
         cpt += 1
 
     result = datasetService.create(Views.stack([o.view(FloatType()) for o in op]))
     result = ops.convert().uint16(result)
 
     # Create the ImagePlus from the result converted back to uint16
-    name = metadata['ChNames'][channel-1]
+    name = metadata['Summary']['ChNames'][channel-1]
     imp_title = "Processed Tiles (channel: %i-%s, z: %i)" % (channel, name, z)
     img_tiles_corr = IL.wrap(result, imp_title)
     # Convert axis to Time and se the resolution
@@ -335,7 +444,7 @@ def run_processing(channel, zimg, zwhite, run_basic=True, optimizeXY=True, savep
     if not os.path.exists(savepath):
        os.makedirs(savepath)
 
-    chname = metadata['ChNames'][channel-1]
+    chname = metadata['Summary']['ChNames'][channel-1]
     # A folder for the processed tiles
     procfolder = savepath+'/DW_processed_tiles_%s_z%i_whitez%i' % (chname, zimg, zwhite)
     if not os.path.exists(procfolder):
@@ -401,8 +510,9 @@ def run_processing(channel, zimg, zwhite, run_basic=True, optimizeXY=True, savep
     proctiles.close()
 
 def run_gui(default_values=None):
-    channels = metadata['ChNames']
-    dataZmax = metadata['Slices']
+
+    channels = metadata['Summary']['ChNames']
+    dataZmax = metadata['Summary']['Slices']
     whiteZmax = white.getNSlices()
 
     if default_values is None:
@@ -475,30 +585,32 @@ def run_gui(default_values=None):
        white.close()
 
 
-def extract_tile_positions(metadata, selected_roi, tilesondisk, tilename='tile_'):
+def getMMtiles(metadata, selected_roi, tilesondisk):
     """
-    Extraction des positions des tuiles depuis les métadonnées
-    pour les rendres compatible avec le plugin de stitching.
+    Extract tiles position from micromanager metadata
+
+    tilesondisk, list of tiles name found on disk
     """
     
-    IJ.log("Extract tiles coordinates from metadata")
-    # Extract pixel size
-    pix2um = metadata['PixelSize_um']
-    # The roi width
-    roiwidth = metadata['Width']
+    mmversion = getMMversion(metadata)
+
+    mmsum = metadata['Summary']
+    xt = None
+    yt = None
     
+    if mmversion == 1:
         # Positions
-    xt = [None] * len(metadata['InitialPositionList'])
-    yt = [None] * len(metadata['InitialPositionList'])
+        xt = [None] * len(mmsum['InitialPositionList'])
+        yt = [None] * len(mmsum['InitialPositionList'])
         test_nameSE = 'roi'+selected_roi
         test_nameGC = '%s-Pos_' % selected_roi
 
-    uniquecols = set(a['GridColumnIndex'] for a in metadata['InitialPositionList'])
-    uniquerows = set(a['GridRowIndex'] for a in metadata['InitialPositionList'])
+        uniquecols = set(a['GridColumnIndex'] for a in mmsum['InitialPositionList'])
+        uniquerows = set(a['GridRowIndex'] for a in mmsum['InitialPositionList'])
         Nrow = max(uniquerows) + 1
         Ncol = max(uniquecols) + 1
 
-    for i, pl in enumerate(metadata['InitialPositionList']):
+        for i, pl in enumerate(mmsum['InitialPositionList']):
             if pl['Label'] in tilesondisk:
                 # print(pl['Label'])
                 if test_nameSE in pl['Label']:
@@ -518,11 +630,77 @@ def extract_tile_positions(metadata, selected_roi, tilesondisk, tilename='tile_'
             else:
                 IJ.log("[WARNING]: File %s listed on metadata does not exist on disk" % pl['Label'])
 
+    else:
+        # for MM2
+        xt = [None] * len(mmsum['StagePositions'])
+        yt = [None] * len(mmsum['StagePositions'])
+        test_nameSE = 'roi'+selected_roi
+        test_nameGC = '%s-Pos_' % selected_roi
+
+        uniquecols = set(a['GridCol'] for a in mmsum['StagePositions'])
+        uniquerows = set(a['GridRow'] for a in mmsum['StagePositions'])
+        Nrow = max(uniquerows) + 1
+        Ncol = max(uniquecols) + 1
+
+        for i, pl in enumerate(mmsum['StagePositions']):
+            if pl['Label'] in tilesondisk:
+                # print(pl['Label'])
+                if test_nameSE in pl['Label']:
+                    ind = i
+
+                if test_nameGC in pl['Label']:
+                    # Compute the globel indice ind = col*(Nrow)+row
+                    # as the order is not the same as the one loaded by bio-format reader
+                    c = pl['GridCol']
+                    r = pl['GridRow']
+                    ind = c * Nrow + r
+                    # print(ind, pl['Label'])
+
+                xt[ind] = pl['DevicePositions'][0]['Position_um'][0]
+                yt[ind] = pl['DevicePositions'][0]['Position_um'][1]
+
+            else:
+                IJ.log("[WARNING]: File %s listed on metadata does not exist on disk" % pl['Label'])
+
+    return xt, yt
+    
+def extract_tile_positions(metadata, selected_roi, tilesondisk, tilename='tile_'):
+    """
+    Extraction des positions des tuiles depuis les métadonnées
+    pour les rendres compatible avec le plugin de stitching.
+    """
+    
+    IJ.log("Extract tiles coordinates from metadata")
+    # Extract pixel size
+    pix2um = getScale(metadata)
+    # The roi width
+    roiwidth, _ = getWidthHeight(metadata)
+
+    # Is the microscope use Andor camera ?
+    i = getFirstMetadata(metadata)
+    if metadata[metadata.keys()[i]]['Camera'] == 'Andor':
+        ACAMAS = True
+    else:
+        ACAMAS = False
+        
+    xt, yt = getMMtiles(metadata, selected_roi, tilesondisk)
+    
+    if ACAMAS:
+        xtmp = yt
+        ytmp = xt
+        xt = xtmp
+        yt = ytmp
+        #xt = [-x for x in xtmp]
+        #yt = [-y for y in ytmp]        
+        
     # Remove None values in xt and yt (when the tile is not on disk)
     igood = [i for i,v in enumerate(xt) if v is not None]
     xt = [xt[i] for i in igood]
     yt = [yt[i] for i in igood]
 
+    if ACAMAS:
+        xt = [-x for x in xt]
+        
     xt = map(lambda x: x/pix2um, xt)
     yt = map(lambda x: x/pix2um, yt)
 
@@ -532,6 +710,7 @@ def extract_tile_positions(metadata, selected_roi, tilesondisk, tilename='tile_'
     xt = [x-xt[0] for x in xt]
     yt = [y-yt[0] for y in yt]
 
+        
     outstr = []
     for i in range(len(yt)):
         outstr += ['%s%04i.tif; ; (%0.2f, %0.2f)'%(tilename, i+1, yt[i], -xt[i])]
@@ -563,7 +742,11 @@ if __name__ in ['__builtin__','__main__']:
     hstack, metadata, rootdir, selectedroi, file_names = load_image_tiles_stack(fdata)
     hstack.show()
 
-    roi = metadata['ROI']
+	# GET the version of MM
+    mmversion = getMMversion(metadata)
+
+    # LOAD ROI info
+    roi = getRoi(metadata, mmversion)
     IJ.log('ROI: %s' % str(roi))
     
     # Make a list of tiles names