Commit 51cb03a6 authored by Georgi Tushev's avatar Georgi Tushev
Browse files

use pyvips buffer

parent c825bec2
.DS_Store
venv/
\ No newline at end of file
# PyStitcher
# CziStitcher
Reads a list of CZI files and stitches based on viewer coordinates
......@@ -8,3 +8,14 @@ Reads a list of CZI files and stitches based on viewer coordinates
[aicspylibczi](https://pypi.org/project/aicspylibczi/)
[numpy](https://pypi.org/project/numpy/)
## Installation
```
$ git clone https://gitlab.mpcdf.mpg.de/mpibr/scic/pystitcher.git
$ virtualenv -p $(which python3) venv
$ source venv/bin/activate
$ pip install --upgrade pip
$ pip install -r requirements.txt
$ python czistitcher.py /<path_to_czi_files>/image_name.queue
```
......@@ -34,14 +34,30 @@ def get_scene_grid(czi_queue: list) -> numpy.array:
return scene_grid
def get_scene_dims(czi_scene: numpy.array, dim: chr) -> int:
if dim == 'X':
return czi_scene[:, 0].max() + czi_scene[:, 2].max()
elif dim == 'Y':
return czi_scene[:, 1].max() + czi_scene[:, 3].max()
elif dim == 'C':
return czi_scene[:, -1].max()
else:
return -1
# def stitch_scene(czi_queue: list, czi_scene: numpy.array) -> numpy.ndarray:
# """stitch czi images in a single stack, after max-Z-projection"""
# stack_size_x = czi_scene[:, 0].max()
# stack_size_y = czi_scene[:, 1].max()
# stack_size_c = czi_scene[:, -1].max()
# stack_size_x = get_scene_dims(czi_scene, 'X')
# stack_size_y = get_scene_dims(czi_scene, 'Y')
# stack_size_c = get_scene_dims(czi_scene, 'C')
# stack = numpy.zeros((stack_size_x, stack_size_y, stack_size_c), dtype=numpy.uint16)
# for l in range(0, len(czi_queue)):
# print(czi_queue[l])
# czi_obj = aicspylibczi.CziFile(czi_queue[l])
# img, shp = czi_obj.read_image()
# img = numpy.squeeze(img)
......@@ -51,23 +67,19 @@ def get_scene_grid(czi_queue: list) -> numpy.array:
# scene_y = czi_scene[l, 1]
# scene_w = czi_scene[l, 2]
# scene_h = czi_scene[l, 3]
# print(czi_queue[l], scene_x, scene_y, scene_w, scene_h, img.shape)
# stack[scene_x:(scene_x + scene_w), scene_y:(scene_y + scene_h), :] = img
# break
# return stack
def get_scene_dims(czi_scene: numpy.array, dim: chr) -> int:
if dim == 'X':
return czi_scene[:, 0].max()
elif dim == 'Y':
return czi_scene[:, 1].max()
# #break
elif dim == 'C':
return czi_scene[:, -1].max()
# return stack
else:
return -1
# # print("convert to pyvips")
# # images = []
# # for c in range(0, stack_size_c):
# # images.append(pyvips.Image.new_from_array(stack[:, :, c]))
# # image = pyvips.Image.arrayjoin(images, across = 1)
# # image.set_type(pyvips.GValue.gint_type, "page-height", stack_size_y)
# # return image
def stitch_scene(czi_queue: list, czi_scene: numpy.array) -> pyvips.Image:
......@@ -85,15 +97,15 @@ def stitch_scene(czi_queue: list, czi_scene: numpy.array) -> pyvips.Image:
print(czi_queue[l])
czi_obj = aicspylibczi.CziFile(czi_queue[l])
img, shp = czi_obj.read_image()
img = numpy.squeeze(img)
img = numpy.max(img, axis=1)
img = numpy.transpose(img)
img = numpy.squeeze(img) #CZYX
img = numpy.max(img, axis=1) #CYX
scene_x = czi_scene[l, 0]
scene_y = czi_scene[l, 1]
for c in range(0, stack_size_c):
frame = pyvips.Image.new_from_array(img[:, :, c])
stack[0] = stack[0].insert(frame, scene_x, scene_y)
break
# pyvips will flip YX to XY with new_from_array()
frame = pyvips.Image.new_from_array(img[c, :, :])
stack[c] = stack[c].insert(frame, scene_x, scene_y)
#break
image = pyvips.Image.arrayjoin(stack, across = 1)
image.set_type(pyvips.GValue.gint_type, "page-height", stack_size_y)
......@@ -118,21 +130,26 @@ if __name__ == "__main__":
czi_queue = get_files_per_queue(sys.argv[1])
czi_scene = get_scene_grid(czi_queue)
image = stitch_scene(czi_queue, czi_scene)
metadata = get_meta_info(czi_queue[0])
metadata = replace_tag(metadata, "SizeX", get_scene_dims(czi_scene, 'X'))
metadata = replace_tag(metadata, "SizeY", get_scene_dims(czi_scene, 'Y'))
metadata = replace_tag(metadata, "SizeC", get_scene_dims(czi_scene, 'C'))
metadata = replace_tag(metadata, "SizeZ", 1)
metadata = re.sub("<\?xml version=\"1.0\"\?>\n","", metadata)
metadata = re.sub("\ +<SizeM>\d+</SizeM>\n","", metadata)
# metadata = get_meta_info(czi_queue[0])
# metadata = replace_tag(metadata, "SizeX", get_scene_dims(czi_scene, 'X'))
# metadata = replace_tag(metadata, "SizeY", get_scene_dims(czi_scene, 'Y'))
# metadata = replace_tag(metadata, "SizeC", get_scene_dims(czi_scene, 'C'))
# metadata = replace_tag(metadata, "SizeZ", 1)
# metadata = re.sub("<\?xml version=\"1.0\"\?>\n","", metadata)
# metadata = re.sub("\ +<SizeM>\d+</SizeM>\n","", metadata)
print("write tif ...")
outPath = os.path.dirname(czi_queue[0])
outName = os.path.basename(czi_queue[0])
outName = os.path.splitext(outName)[0] + "_stitched.tif"
outFile = os.path.join(outPath, outName)
image_bands = get_scene_dims(czi_scene, 'C')
image_width = get_scene_dims(czi_scene, 'X')
image_height = get_scene_dims(czi_scene, 'Y')
image.set_type(pyvips.GValue.gstr_type, "image-description",
f"""<?xml version="1.0" encoding="UTF-8"?>
<OME xmlns="http://www.openmicroscopy.org/Schemas/OME/2016-06"
......@@ -140,8 +157,16 @@ if __name__ == "__main__":
xsi:schemaLocation="http://www.openmicroscopy.org/Schemas/OME/2016-06 http://www.openmicroscopy.org/Schemas/OME/2016-06/ome.xsd">
<Image ID="Image:0">
<!-- Minimum required fields about image dimensions -->
<Pixels DimensionOrder="XYCZT"
<Pixels
BigEndian="false"
DimensionOrder="XYCZT"
ID="Pixels:0"
Interleaved="true"
PhysicalSizeX="0.10378322801247437"
PhysicalSizeXUnit="µm"
PhysicalSizeY="0.10378322801247437"
PhysicalSizeYUnit="µm"
SignificantBits="16"
SizeC="{image_bands}"
SizeT="1"
SizeX="{image_width}"
......@@ -152,104 +177,5 @@ if __name__ == "__main__":
</Image>
</OME>""")
image.set_type(pyvips.GValue.gstr_type, "image-description", metadata)
image.write_to_file("/Users/tushevg/Desktop/stitchme/stack.tif", pyramid=True, subifd=True, tile=True, compression="none")
# def ExportNDImage(image, fileName, xmlstr):
# image_width = image.shape[0]
# image_height = image.shape[1]
# image_bands = image.shape[2]
# images = []
# for c in range(0, image_bands):
# images.append(pyvips.Image.new_from_array(image[:,:,c]))
# out = pyvips.Image.arrayjoin(images, across=1)
# out.set_type(pyvips.GValue.gint_type, "page-height", images[0].height)
# out.set_type(pyvips.GValue.gstr_type, "image-description",
# f"""<?xml version="1.0" encoding="UTF-8"?>
# <OME xmlns="http://www.openmicroscopy.org/Schemas/OME/2016-06"
# xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance"
# xsi:schemaLocation="http://www.openmicroscopy.org/Schemas/OME/2016-06 http://www.openmicroscopy.org/Schemas/OME/2016-06/ome.xsd">
# <Image ID="Image:0">
# <!-- Minimum required fields about image dimensions -->
# <Pixels DimensionOrder="XYCZT"
# ID="Pixels:0"
# SizeC="{image_bands}"
# SizeT="1"
# SizeX="{image_width}"
# SizeY="{image_height}"
# SizeZ="1"
# Type="uint16">
# </Pixels>
# </Image>
# </OME>""")
# #out.set_type(pyvips.GValue.gstr_type, "image-description", xmlstr)
# out.write_to_file(fileName, pyramid=True, subifd=True, tile=True, compression="none")
# def ExportPyVipsImage(images, fileName, image_width, image_height, image_bands):
# out = pyvips.Image.arrayjoin(images, across=1)
# out.set_type(pyvips.GValue.gint_type, "page-height", images[0].height)
# out.set_type(pyvips.GValue.gstr_type, "image-description",
# f"""<?xml version="1.0" encoding="UTF-8"?>
# <OME xmlns="http://www.openmicroscopy.org/Schemas/OME/2016-06"
# xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance"
# xsi:schemaLocation="http://www.openmicroscopy.org/Schemas/OME/2016-06 http://www.openmicroscopy.org/Schemas/OME/2016-06/ome.xsd">
# <Image ID="Image:0">
# <!-- Minimum required fields about image dimensions -->
# <Pixels DimensionOrder="XYCZT"
# ID="Pixels:0"
# SizeC="{image_bands}"
# SizeT="1"
# SizeX="{image_width}"
# SizeY="{image_height}"
# SizeZ="1"
# Type="uint16">
# </Pixels>
# </Image>
# </OME>""")
# #out.set_type(pyvips.GValue.gstr_type, "image-description", xmlstr)
# out.write_to_file(fileName, pyramid=True, subifd=True, tile=True, compression="none")
# if __name__ == "__main__":
# fileQueue = "/Users/tushevg/Desktop/projects/InstinctSex/20220502_AFM_15431_I/20220502_AFM_15431_1d_I.queue"
# queue = GetFilesInQueue(fileQueue)
# grid = FindGridPoints(queue)
# ndims = FindNDims(queue)
# # calculate field of view
# image_width = (grid[:,0] + grid[:,2]).max() + 1
# image_height = (grid[:,1] + grid[:,3]).max() + 1
# image_bands = ndims['C']
# # allocate stitched image
# images = AllocateImages(image_width, image_height, image_bands)
# #ndimage = numpy.zeros((image_width, image_height, image_bands), dtype="uint16")
# #print(ndimage.shape)
# # fill
# k = 0
# for file in queue:
# czi = CziFile(file)
# data = czi.read_image()
# img = numpy.squeeze(data[0])
# img = numpy.max(img, axis = 1)
# #img = numpy.moveaxis(img, [0,1,2], [2,1,0])
# #xmlstr = ElementTree.tostring(czi.meta, encoding='unicode', method='xml')
# print(file, img.shape)
# xpos = grid[k,0]
# ypos = grid[k,1]
# #ndimage[xpos:(xpos+img.shape[0]), ypos:(ypos+img.shape[1]), :] = img
# for c in range(0, ndims['C']):
# slide = pyvips.Image.new_from_array(img[c,:,:])
# images[c] = images[c].insert(slide, xpos, ypos)
# #break
# k = k + 1
# # export
# print("export tif")
# #ExportNDImage(ndimage, "stitched.tif", xmlstr)
# ExportPyVipsImage(images, "stack.tif", image_width, image_height, image_bands)
#image.set_type(pyvips.GValue.gstr_type, "image-description", metadata)
image.write_to_file(outFile, pyramid=True, subifd=True, tile=True, compression="none")
\ No newline at end of file
aicspylibczi==3.0.5
cffi==1.15.1
numpy==1.23.1
pycparser==2.21
pyvips==2.2.1
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment