Extracting patches using sliding window?

Hi 

Does the following code describe that data stored in numpy_array in. pkl format the author of the code try to extract patches using the sliding window function, I have following queries if someone can answer?
1- As this data is numpy_array and extracted patches are in numpy_array as well I want to ask is that possible to create images from it?
2- The data is of low dose CT images so I think the author extracted 3D patches from it?
3-To recreate it will be a 3D image (axial/Sagittal/Coronal)?
4-Size of extracted patches are 64*64?


Code:

import os

import joblib

import cv2

import numpy as np

def window_leveling(x, c=50, w=400, ymin=0, ymax=255):

    output = x.copy()

    wmin = (c - 0.5) - ((w - 1) / 2)

    wmax = (c - 0.5) + ((w - 1) / 2)

    output = ((wmin < x) & (x <= wmax)) * ((x - (c - 0.5)) / (w - 1) + 0.5) * (ymax - ymin) + ymin

    output[x <= wmin] = ymin

    output[x > wmax] = ymax

    return output.astype(np.uint8)

def main():

    path = 'cross_validation path'

    out_path = 'saved patches'

    emr_feature_number = [1, 6, 9, 11, 14, 15, 16]

    ct_slice_num = 80

    default_value = 0

    for i in range(5):

        ct_data = joblib.load(f'{path}/cross_validation_ct_threshold_{i+1}.pkl')

        voxel_list = joblib.load(f'{path}/cross_validation_voxel_list_postprocess_{i+1}.pkl')

        emr_data = joblib.load(f'{path}/emr_data_{i+1}.pkl')

       

        sub_path = f'{out_path}/fold_{i+1}/'

        os.makedirs(f'{sub_path}')

       

        for voxel_info in voxel_list:

            image_index, z, y, x, answer = voxel_info

            emr = np.array(emr_data[image_index])[emr_feature_number]

            ct = window_leveling(ct_data[image_index])

            start_index = (ct_slice_num - ct.shape[0]) //2

            np_extended = np.full((ct_slice_num, 512, 512), default_value, dtype=np.uint8)

            np_extended[start_index:start_index + ct.shape[0], :, :] = ct

            np_padding = np.pad(np_extended, ((10, 10), (65, 65), (65, 65)), 'constant', constant_values=default_value)

           

            ct_z = z + 10

            ct_y = y + 65

            ct_x = x + 65

           

            axial =  np_padding[ct_z, ct_y - 32:ct_y + 33, ct_x - 32:ct_x + 33]

            sagittal = np_padding[ct_z - 10:ct_z + 11, ct_y - 32:ct_y + 33, ct_x]

            coronal = np_padding[ct_z - 10:ct_z + 11, ct_y, ct_x - 32:ct_x + 33]

           

            data = [axial, sagittal, coronal, emr, int(answer)]

            save_path = f'{sub_path}/fold_{i+1}_{image_index:02d}_{answer}_{z:02d}_{y:03d}_{x:03d}.pkl'

            joblib.dump(data, save_path)

if __name__ == '__main__':

    main()

1. Is it possible to create images from it?

Yes, you can use the library matplotlib.pyplot to plot numpy arrays. If you have 2d images you can use something like imshow.

https://matplotlib.org/3.5.0/tutorials/introductory/images.html#sphx-glr-tutorials-introductory-images-py

If you have 3d data, there are some ways to visualize it, but if the data is dense, it will be difficult to understand what is going on.

https://matplotlib.org/2.0.2/mpl_toolkits/mplot3d/tutorial.html


2. Did the author extract 3d patches from it?

I believe the author is extrading three 2d patches. If you see how he extracts the slices there is always one dimension that does not “move”. For the axial dim ct_z remains constant, for sagittal is ct_x and for coronal is ct_y.

What this means is that for a given voxel you have three 2d images, one for each plane with the desired voxel in the middle.


3. Will it be a 3D image?

As we saw in the last question, it actually return three 2D images. If you wanted to get a 3D image you could do something like this.

scan_3d = np_padding[ct_z - 10:ct_z + 11, ct_y - 32:ct_y + 33, ct_x - 32:ct_x + 33]

However, as we saw in the first answer, it might be difficult to plot this 3D data in an understandable way.

Maybe one of these tutorials might helps you plotting the 3d data.


4. Size of extracted patches are 64*64?

Sizes:

  • Axial = Y: 64, X: 64
  • Sagittal = Z: 20, Y: 64
  • Coronal = Z: 20, X: 64

Hope this helps :smile:

thank you, you explain everything very well :slight_smile:

1 Like

Sorry to bother you again but as you can see these are 2D images of axial/sagittal and coronal planes.
I want to implement 3D CNN model for classification and further I want to implement localization techniques as well for automatic calcium scoring.

This is an example of 2D CNN model implement on extract patches

but I am little confuse and following is my question if you can guide me

How can I change this model into 3D CNN. As this is implemented on all X/y/z planes as I follow some articles but there is no example of all planes they just implement it on Axial images?

Thank you in advance

If you want to analyze your data in 3D I would really recommend this book written by @tom et. al.

In the first part they teach the most important aspects of pytorch.

The second part is completely dedicated to lung cancer detection in 3D.

Since this is a big topic, I really thing this book will help you more to understand what is happening.

You can read it for free in the link that I gave you.

1 Like