Setup


In [1]:
%matplotlib inline


/usr/local/python-2.7.6/lib/python2.7/site-packages/matplotlib-1.5.1-py2.7-linux-x86_64.egg/matplotlib/font_manager.py:273: UserWarning: Matplotlib is building the font cache using fc-list. This may take a moment.
  warnings.warn('Matplotlib is building the font cache using fc-list. This may take a moment.')

In [2]:
from numpy import array
import matplotlib.pyplot as plt
import seaborn as sns


/usr/local/python-2.7.6/lib/python2.7/site-packages/matplotlib-1.5.1-py2.7-linux-x86_64.egg/matplotlib/__init__.py:872: UserWarning: axes.color_cycle is deprecated and replaced with axes.prop_cycle; please use the latter.
  warnings.warn(self.msg_depr % (key, alt_key))

In [3]:
import thunder
from showit import image, tile
import matplotlib.animation as animation

In [4]:
from os.path import join, exists
from os import mkdir, makedirs

In [5]:
from numpy import save

In [6]:
from skimage.io import imsave, imread

In [7]:
from regional import many
from numpy import random

In [8]:
def normalize(oim):
    # normalizes 3D image across first axis
    assert oim.ndim == 3
    
    means = oim.mean(axis=(1, 2), dtype='float32')
    return array([oim[i]/means[i]/4 for i in range(oim.shape[0])]).clip(0, 1)

Set directory and session information


In [9]:
directory = '/tier2/freeman/Nick/lfov.calibration'

In [10]:
key = '2016-04-16-L23'
name = 'anm-0330549'

In [11]:
path = join(directory, 'reprocessed', name, key)
print exists(path)


True

Create savepath


In [12]:
savepath = join(path, 'sources')
if not exists(savepath):
    makedirs(savepath)

Load data


In [15]:
data = thunder.images.frombinary(join(path, 'patch'), engine=sc, npartitions = 1000)

In [21]:
data = data[:,:200,:200]

In [23]:
data.compute()

In [16]:
from numpy import arange, polyfit, polyval

In [24]:
def detrend(y, order=5):
        """
        Detrend series data with linear or nonlinear detrending.
        Preserve intercept so that subsequent operations can adjust the baseline.
        Parameters
        ----------
        method : str, optional, default = 'linear'
            Detrending method
        order : int, optional, default = 5
            Order of polynomial, for non-linear detrending only
        """

        x = arange(len(y))
        p = polyfit(x, y, order)
        p[-1] = 0
        yy = polyval(p, x)
        return y - yy

In [25]:
detrended = data.map_as_series(lambda x: detrend(x))

In [26]:
detrended.cache()
detrended.compute()

In [29]:
#detrended.tobinary(join(path, 'patch'))

In [30]:
detrended


Out[30]:
Images
mode: spark
dtype: float64
shape: (8279, 200, 200)

In [32]:
mean = detrended.mean().toarray()

In [34]:
image(mean, clim=(0, 3.5*mean.mean()), size=12)


Out[34]:
<matplotlib.image.AxesImage at 0x5d9dc90>

In [57]:
meanG = detrended.uniform_filter(2).mean().toarray()

In [58]:
image(meanG, clim=(0, 3.5*meanG.mean()), size=12)


Out[58]:
<matplotlib.image.AxesImage at 0x73f5450>

compute local corr


In [16]:
localcorr2 = data.localcorr(2).astype('float32')

In [18]:
image(localcorr2, clim=(0, 2*localcorr2.mean()), size=12)


Out[18]:
<matplotlib.image.AxesImage at 0x6b57590>

In [ ]:


In [35]:
localcorr4 = data.localcorr(4).astype('float32')

In [36]:
image(localcorr4, clim=(0, 3.5*localcorr4.mean()), size=12)


Out[36]:
<matplotlib.image.AxesImage at 0x6abae10>

In [38]:
localcorr8 = detrended.localcorr(8).astype('float32')

In [39]:
image(localcorr8, clim=(0, 3.5*localcorr8.mean()), size=12)


Out[39]:
<matplotlib.image.AxesImage at 0x73e8a10>

In [59]:
localcorr1G8 = detrended.uniform_filter(2).localcorr(8).astype('float32')

In [60]:
image(localcorr1G8, clim=(0, 3.5*localcorr1G8.mean()), size=12)


Out[60]:
<matplotlib.image.AxesImage at 0x6a8ee90>

In [52]:
localcorr16 = detrended.localcorr(16).astype('float32')

In [53]:
image(localcorr16, clim=(0, 3.5*localcorr16.mean()), size=12)


Out[53]:
<matplotlib.image.AxesImage at 0x82c79d0>

In [ ]:


In [ ]:


In [ ]:


In [ ]:


In [ ]:


In [ ]:


In [17]:
data = imread(join(path, 'summary', 'localcorr.tif'), plugin='tifffile')
#data = data.astype('float32')/255

In [16]:
#norm = normalize(data)

In [32]:
img = data[:200,:200]

In [ ]:


In [33]:
image(img, size=12)


Out[33]:
<matplotlib.image.AxesImage at 0x10f04f850>

Blob finding


In [19]:
from skimage.filters import gaussian_filter
from skimage.feature import peak_local_max

In [20]:
from scipy.optimize import curve_fit
from scipy.stats import multivariate_normal

In [21]:
def inside(im, center, window):
    inside = [(center[i]-window[i] >= 0) & (center[i]+window[i] <= im.shape[i]) for i in range(0,2)]
    inside = np.all(inside)
    return inside

In [22]:
def volume(im, center, window):
    if np.all(inside(im, center, window)):
        volume = im[(center[0]-window[0]):(center[0]+window[0]), (center[1]-window[1]):(center[1]+window[1])]
        volume = volume.astype('float64')
        baseline = volume[[0,-1],[0,-1]].mean()
        volume = volume - baseline
        volume = volume/volume.max()
        return volume

In [24]:
imgT = thunder.images.fromarray(localcorr2)

In [25]:
smoothed = imgT.uniform_filter(4)

In [26]:
image(smoothed, size=12)


Out[26]:
<matplotlib.image.AxesImage at 0x74761d0>

In [27]:
import numpy as np

In [28]:
threshold_abs = np.percentile(smoothed.toarray().flatten(),70)
print threshold_abs
threshold_rel = np.percentile(smoothed.toarray().flatten(),80) - threshold_abs
print threshold_rel


0.535475963354
0.00891479849815

In [29]:
from skimage.feature import blob_log, blob_doh, blob_dog

In [31]:
centers = blob_log(localcorr2, min_sigma=.2, max_sigma=4, num_sigma=20, threshold=.05)
len(centers)


Out[31]:
145

In [215]:
centers = blob_dog(localcorr16, min_sigma=.3, max_sigma=3, sigma_ratio=1.05, threshold=.01)
len(centers)


Out[215]:
66

In [140]:
def findBeads(img):
    centers = peak_local_max(img, min_distance=1, threshold_rel=threshold_rel, threshold_abs=threshold_abs, exclude_border=False)
    return centers

In [141]:
centers = findBeads(smoothed.toarray())
len(centers)


Out[141]:
418

In [34]:
image(localcorr2, size=12)
plt.plot(centers[:,1], centers[:,0], 'g.', ms=20);
plt.xlim([0, localcorr2.shape[1]])
plt.ylim([localcorr2.shape[0], 0])


Out[34]:
(200, 0)

In [219]:
[centers[20,1], centers[20,0]]


Out[219]:
[33, 50]

In [239]:
blurred = detrended.uniform_filter(16)

In [242]:
trace = blurred[:,50,33].toarray()

In [243]:
plt.plot(trace)


Out[243]:
[<matplotlib.lines.Line2D at 0xcc40ad0>]

In [232]:
records = detrended.toseries()

In [245]:
corIndvB = records.correlate(trace).toarray()

In [244]:
corIndv = X

In [240]:
image(corIndv, size=12)


Out[240]:
<matplotlib.image.AxesImage at 0x92e7290>

In [246]:
image(corIndvB, size=12)


Out[246]:
<matplotlib.image.AxesImage at 0x16aefa10>

In [ ]:


In [ ]:


In [ ]:


In [48]:
def getCenters(im, options):
    window = [options['windowUm'][0]*options['pxPerUmAx'], options['windowUm'][1]*options['pxPerUmLat']]
    window = [round(x) for x in window]
    centers, smoothed = findBeads(im, window)
    print len(centers)
    beads = [volume(im, x, window) for x in centers]
    maxima = [im[x[0], x[1]] for x in centers]
    return beads, maxima, centers, smoothed

In [49]:
FOVumLat = 590.0
FOVpxLat = 512.0 # 512
pxPerUmLat = FOVpxLat/FOVumLat
pxPerUmAx = .05 # 2.0
wavelength = 970.0
NA = 0.6
windowUm = [2, 2]
options = {'FOVumLat':FOVumLat, 'FOVpxLat':FOVpxLat, 'pxPerUmLat':FOVpxLat/FOVumLat, 'pxPerUmAx':pxPerUmAx, 'wavelength':970.0, 'NA':0.6, 'windowUm':windowUm}

In [50]:
import numpy as np

In [51]:
beads, maxima, centers, smoothed = getCenters(patch, options)


942
/Users/sofroniewn/anaconda/lib/python2.7/site-packages/IPython/kernel/__main__.py:3: VisibleDeprecationWarning: using a non-integer number instead of an integer will result in an error in the future
  app.launch_new_instance()
---------------------------------------------------------------------------
IndexError                                Traceback (most recent call last)
<ipython-input-51-656fda8f7292> in <module>()
----> 1 beads, maxima, centers, smoothed = getCenters(patch, options)

<ipython-input-48-a9ac1256a193> in getCenters(im, options)
      4     centers, smoothed = findBeads(im, window)
      5     print len(centers)
----> 6     beads = [volume(im, x, window) for x in centers]
      7     maxima = [im[x[0], x[1]] for x in centers]
      8     return beads, maxima, centers, smoothed

<ipython-input-47-67fb59aa567f> in volume(im, center, window)
      3         volume = im[(center[0]-window[0]):(center[0]+window[0]), (center[1]-window[1]):(center[1]+window[1])]
      4         volume = volume.astype('float64')
----> 5         baseline = volume[[0,-1],[0,-1]].mean()
      6         volume = volume - baseline
      7         volume = volume/volume.max()

IndexError: index 0 is out of bounds for axis 0 with size 0

In [ ]:


In [ ]:


In [18]:
img = data

In [ ]:


In [ ]:


In [ ]:


In [ ]:

Select plane (if necessary)


In [56]:
plane = 0
img = norm[plane]

Load sources


In [205]:
with open(join(path, 'sources', 'sources-%04d.tif' % plane)) as fid:
    sources = many([x['coordinates'] for x in json.load(fid)])


---------------------------------------------------------------------------
IOError                                   Traceback (most recent call last)
<ipython-input-205-5a0ac918c574> in <module>()
----> 1 with open(join(path, 'sources', 'sources-%04d.tif' % plane)) as fid:
      2     sources = many([x['coordinates'] for x in json.load(fid)])

IOError: [Errno 2] No such file or directory: '/volumes/freeman/Nick/lfov.calibration/reprocessed/anm-0326049/2016-04-01/sources/sources-0000.tif'

Show image


In [19]:
viz = lgn.imagepoly(img)
viz


---------------------------------------------------------------------------
ConnectionError                           Traceback (most recent call last)
<ipython-input-19-83d7fe3aa1e0> in <module>()
----> 1 viz = lgn.imagepoly(img)
      2 viz

<string> in <lambda>(self, imagedata, coordinates, xy)

/groups/svoboda/home/sofroniewn/lightning-python/lightning/types/decorators.pyc in plotter(self, *args, **kwargs)
     10     def plotter(self, *args, **kwargs):
     11         if not hasattr(self, 'session'):
---> 12             self.create_session()
     13         viz = VizType.baseplot(self.session, VizType._name, *args, **kwargs)
     14         self.session.visualizations.append(viz)

/groups/svoboda/home/sofroniewn/lightning-python/lightning/main.pyc in create_session(self, name)
     83         will be "Session No." with the number automatically generated.
     84         """
---> 85         self.session = Session.create(self.host, name=name, auth=self.auth)
     86         return self.session
     87 

/groups/svoboda/home/sofroniewn/lightning-python/lightning/session.pyc in create(cls, host, name, auth)
     46         headers = {'Content-type': 'application/json', 'Accept': 'text/plain'}
     47 
---> 48         r = requests.post(url, data=json.dumps(payload), headers=headers, auth=auth)
     49         return cls(host=host, json=r.json(), auth=auth)

/usr/local/python-2.7.6/lib/python2.7/site-packages/requests-2.9.1-py2.7.egg/requests/api.pyc in post(url, data, json, **kwargs)
    105     """
    106 
--> 107     return request('post', url, data=data, json=json, **kwargs)
    108 
    109 

/usr/local/python-2.7.6/lib/python2.7/site-packages/requests-2.9.1-py2.7.egg/requests/api.pyc in request(method, url, **kwargs)
     51     # cases, and look like a memory leak in others.
     52     with sessions.Session() as session:
---> 53         return session.request(method=method, url=url, **kwargs)
     54 
     55 

/usr/local/python-2.7.6/lib/python2.7/site-packages/requests-2.9.1-py2.7.egg/requests/sessions.pyc in request(self, method, url, params, data, headers, cookies, files, auth, timeout, allow_redirects, proxies, hooks, stream, verify, cert, json)
    466         }
    467         send_kwargs.update(settings)
--> 468         resp = self.send(prep, **send_kwargs)
    469 
    470         return resp

/usr/local/python-2.7.6/lib/python2.7/site-packages/requests-2.9.1-py2.7.egg/requests/sessions.pyc in send(self, request, **kwargs)
    574 
    575         # Send the request
--> 576         r = adapter.send(request, **kwargs)
    577 
    578         # Total elapsed time of the request (approximately)

/usr/local/python-2.7.6/lib/python2.7/site-packages/requests-2.9.1-py2.7.egg/requests/adapters.pyc in send(self, request, stream, timeout, verify, cert, proxies)
    435                 raise RetryError(e, request=request)
    436 
--> 437             raise ConnectionError(e, request=request)
    438 
    439         except ClosedPoolError as e:

ConnectionError: HTTPConnectionPool(host='kafka1.int.janelia.org', port=3000): Max retries exceeded with url: /sessions/ (Caused by NewConnectionError('<requests.packages.urllib3.connection.HTTPConnection object at 0x4ca0190>: Failed to establish a new connection: [Errno 111] Connection refused',))

In [ ]:
X = viz.points(z=plane)

In [ ]:
X

In [ ]:


In [58]:
sources = many(viz.points(z=plane))

Show sources


In [65]:
from numpy import maximum, tile

In [77]:
base.shape


Out[77]:
(512, 512, 4, 3)

In [74]:
sources[0]


Out[74]:
region
center: [254.6043956043956, 282.5824175824176, 0.0]
bbox: [[250, 259], [277, 288], [0, 0]]

In [78]:
base = tile(norm,(3,1,1,1)).transpose(1,2,3,0).transpose(1,2,0,3)
masks = sources.mask((512,512,4), background='black', fill='blue', stroke='orange')
blend = maximum(base, masks)

In [81]:
fig = plt.figure(figsize=[12,12])
ax = plt.axes()
im = image(blend[:,:,0], ax=ax)
plt.xlim([0, blend.shape[1]]);
plt.ylim([blend.shape[0], 0]);
#for s in range(sources.count):
#    plt.annotate(s=str(s), xy=[sources.center[s][1],sources.center[s][0]], color='w', size = 20);



In [82]:
imsave(join(path, 'sources', 'sources.tif'), (255*blend).astype('uint8'), plugin='tifffile', photometric='rgb')
#imsave(join(path, 'sources', 'sources.tif'), (255*blend).astype('uint8'), plugin='tifffile', photometric='rgb')


---------------------------------------------------------------------------
IOError                                   Traceback (most recent call last)
<ipython-input-82-b30365c4f171> in <module>()
----> 1 imsave(join(path, 'sources', 'sources.tif'), (255*blend).astype('uint8'), plugin='tifffile', photometric='rgb')
      2 #imsave(join(path, 'sources', 'sources.tif'), (255*blend).astype('uint8'), plugin='tifffile', photometric='rgb')

/Users/sofroniewn/anaconda/lib/python2.7/site-packages/skimage/io/_io.pyc in imsave(fname, arr, plugin, **plugin_args)
    158     if is_low_contrast(arr):
    159         warnings.warn('%s is a low contrast image' % fname)
--> 160     return call_plugin('imsave', fname, arr, plugin=plugin, **plugin_args)
    161 
    162 

/Users/sofroniewn/anaconda/lib/python2.7/site-packages/skimage/io/manage_plugins.pyc in call_plugin(kind, *args, **kwargs)
    205                                (plugin, kind))
    206 
--> 207     return func(*args, **kwargs)
    208 
    209 

/Users/sofroniewn/anaconda/lib/python2.7/site-packages/tifffile/tifffile.pyc in imsave(filename, data, **kwargs)
    283 
    284     with TiffWriter(filename, **tifargs) as tif:
--> 285         tif.save(data, **kwargs)
    286 
    287 

/Users/sofroniewn/anaconda/lib/python2.7/site-packages/tifffile/tifffile.pyc in __exit__(self, exc_type, exc_value, traceback)
   1057 
   1058     def __exit__(self, exc_type, exc_value, traceback):
-> 1059         self.close()
   1060 
   1061 

/Users/sofroniewn/anaconda/lib/python2.7/site-packages/tifffile/tifffile.pyc in close(self, truncate)
   1049         """Write remaining pages (if not truncate) and close file handle."""
   1050         if not truncate:
-> 1051             self._write_remaining_pages()
   1052         self._write_image_description()
   1053         self._fh.close()

/Users/sofroniewn/anaconda/lib/python2.7/site-packages/tifffile/tifffile.pyc in _write_remaining_pages(self)
    973             tag_offset = fh.tell()
    974             fh.write(tag_bytes)
--> 975             self._ifd_offset = fh.tell()
    976             fh.write(pack(offset_format, 0))  # offset to next IFD
    977 

IOError: [Errno 57] Socket is not connected

In [203]:
imsave(join(path, 'sources', 'sources-%04d.tif' % plane), (255*blend).astype('uint8'), plugin='tifffile', photometric='rgb')


---------------------------------------------------------------------------
IOError                                   Traceback (most recent call last)
<ipython-input-203-642bd4be944b> in <module>()
----> 1 imsave(join(path, 'sources', 'sources-%04d.tif' % plane), (255*blend).astype('uint8'), plugin='tifffile', photometric='rgb')

/Users/sofroniewn/anaconda/lib/python2.7/site-packages/skimage/io/_io.pyc in imsave(fname, arr, plugin, **plugin_args)
    158     if is_low_contrast(arr):
    159         warnings.warn('%s is a low contrast image' % fname)
--> 160     return call_plugin('imsave', fname, arr, plugin=plugin, **plugin_args)
    161 
    162 

/Users/sofroniewn/anaconda/lib/python2.7/site-packages/skimage/io/manage_plugins.pyc in call_plugin(kind, *args, **kwargs)
    205                                (plugin, kind))
    206 
--> 207     return func(*args, **kwargs)
    208 
    209 

/Users/sofroniewn/anaconda/lib/python2.7/site-packages/tifffile/tifffile.pyc in imsave(filename, data, **kwargs)
    282         tifargs['bigtiff'] = True
    283 
--> 284     with TiffWriter(filename, **tifargs) as tif:
    285         tif.save(data, **kwargs)
    286 

/Users/sofroniewn/anaconda/lib/python2.7/site-packages/tifffile/tifffile.pyc in __init__(self, filename, bigtiff, byteorder, software, imagej)
    369         self._tag_offsets = None  # strip or tile offset tag code
    370 
--> 371         self._fh = open(filename, 'wb')
    372         self._fh.write({'<': b'II', '>': b'MM'}[byteorder])
    373 

IOError: [Errno 2] No such file or directory: '/volumes/freeman/Nick/lfov.calibration/reprocessed/anm-0326049/2016-04-01/sources/sources-0000.tif'

In [207]:
###multi dim nature of this ...

In [ ]:
foo = sources.masks((512,512,4), color='blue', base=X)
image(foo[:,:,1,:])

In [227]:
from numpy import concatenate

In [231]:
x = many([[concatenate((x, [plane])) for x in roi.coordinates] for roi in sources])

In [235]:
if norm.ndim == 2:
    base = tile(norm,(3,1,1)).transpose(1,2,0)
    masks = sources.mask(norm.shape, background='black', fill='blue', stroke='orange')
    blend = maximum(base, masks)
else:
    base = tile(norm,(3,1,1,1)).transpose(1,2,3,0)
    masks = [sources.mask(norm.shape, background='black', fill='blue', stroke='orange')


---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-235-4058578af55b> in <module>()
      4     base = tile(norm,(3,1,1,1)).transpose(1,2,3,0)
      5 
----> 6 masks = x.mask(norm.shape, background='black', fill='blue', stroke='orange')
      7 blend = maximum(base, masks)

/Users/sofroniewn/anaconda/lib/python2.7/site-packages/regional/regional.pyc in mask(self, dims, base, fill, stroke, background, cmap, value)
    410         minbound = asarray([b[0:2] for b in self.bbox]).min(axis=0)
    411         maxbound = asarray([b[2:] for b in self.bbox]).max(axis=0)
--> 412         extent = maxbound - minbound + 1
    413 
    414         if dims is None and base is None:

ValueError: operands could not be broadcast together with shapes (4,) (2,) 

In [ ]:
def mask(sources, shape, **kwargs):
    if len(shape) == 3:
        
    else:
        return sources.mask(shape, **kwargs)

In [1]:
import regional

In [2]:
regional.__version__


Out[2]:
'1.0.1'

In [19]:
a = [1, 18, 35]
b = [2, 20, 40]

In [5]:
import numpy as np

In [20]:
x = np.stack((a, b), axis = 1)

In [21]:
x


Out[21]:
array([[ 1,  2],
       [18, 20],
       [35, 40]])

In [23]:
x[:,0]


Out[23]:
array([ 1, 18, 35])

Save sources


In [ ]:

ImageJ


In [21]:
from sima.misc.imagej import read_imagej_roi_zip

In [25]:
X = read_imagej_roi_zip('/Users/sofroniewn/github/regional/RoiSet.zip')

In [33]:
X[0]


Out[33]:
{u'label': u'0056-0064-0067', u'polygons': array([[ 39.,  79.,  55.],
        [ 39.,  80.,  55.],
        [ 38.,  80.,  55.],
        [ 38.,  81.,  55.],
        [ 38.,  82.,  55.],
        [ 38.,  83.,  55.],
        [ 38.,  84.,  55.],
        [ 38.,  85.,  55.],
        [ 38.,  86.,  55.],
        [ 38.,  87.,  55.],
        [ 38.,  88.,  55.],
        [ 38.,  89.,  55.],
        [ 39.,  89.,  55.],
        [ 40.,  89.,  55.],
        [ 41.,  90.,  55.],
        [ 42.,  90.,  55.],
        [ 44.,  90.,  55.],
        [ 45.,  90.,  55.],
        [ 46.,  90.,  55.],
        [ 48.,  90.,  55.],
        [ 49.,  89.,  55.],
        [ 50.,  88.,  55.],
        [ 50.,  87.,  55.],
        [ 50.,  86.,  55.],
        [ 49.,  85.,  55.],
        [ 48.,  84.,  55.],
        [ 48.,  83.,  55.],
        [ 47.,  82.,  55.],
        [ 46.,  81.,  55.],
        [ 46.,  80.,  55.],
        [ 45.,  80.,  55.],
        [ 44.,  80.,  55.],
        [ 43.,  80.,  55.],
        [ 42.,  80.,  55.],
        [ 41.,  80.,  55.],
        [ 40.,  80.,  55.]])}

In [38]:
Y[0]


Out[38]:
array([[78, 38],
       [79, 38],
       [79, 37],
       [80, 37],
       [81, 37],
       [82, 37],
       [83, 37],
       [84, 37],
       [85, 37],
       [86, 37],
       [87, 37],
       [88, 37],
       [88, 38],
       [88, 39],
       [89, 40],
       [89, 41],
       [89, 43],
       [89, 44],
       [89, 45],
       [89, 47],
       [88, 48],
       [87, 49],
       [86, 49],
       [85, 49],
       [84, 48],
       [83, 47],
       [82, 47],
       [81, 46],
       [80, 45],
       [79, 45],
       [79, 44],
       [79, 43],
       [79, 42],
       [79, 41],
       [79, 40],
       [79, 39]], dtype=int16)

In [35]:
Y = read_roi_zip('/Users/sofroniewn/github/regional/RoiSet.zip')

In [34]:
# Copyright: Luis Pedro Coelho <luis@luispedro.org>, 2012
# License: MIT
import numpy as np

def read_roi(fileobj):
    '''
    points = read_roi(fileobj)
    Read ImageJ's ROI format
    '''
# This is based on:
# http://rsbweb.nih.gov/ij/developer/source/ij/io/RoiDecoder.java.html
# http://rsbweb.nih.gov/ij/developer/source/ij/io/RoiEncoder.java.html


    SPLINE_FIT = 1
    DOUBLE_HEADED = 2
    OUTLINE = 4
    OVERLAY_LABELS = 8
    OVERLAY_NAMES = 16
    OVERLAY_BACKGROUNDS = 32
    OVERLAY_BOLD = 64
    SUB_PIXEL_RESOLUTION = 128
    DRAW_OFFSET = 256


    pos = [4]
    def get8():
        pos[0] += 1
        s = fileobj.read(1)
        if not s:
            raise IOError('readroi: Unexpected EOF')
        return ord(s)

    def get16():
        b0 = get8()
        b1 = get8()
        return (b0 << 8) | b1

    def get32():
        s0 = get16()
        s1 = get16()
        return (s0 << 16) | s1

    def getfloat():
        v = np.int32(get32())
        return v.view(np.float32)

    magic = fileobj.read(4)
    if magic != 'Iout':
        raise IOError('Magic number not found')
    version = get16()

    # It seems that the roi type field occupies 2 Bytes, but only one is used
    roi_type = get8()
    # Discard second Byte:
    get8()

    if not (0 <= roi_type < 11):
        raise ValueError('roireader: ROI type %s not supported' % roi_type)

    if roi_type != 7:
        raise ValueError('roireader: ROI type %s not supported (!= 7)' % roi_type)

    top = get16()
    left = get16()
    bottom = get16()
    right = get16()
    n_coordinates = get16()

    x1 = getfloat() 
    y1 = getfloat() 
    x2 = getfloat() 
    y2 = getfloat()
    stroke_width = get16()
    shape_roi_size = get32()
    stroke_color = get32()
    fill_color = get32()
    subtype = get16()
    if subtype != 0:
        raise ValueError('roireader: ROI subtype %s not supported (!= 0)' % subtype)
    options = get16()
    arrow_style = get8()
    arrow_head_size = get8()
    rect_arc_size = get16()
    position = get32()
    header2offset = get32()

    if options & SUB_PIXEL_RESOLUTION:
        getc = getfloat
        points = np.empty((n_coordinates, 2), dtype=np.float32)
    else:
        getc = get16
        points = np.empty((n_coordinates, 2), dtype=np.int16)
    points[:,1] = [getc() for i in xrange(n_coordinates)]
    points[:,0] = [getc() for i in xrange(n_coordinates)]
    points[:,1] += left
    points[:,0] += top
    points -= 1
    return points

def read_roi_zip(fname):
    import zipfile
    with zipfile.ZipFile(fname) as zf:
        return [read_roi(zf.open(n))
                    for n in zf.namelist()]

In [ ]: