Merge pull request #76 from WeisongZhao/master

Temporal color-code and bleach correction

Author: yxdragon

imagepy/core/util/fileio.py

@@ -8,7 +8,8 @@
 def show_img(img, title):
     if isinstance(img, list):
         return IPy.show_img(img, title)
-    if img.dtype!=np.uint8 and img.ndim>2 and img.shape[2]!=3:
+    # if img.dtype!=np.uint8 and img.ndim>2 and img.shape[2]!=3:
+    if img.ndim>2 and img.shape[2]!=3:
         return IPy.show_img(img, title)
     if img.dtype==np.uint8 and img.ndim==3 and img.shape[2]==4:
         img = img[:,:,:3].copy()

imagepy/menus/Image/Adjust/__init__.py

@@ -1 +1 @@
-catlog = ['threshold_plg','graystairs_plg','brightcons_plg','curve_plg','-','colorbalance_plg','colorstairs_plg','-','histogram_plgs']
+catlog = ['threshold_plg','graystairs_plg','brightcons_plg','curve_plg','-','bleachCorrection_plg','enhanceContrast_plg','normalize_plg','-','colorbalance_plg','colorstairs_plg','-','histogram_plgs']

imagepy/menus/Image/Adjust/bleachCorrection_plg.py

@@ -0,0 +1,91 @@
+"""
+Created on Sun Jan 23 11:53:00 2020
+@author: weisong
+"""
+from imagepy.core.engine import Simple
+import numpy as np
+from imagepy import IPy
+from scipy.optimize import curve_fit
+from skimage.exposure import histogram_matching
+import matplotlib.pyplot as plt 
+import pandas as pd
+
+def copy(imgs):
+    if isinstance(imgs, list):
+        return [np.zeros_like(imgs[0])]
+    else: return np.zeros_like(imgs)
+
+def exponential_func(t, ref, k, offset):
+    return ref * np.exp(- k * t) + offset
+
+def simple_ratio(imgs, back=0, inplace=True, out=print):
+    if isinstance(back, int): back=imgs[back]
+    buf = imgs if inplace else copy(imgs)
+    z, (x, y) = len(imgs), imgs[0].shape 
+    values, k0 = np.zeros(z), back.sum()/x/y
+    lim = 255 if imgs[0].dtype.type==np.uint8 else 65535
+    for i in range(z):
+        values[i] = imgs[i].sum()/x/y 
+        np.clip(imgs[i], 0, lim/(k0/values[i]), out=buf[i])
+        np.multiply(buf[i], k0/values[i], out=buf[i], casting='unsafe')
+        out(i, z)
+    return buf, values, k0/values
+
+def exponential_fit(imgs, inplace=True, out=print):
+    buf = imgs if inplace else copy(imgs)
+    z, (x, y) = len(imgs), imgs[0].shape 
+    intensity = [i.sum()/x/y for i in imgs]
+    popt, pcov = curve_fit(exponential_func, np.arange(z), intensity)
+    k0 = exponential_func(0, popt[0], popt[1], popt[2])
+    rst = exponential_func(np.arange(z), popt[0], popt[1], popt[2])
+    lim = 255 if imgs[0].dtype.type==np.uint8 else 65535
+    for i in range(z):
+        np.clip(imgs[i], 0, lim/(rst[i]/k0), out=buf[i])
+        np.multiply(buf[i], rst[i]/k0, out=buf[i], casting='unsafe')
+        out(i, z)
+    return buf, popt, intensity, rst
+
+def histogram_match(imgs, back=0, inplace=True, out=print): 
+    if isinstance(back, int): back=imgs[back]   
+    buf = imgs if inplace else copy(imgs)
+    z, (x, y) = len(imgs), imgs[0].shape 
+    for i in range(z):
+        buf[i] = histogram_matching.match_histograms(imgs[i], back)
+        out(i, z)
+    return buf
+
+def plot(popt, intensity, fitresult):
+    t = np.arange(len(intensity))
+    plt.plot(t, intensity,'r.',label='Experiment')
+    plt.plot(t, fitresult,'k',label=
+        'Exponential fitted curve\n y=a*exp(-bx)+c\n a=%f\n b=%f\n c=%f'%tuple(popt))
+    plt.title('Exponential fitted result')
+    plt.legend()
+    plt.show()
+
+def plot_after(popt, intensity, fitresult):
+    import wx
+    wx.CallAfter(plot, popt, intensity, fitresult)
+
+class Plugin(Simple):
+    title = 'Bleach Correction'
+    note = ['8-bit','16-bit','stack']
+    para = {'method':'Simple Ratio', 'new':True}
+    view = [(list, 'method', ['Simple Ratio','Exponential Fit','Histogram Match'], 
+        str, 'Correction Method',''),
+            (bool, 'new', 'show new window'),
+            ('lab', 'lab', 'Correct intensity based on your current slice!')]
+            
+    def run(self, ips, imgs, para = None):
+        if para['method'] == 'Simple Ratio':
+            rst, value, ratio = simple_ratio(imgs, ips.img, not para['new'], self.progress)
+            body = pd.DataFrame({'Mean value': value, 'Ratio':ratio})
+            IPy.show_table(body, '%s-simple ratio'%ips.title)
+        if para['method'] == 'Exponential Fit':
+            rst, popt, intensity, fitrst = exponential_fit(imgs, not para['new'], self.progress)
+            plot_after(popt, intensity, fitrst)
+            body = {'Intensity':intensity, 'Fit':fitrst}
+            IPy.show_table(pd.DataFrame(body), '%s-exp fit'%ips.title)
+        if para['method'] == 'Histogram Match':
+            rst = histogram_match(imgs, ips.img, not para['new'], self.progress)
+        if para['new']:  IPy.show_img(rst, '%s-corrected'%ips.title)

imagepy/menus/Image/Adjust/enhanceContrast_plg.py

@@ -0,0 +1,17 @@
+"""
+Created on Sun Jan 25 17:00:00 2020
+@author: weisong
+"""
+from imagepy.core.engine import Filter
+from skimage import exposure
+import numpy as np
+
+class Plugin(Filter):
+    title = 'Enhance contrast'
+    note = ['all', 'auto_msk', 'auto_snap','preview']
+    para = {'percentage': 0.3}
+    view = [(float, 'percentage', (0,100), 4, 'Saturated pixels', '%')]
+    
+    def run(self, ips, snap, img, para = None):
+        up, down = np.percentile(snap, (0, 100 - para['percentage']))
+        return exposure.rescale_intensity(snap, in_range=(up, down))

imagepy/menus/Image/Adjust/normalize_plg.py

@@ -0,0 +1,33 @@
+"""
+Created on Sun Jan 25 9:00:00 2020
+@author: weisong
+"""
+from imagepy.core.engine import Simple
+import numpy as np
+from imagepy import IPy
+
+class Plugin(Simple):
+    title = 'Normalize'
+    note = ['8-bit','16-bit','float','stack']
+    para = {'if3d': False, 'Sb':False}
+    view = [(bool, 'if3d', '3D stack'),
+    (bool, 'Sb', 'Subtract background')]
+
+    def run(self, ips, imgs, para = None):
+        imgs_= np.array(imgs).astype('float64').transpose(1,2,0)
+        if para['if3d']:
+            if para['Sb']:
+                imgs_ -= imgs_.min()
+            imgs_ = imgs_ / imgs_.max()
+        else:
+            if para['Sb']:
+                for i in range(z):
+                    imgs_[:,:,i] -= imgs_[:,:,i].min()
+            for i in range(z):
+                imgs_[:,:,i] = imgs_[:,:,i] / imgs_[:,:,i].max()        
+        if imgs.dtype == np.uint8:
+            imgs[:] = (255*imgs_).astype(imgs.dtype).transpose(2,0,1)
+        elif imgs.dtype == np.uint16:
+            imgs[:] = (65535*imgs_).astype(imgs.dtype).transpose(2,0,1)
+        else:
+            imgs[:] = (imgs_).astype(imgs.dtype).transpose(2,0,1)

imagepy/menus/Image/duplicate_plg.py

@@ -14,12 +14,18 @@ class Duplicate(Simple):
     title = 'Duplicate'
     note = ['all']
 
-    para = {'name':'Undefined','stack':True}
+    para = {'name':'Undefined','start':1,'end':2,'stack':True}
 
     def load(self, ips):
+        self.slength = len(ips.imgs)
         self.para['name'] = ips.title+'-copy'
+        self.para['end'] = self.slength
         self.view = [(str, 'name', 'Name', '')]
         if ips.get_nslices()>1:
+            self.view.append((int, 'start', 
+                (1,self.slength),0,'Start slice','1~%d'%self.slength))
+            self.view.append((int, 'end', 
+                (1,self.slength),0,'End slice','1~%d'%self.slength))
             self.view.append((bool, 'stack', 'duplicate stack'))
         return True
     #process
@@ -42,6 +48,7 @@ def run(self, ips, imgs, para = None):
                     ipsd.backimg = ips.backimg[sr, sc]
                 '''
         elif ips.get_nslices()>1 and self.para['stack']:
+            imgs=imgs[para['start']-1: para['end']]
             if ips.roi == None:
                 if ips.is3d:imgs=imgs.copy()
                 else:imgs = [i.copy() for i in imgs]

imagepy/menus/Plugins/__init__.py

@@ -1 +1 @@
-catlog = ['New', 'Macros', 'Manager', '-', 'Install', 'Contribute', 'update_plg', '-', 'StackReg', 'Games', 'screencap_plg']
+catlog = ['New', 'Macros', 'Manager', '-', 'Install', 'Contribute', 'update_plg', '-', 'temporal_plg','StackReg', 'Games', 'screencap_plg']

imagepy/menus/Plugins/temporal_plg.py

@@ -0,0 +1,42 @@
+"""
+Created on Sun Jan 22 12:56:00 2020
+@author: weisong
+"""
+from imagepy.core.engine import Simple
+import numpy as np
+from imagepy.core.manager import ColorManager
+from imagepy import IPy
+
+def color_code(img, lut):
+    idx = np.linspace(0,255,len(img)).astype(int)
+    cs = lut[idx].astype(np.uint32)
+    buf = np.zeros(img[0].shape+(3,), dtype=np.uint32)
+    for im, c in zip(img, cs):
+        buf += im.reshape(im.shape+(-1,)) * c
+    k = 255/buf.max(axis=(0,1)).reshape((1,1,3))
+    return (buf * k).astype(np.uint8)
+
+class Plugin(Simple):
+    title = 'Temporal color-code'
+    note = ['all', 'stack']
+    para = {'LUT':'Jet',
+            'Start image':1,
+            'End image': 2,
+            'Creatbar':True}
+    def load(self, ips): 
+        self.slength = len(ips.imgs)
+        self.para['End image'] = self.slength
+        self.view = [(list, 'LUT', list(ColorManager.luts.keys()), str, 'LUT',''),
+            (int, 'Start image', (1,self.slength),0,'Start image','1~%d'%self.slength),
+            (int, 'End image', (2,self.slength),0,'End image','start~%d'%self.slength),
+            (bool, 'Creatbar', 'Creat time color scale bar')]
+        return True
+
+    def run(self, ips, imgs, para = None):
+        cmap = ColorManager.luts[ para['LUT'] ]
+        imglut = color_code(imgs[para['Start image']-1: para['End image']], cmap)
+        IPy.show_img([imglut],'Color-coded %s'%ips.title)
+        if para['Creatbar']:
+            cmapshow = np.ones([32,256,3])*cmap
+            IPy.show_img([cmapshow.astype('uint8')],'Color bar')
+