|
| 1 | +import itk |
| 2 | +from scipy.signal import find_peaks |
| 3 | +from skimage.morphology import disk |
| 4 | +from skimage.morphology import dilation |
| 5 | +from skimage.feature import match_template |
| 6 | +import skvideo.io |
| 7 | +from itkpocus.util import get_framerate |
| 8 | +import skimage.filters |
| 9 | +import numpy as np |
| 10 | +import os |
| 11 | + |
| 12 | +def load_and_preprocess_video(fp, version=None): |
| 13 | + ''' |
| 14 | + Loads and preprocesses a Sonosite video. |
| 15 | + |
| 16 | + Parameters |
| 17 | + ---------- |
| 18 | + fp : str |
| 19 | + filepath to video (e.g. .mp4) |
| 20 | + version : optional |
| 21 | + Reserved for future use. |
| 22 | + |
| 23 | + Returns |
| 24 | + ------- |
| 25 | + img : itk.Image[itk.F,3] |
| 26 | + Floating point image with physical spacing set (3rd dimension is framerate), origin is at upper-left of image, and intensities are between 0.0 and 1.0 |
| 27 | + meta : dict |
| 28 | + Meta data (includes spacing and crop) |
| 29 | + ''' |
| 30 | + |
| 31 | + npvid_rgb = skvideo.io.vread(fp) |
| 32 | + return preprocess_video(npvid_rgb, framerate=get_framerate(skvideo.io.ffprobe(fp))) |
| 33 | + |
| 34 | +def load_and_preprocess_image(fp, version=None): |
| 35 | + ''' |
| 36 | + Loads and preprocesses a Sonosite image. |
| 37 | + |
| 38 | + Parameters |
| 39 | + ---------- |
| 40 | + fp : str |
| 41 | + filepath to image |
| 42 | + version : optional |
| 43 | + Reserved for future use. |
| 44 | + |
| 45 | + Returns |
| 46 | + ------- |
| 47 | + img : itk.Image[itk.F,2] |
| 48 | + Floating point image with physical spacing set (3rd dimension is framerate), origin is at upper-left of image, and intensities are between 0.0 and 1.0 |
| 49 | + meta : dict |
| 50 | + Meta data (includes spacing and crop) |
| 51 | + ''' |
| 52 | + npimg_rgb = itk.array_from_image(itk.imread(fp)) |
| 53 | + return preprocess_image(npimg_rgb) |
| 54 | + |
| 55 | +def _find_spacing(npimg): |
| 56 | + ''' |
| 57 | + Finds the white ticks on the right of the image and calculates their spacing. |
| 58 | + |
| 59 | + Parameters |
| 60 | + ----------- |
| 61 | + npimg (numpy array): 0-255 numpy array MxN |
| 62 | + |
| 63 | + Returns |
| 64 | + ------- |
| 65 | + |
| 66 | + spacing : float |
| 67 | + mm per pixel |
| 68 | + ''' |
| 69 | + rulerpos = [npimg.shape[1]-6, npimg.shape[1]] |
| 70 | + rulerimg = npimg[:,rulerpos[0]:rulerpos[1]].copy() |
| 71 | + |
| 72 | + ruler1d = np.max(rulerimg, axis=1) |
| 73 | + peaks, properties = find_peaks(ruler1d, height=[240,255], distance=5) |
| 74 | + rulerdiff = peaks[1:-1] - peaks[:-2] |
| 75 | + if np.std(rulerdiff) >= 2.0: |
| 76 | + raise ValueError("Error detecting ruler, peaks not equally-spaced {}".format(peaks)) |
| 77 | + pixelsize = 5.0 / np.mean(rulerdiff) # in mm, 0.5 cm ruler breaks |
| 78 | + |
| 79 | + return pixelsize |
| 80 | + |
| 81 | +def _normalize(npimg, npimgrgb): |
| 82 | + ''' |
| 83 | + A bunch of pixel-hacking to find the overlay elements in the image. Returns the overlay (necessary) |
| 84 | + for cropping later on and the image with the overlay elements in-filled (median filter at overlay pixels). |
| 85 | + |
| 86 | + Returns |
| 87 | + ------- |
| 88 | + npnorm : ndarray |
| 89 | + MxN normalized image to be cropped later (median filtered hud elements) |
| 90 | + npmasked : ndarray |
| 91 | + MxN hud image that is input to cropping algorithm |
| 92 | + ''' |
| 93 | + |
| 94 | + # try to detect hud elements by their color (ultrasound should be grey) and pure whiteness |
| 95 | + npcolor = np.logical_not(np.logical_and(npimgrgb[:,:,0] == npimgrgb[:,:,1], npimgrgb[:,:,1] == npimgrgb[:,:,2])) |
| 96 | + npwhite = npimg > 235 |
| 97 | + npblack = npimg < 10 # tries to get rid of black noise |
| 98 | + nphud = np.logical_or(npcolor, npwhite, npblack) |
| 99 | + nphud2 = dilation(nphud, disk(4)) |
| 100 | + npmasked = npimg.copy() |
| 101 | + npmasked[nphud2] = 0 |
| 102 | + |
| 103 | + # now we have a cropped image, need to get rid of the annotation marks |
| 104 | + nphud3 = dilation(nphud, disk(2)) |
| 105 | + |
| 106 | + npmed = skimage.filters.median(npimg, disk(4)) |
| 107 | + npnorm = npimg.copy() |
| 108 | + npnorm[nphud3] = npmed[nphud3] |
| 109 | + |
| 110 | + return npnorm, npmasked |
| 111 | + |
| 112 | + |
| 113 | +def _crop_image(npimg, crop): |
| 114 | + ''' |
| 115 | + Crops an ndarray. |
| 116 | + |
| 117 | + Parameters |
| 118 | + ---------- |
| 119 | + npimg : ndarray |
| 120 | + MxN |
| 121 | + crop : ndarray |
| 122 | + 2x2 [[ymin, ymax], [xmin, xmax]] where y=row and x=column |
| 123 | + |
| 124 | + Returns |
| 125 | + ------- |
| 126 | + ndarray |
| 127 | + ''' |
| 128 | + return npimg[crop[0,0]:crop[0,1]+1, crop[1,0]:crop[1,1]+1] |
| 129 | + |
| 130 | +def _find_crop(npnorm, npmasked): |
| 131 | + ''' |
| 132 | + A bunch of funkiness with removing the SonoSite overlay. Overlay elements can overlap the ultrasound image |
| 133 | + and jpeg compression can add noise. This works by a fixed crop in the y-direction and then a search on the |
| 134 | + overlay-removed image to find the x-direction boundaries. |
| 135 | + |
| 136 | + Returns |
| 137 | + ------- |
| 138 | + ndarray |
| 139 | + 2x2 [[ymin, ymax], [xmin, xmax]] where y=row, x=column |
| 140 | + ''' |
| 141 | + |
| 142 | + def array_crop(arr, threshold): |
| 143 | + c1 = 0 |
| 144 | + c2 = len(arr)-1 |
| 145 | + while (c1 <= c2 and arr[c1] < threshold): |
| 146 | + c1 += 1 |
| 147 | + while (c2 > c1 and arr[c2] < threshold): |
| 148 | + c2 -= 1 |
| 149 | + return np.array([c1, c2]) |
| 150 | + |
| 151 | + # oy, hard-code cropping, then looking for image intensity in the x-axis |
| 152 | + # zooming seems to stretch the image in x, but y is always maxed |
| 153 | + # confounding is grey-scale overlay elements and jpeg compression artifacts |
| 154 | + crop1 = np.array([[50, 620], [180, 830]]) |
| 155 | + npcrop1 = _crop_image(npmasked, crop1) |
| 156 | + xintensity = np.sum(npcrop1, axis=0) |
| 157 | + yintensity = np.sum(npcrop1, axis=1) |
| 158 | + xcrop = array_crop(xintensity, 2000) |
| 159 | + crop2 = np.array([crop1[0,:], [crop1[1,0] + xcrop[0], crop1[1,1]+xcrop[1]-npcrop1.shape[1]]]) |
| 160 | + |
| 161 | + return crop2 |
| 162 | + |
| 163 | +def preprocess_image(npimg, version=None): |
| 164 | + ''' |
| 165 | + Loads and preprocesses a Sonosite image. |
| 166 | + |
| 167 | + Parameters |
| 168 | + ---------- |
| 169 | + npimg : ndarray |
| 170 | + MxNx3 image array |
| 171 | + version : optional |
| 172 | + Reserved for future use. |
| 173 | +
|
| 174 | + Returns |
| 175 | + ------- |
| 176 | + img : itk.Image[itk.F,2] |
| 177 | + Floating point image with physical spacing set, origin is at upper-left of image, and intensities are between 0.0 and 1.0 |
| 178 | + meta : dict |
| 179 | + Meta data (includes spacing and crop) |
| 180 | + ''' |
| 181 | + |
| 182 | + npimg_bw = npimg[:,:,0].squeeze() |
| 183 | + |
| 184 | + spacing = _find_spacing(npimg_bw) |
| 185 | + npnorm, npmasked = _normalize(npimg_bw, npimg) |
| 186 | + crop = _find_crop(npnorm, npmasked) |
| 187 | + npnormcrop = _crop_image(npnorm, crop) # this is the us image for analysis |
| 188 | + |
| 189 | + img = itk.image_from_array((npnormcrop / 255.0).astype('float32')) |
| 190 | + img.SetSpacing(np.array([spacing, spacing])) |
| 191 | + return img, {'spacing' : np.array([spacing, spacing]), 'crop' : crop} |
| 192 | + |
| 193 | +def preprocess_video(npvid, framerate=1, version=None): |
| 194 | + ''' |
| 195 | + Loads and preprocesses a Sonosite video. |
| 196 | + |
| 197 | + Parameters |
| 198 | + ---------- |
| 199 | + ndarray : ndarray |
| 200 | + FxMxNxRGB |
| 201 | + version : optional |
| 202 | + Reserved for future use. |
| 203 | + |
| 204 | + Returns |
| 205 | + ------- |
| 206 | + img : itk.Image[itk.F,3] |
| 207 | + Floating point image with physical spacing set (3rd dimension is framerate), origin is at upper-left of image, and intensities are between 0.0 and 1.0 |
| 208 | + meta : dict |
| 209 | + Meta data (includes spacing and crop) |
| 210 | + ''' |
| 211 | + frm_cnt = npvid.shape[0] |
| 212 | + frm1rgb = npvid[0,:,:,:].squeeze() |
| 213 | + frm1 = frm1rgb[:,:,0].squeeze() |
| 214 | + spacing = _find_spacing(frm1) |
| 215 | + tmpnorm, tmpmasked = _normalize(frm1, frm1rgb) |
| 216 | + crop = _find_crop(tmpnorm, tmpmasked) |
| 217 | + tmpnormcrop = _crop_image(tmpnorm, crop) |
| 218 | + npnormvid = np.zeros([frm_cnt, tmpnormcrop.shape[0], tmpnormcrop.shape[1]]) |
| 219 | + |
| 220 | + for i in range(frm_cnt): |
| 221 | + frmrgb = npvid[i,:,:,:].squeeze() |
| 222 | + frm = frmrgb[:,:,0].squeeze() # just pick a channel for greyscale |
| 223 | + frmnorm, frmmasked = _normalize(frm, frmrgb) |
| 224 | + npnormvid[i,:,:] = _crop_image(frmnorm, crop) |
| 225 | + |
| 226 | + img = itk.image_from_array((npnormvid / 255.0).astype('float32')) |
| 227 | + tmp = np.array([spacing, spacing, framerate]) |
| 228 | + img.SetSpacing(tmp) |
| 229 | + return (npnormvid / 255.0).astype('float32'), {'spacing' : tmp, 'crop' : crop} |
0 commit comments