How to merge tiles obtained using openslide-python

I am trying to combine tiles in the correct order so they end up as the same whole slide image (.svs file).

The .svs file is read from a filepath according to the function beloew:

def open_slide(filepath = None):
    try:
        slide = openslide.open_slide(filepath)
       
    except OpenSlideError as o:
        print("Error" + str(o))
        
        slide = None
    except FileNotFoundError as f:
        print("Error" +  str(f))
        slide = None
    return slide

In the picture below I am trying the merge the tiles I got using openslide-python's DeepZoom generator (see code snippet below)

def create_tile_generator(slide, tile_size, overlap):
     gen = DeepZoomGenerator(slide, tile_size=tile_size, overlap=overlap, limit_bounds=False)

This is how I split the .svs into its tiles:

def split_wsi_to_tiles(wsi_path = None):
    print("splitting wsi into tiles")
    tile_indices = process_slide(slide_num = SLIDE_NUM , filepath= wsi_path, tile_size = TILE_SIZE, overlap = OVERLAP)
    i = 0
    tile_indices_savepath  = os.path.join(os.getcwd(),"saved","tile_indices")
    save_file(filepath = tile_indices_savepath,filename=name,file= tile_indices)
    for ti in tile_indices:
        suffix =  str(i)
        (slide_num,tile) = process_tile_index(tile_index =ti,filepath = svs_path )
        tile = cv2.cvtColor(tile, cv2.COLOR_BGR2RGB)
        cv2.imwrite(save_path + suffix + ext,tile)
        i = i + 1
    print("done splitting wsi into tiles")
    return tile_indices_savepath

The helper functions process_slide and process_tile_index are given below

def process_slide(slide_num =1 , filepath= None, tile_size = 256, overlap = 0):
    slide = open_slide(filepath = filepath)
    generator = create_tile_generator(slide, tile_size, overlap)
    zoom_level = get_40x_zoom_level(slide, generator)
    print("zoom level set to " + str(zoom_level))
    cols, rows = generator.level_tiles[zoom_level - 1]
    tile_indices = [(slide_num, tile_size, overlap, zoom_level, col, row)
                  for col in range(cols) for row in range(rows)]
    return tile_indices

def process_tile_index(tile_index=None,filepath= None):
        slide_num, tile_size, overlap, zoom_level, col, row = tile_index
        slide = open_slide(filepath = filepath)
        generator = create_tile_generator(slide, tile_size, overlap)
        tile = np.asarray(generator.get_tile(zoom_level, (col, row)))
        return (slide_num, tile)

The get_40x_zoom_level function description:

def get_40x_zoom_level(slide, generator):
    global level
    highest_zoom_level = generator.level_count - 1  # 0-based indexing
    try:
        mag = int(slide.properties[openslide.PROPERTY_NAME_OBJECTIVE_POWER])
        
        offset = math.floor((mag / 40) / 2)
        level = highest_zoom_level - offset
    except (ValueError, KeyError) as e:
        level = highest_zoom_level
    print("zoom level set at " +str(level) )
    save_file(filepath= os.path.join(os.getcwd(),"saved"),filename = "level.pickle",file = level)
    return level

This is how I try to merge the tiles back to its whole slide image (not necessarily in .svs format but the same image):

def merge_tiles_to_wsi(tile_path= None,wsi_path = None):
    print("merging tiles into wsi")
    tile_indices = load_file(filepath = tile_indices_savepath,filename = name)
    slide = open_slide(filepath = wsi_path)

    level = load_file(filepath= os.path.join(os.getcwd(),"saved"),filename = "level.pickle")
    generator = create_tile_generator(slide, TILE_SIZE, OVERLAP)
    slide_dims = generator.level_dimensions[level]
    row_size = slide_dims[0]
    col_size = slide_dims[1]
    channel_size = 3
    slide_shape = (row_size,col_size,channel_size)
    print("shape of slide is " + str(slide_shape))
    wsi = np.zeros(slide_shape)
    for ti in tile_indices:
        slide_num, tile_size, overlap, zoom_level, col, row  = ti
        generator = create_tile_generator(slide, tile_size, overlap)
        tile = np.asarray(generator.get_tile(zoom_level, (col, row)))
        
        row_length = tile.shape[0]
        col_length = tile.shape[1]
        row_end = row + row_length
        col_end = col + col_length
        print("col: " + str(col) + " row: " + str(row) + str(wsi[row:row_end,col:col_end].shape) + " " + str(tile.shape))
        wsi[row:row_end,col:col_end] = tile
        # view_image(img= wsi)
    print("merging tiles into wsi")

    return wsi

Here is what the final output looks like out.png

libvips can do this merge and join for you. You can call it from pyvips, the Python binding.

To load an svs image and split it into tiles you can write:

import pyvips

image = pyvips.Image.new_from_file("my-slide.svs")
image.dzsave("my-deepzoom")

And it'll write my-deepzoom.dzi and a directory, my-deepzoom_files, containing all the tiles. There are a lot of parameters you can adjust, see the chapter in the docs:

https://libvips.github.io/libvips/API/current/Making-image-pyramids.md.html

It's very fast and can make pyramids of any size on even modest hardware.

You can recombine tiles to form images with arrayjoin. You give it a list of images in row-major order and set across to the number of images per row. For example:

import pyvips

tiles = [pyvips.Image.new_from_file(f"{x}_{y}.jpeg", access="sequential")
         for y in range(height) for x in range(width)] 
image = pyvips.Image.arrayjoin(tiles, across=width)
image.write_to_file("huge.tif", compression="jpeg", tile=True)

It's very fast and can join extremely large arrays of images.