Generate CSV file
import csv
import random
import numpy as np
PI = np.pi
csvName = ‘/home/kong/400G/new_detect/augTable.csv’
# 71
with open(csvName, ‘wb’) as f:
csvwriter = csv.writer(f)
firstLine = [‘augID’]
firstLine.append(‘rotateAngleX’)
firstLine.append(‘rotateAngleY’)
firstLine.append(‘rotateAngleZ’)
firstLine.append(‘displacementX’)
firstLine.append(‘displacementY’)
firstLine.append(‘displacementZ’)
firstLine.append(‘zoom’)
csvwriter.writerow(firstLine)
for i in range(71):
if i < 20:
line = [i]
line.extend([0, 0, 0, 0, 0, 0, 1])
csvwriter.writerow(line)
else:
line = [i]
line.append(‘{:.2f}’.format(random.random() PI / 6))
line.append(‘{:.2f}’.format(random.random() PI / 6))
line.append(‘{:.2f}’.format(random.random() PI / 6))
line.append(random.randint(-5, 5))
line.append(random.randint(-5, 5))
line.append(random.randint(-5, 5))
line.append(‘{:.2f}’.format(1 + random.random() 0.2 - 0.1))
csvwriter.writerow(line)
Generate Binary files:
import csv
import os
import numpy as np
import SimpleITK as sitk
import time
import sys
import pprocess
import cv2
from scipy.ndimage.interpolation import zoom
CUT_SIZE = 64
SAVE_SIZE = 40
def createPatchCoodinates(centralCoor, size):
topMostPoint = np.array(
[centralCoor[0] - size / 2, centralCoor[1] - size / 2, centralCoor[2] - size / 2])
pointOne = np.array([centralCoor[0] + size / 2, centralCoor[1] -
size / 2, centralCoor[2] - size / 2])
pointTwo = np.array([centralCoor[0] - size / 2, centralCoor[1] +
size / 2, centralCoor[2] - size / 2])
pointThree = np.array(
[centralCoor[0] - size / 2, centralCoor[1] - size / 2, centralCoor[2] + size / 2])
coordinates = np.vstack((topMostPoint, pointOne, pointTwo, pointThree))
return coordinates
def translate(coordinates, translation):
translation = np.vstack((translation, translation, translation, translation))
coordinates = coordinates + translation
return coordinates
def zoom_diy(coordinates, zoomCenter, scale):
for i in range(coordinates.shape[0]):
coordinates[i, :] = (coordinates[i, :] - zoomCenter) * scale + zoomCenter
return coordinates
def rotate(coordinates, rotateCenter, rotateAngle):
rotation = sitk.Euler3DTransform()
rotation.SetCenter(rotateCenter)
rotation.SetRotation(rotateAngle[0], rotateAngle[1], rotateAngle[2])
coordinates = np.array(coordinates, np.int32)
for i in range(coordinates.shape[0]):
point = np.array(coordinates, np.float)
point = tuple(point[i, :])
rotatePoint = rotation.TransformPoint(point)
coordinates[i, :] = np.array(rotatePoint)
return coordinates
def pointBoardcastToMatrix(coordinates, patchSize):
topMostPoint, pointOne, pointTwo, pointThree = coordinates[
0], coordinates[1], coordinates[2], coordinates[3]
point1VectorZ = np.linspace(
0, pointOne[0] - topMostPoint[0], patchSize).reshape(patchSize, 1, 1)
point1VectorY = np.linspace(
0, pointOne[1] - topMostPoint[1], patchSize).reshape(patchSize, 1, 1)
point1VectorX = np.linspace(
0, pointOne[2] - topMostPoint[2], patchSize).reshape(patchSize, 1, 1)
point2VectorZ = np.linspace(
0, pointTwo[0] - topMostPoint[0], patchSize).reshape(patchSize, 1)
point2VectorY = np.linspace(
0, pointTwo[1] - topMostPoint[1], patchSize).reshape(patchSize, 1)
point2VectorX = np.linspace(
0, pointTwo[2] - topMostPoint[2], patchSize).reshape(patchSize, 1)
point3VectorZ = np.linspace(0, pointThree[0] - topMostPoint[0], patchSize)
point3VectorY = np.linspace(0, pointThree[1] - topMostPoint[1], patchSize)
point3VectorX = np.linspace(0, pointThree[2] - topMostPoint[2], patchSize)
xCoorMatrix = point1VectorX + point2VectorX + point3VectorX + topMostPoint[2]
yCoorMatrix = point1VectorY + point2VectorY + point3VectorY + topMostPoint[1]
zCoorMatrix = point1VectorZ + point2VectorZ + point3VectorZ + topMostPoint[0]
coordinates = [zCoorMatrix, yCoorMatrix, xCoorMatrix]
return coordinates
def neighbourInterpolate(coordinates, oriImage):
coordinates = np.round(coordinates).astype(np.int16)
zCoorMatrix, yCoorMatrix, xCoorMatrix = coordinates[0], coordinates[1], coordinates[2]
shape = zCoorMatrix.shape
zCoorLine = tuple(zCoorMatrix.reshape(1, -1)[0])
yCoorLine = tuple(yCoorMatrix.reshape(1, -1)[0])
xCoorLine = tuple(xCoorMatrix.reshape(1, -1)[0])
totalCoor = [zCoorLine, yCoorLine, xCoorLine]
interpolatedImage = oriImage[totalCoor].reshape(shape)
return interpolatedImage
def readCSV(filename):
‘’’read lines from a csv file.
‘’’
lines = []
with open(filename, “rb”) as f:
csvreader = csv.reader(f)
for line in csvreader:
lines.append(line)
return lines
def normalizePlanes(npzarray):
maxHU = 400.
minHU = -1000.
npzarray = (npzarray - minHU) / (maxHU - minHU)
npzarray[npzarray > 1] = 1.
npzarray[npzarray < 0] = 0.
return npzarray
def saveImageSlice(image, path):
if not os.path.exists(path):
os.makedirs(path)
for index in np.arange(image.shape[0]):
sliceTemp = image[index, :, :]
cv2.imwrite(path + str(index) + ‘.png’, sliceTemp * 255)
def worldToVoxelCoord(worldCoord, origin, outputSpacing):
stretchedVoxelCoord = np.absolute(worldCoord - origin)
voxelCoord = stretchedVoxelCoord / outputSpacing
return voxelCoord
def view_bar(args, num, total):
rate = 1.0 num / total
rate_num = int(rate 100)
r = ‘\rProcess {} task percentage (%): ‘.format(args) + ‘{}’.format(rate_num)
sys.stdout.write(r)
sys.stdout.flush()
def aug3d(borderedImage, rotateCenter, rotateAngle, displacement, zoomScale, patchSize):
coordinates = createPatchCoodinates(rotateCenter, patchSize)
coordinates = translate(coordinates, displacement)
coordinates = rotate(coordinates, rotateCenter, rotateAngle)
coordinates = zoom_diy(coordinates, rotateCenter, zoomScale)
coordinates = pointBoardcastToMatrix(coordinates, patchSize)
interpolatedImage = neighbourInterpolate(coordinates, borderedImage)
return interpolatedImage
def interpolatefilter(inputpath):
inputimage = sitk.ReadImage(inputpath)
origin = inputimage.GetOrigin()
spacing = inputimage.GetSpacing()
direction = inputimage.GetDirection()
outputspacing = (spacing[0], spacing[0], spacing[0])
size = inputimage.GetSize()
tmp = int(spacing[2] * size[2] / spacing[0])
if tmp % 2 != 0:
tmp = tmp - 1
outputsize = (size[0], size[1], tmp)
resamplefilter = sitk.ResampleImageFilter()
resamplefilter.SetOutputDirection(direction)
resamplefilter.SetSize(outputsize)
resamplefilter.SetOutputOrigin(origin)
resamplefilter.SetOutputSpacing(outputspacing)
outputimage = resamplefilter.Execute(inputimage)
outputsize = outputimage.GetSize()
numpyImage = sitk.GetArrayFromImage(outputimage)
numpyImage = normalizePlanes(numpyImage)
return numpyImage, list(outputsize), spacing
def createImageBorder(numpyImage, outputsize):
BackImage = np.zeros(((800, 800, 800)))
BackImage[400 - int(outputsize[2] / 2):400 + int(outputsize[2] / 2),
400 - int(outputsize[1] / 2):400 + int(outputsize[1] / 2),
400 - int(outputsize[0] / 2):400 + int(outputsize[0] / 2)] = numpyImage
return BackImage
def write_csv(processorID, line_start, line_end):
fBin = open(os.path.join(binOutputPath + ‘/‘ +
‘view{}’.format(processorID) + ‘.bin’), ‘wb’)
uidTemp = -1
for line in xrange(line_start, line_end):
cand = csvLines[line]
rotateAngle = np.array(augLines\[int(cand\[-1\])\]\[1:4\], np.float)\[::-1\]
displacement = np.array(augLines\[int(cand\[-1\])\]\[4:7\], np.int)\[::-1\]
zoomScale = float(augLines\[int(cand\[-1\])\]\[-1\])
if cand\[2\] != uidTemp:
mhdOriginalPath = originalPath + '/subset{}/{}.mhd'.format(cand\[1\], cand\[2\])
interpolatedImage, outputsize, spacing = interpolatefilter(mhdOriginalPath)
BackImage = createImageBorder(interpolatedImage, outputsize)
uidTemp = cand\[2\]
print '\\n{}'.format(cand\[2\])
rotateCenter = np.round(np.array(cand\[5:2:-1\], np.float) * spacing\[::-1\] / spacing\[0\])
rotateCenter += (400 - np.array(outputsize\[::-1\]) / 2)
augImage = aug3d(BackImage, rotateCenter, rotateAngle,
displacement, zoomScale, CUT_SIZE)
augImage = zoom(augImage, float(SAVE\_SIZE) / CUT\_SIZE)
# saveImageSlice(augImage, 'test/')
# imageCompare = BackImage\[rotateCenter\[0\] - 20:rotateCenter\[0\] + 20,
# rotateCenter\[1\] - 20:rotateCenter\[1\] + 20,
# rotateCenter\[2\] - 20:rotateCenter\[2\] + 20\]
# saveImageSlice(imageCompare, 'test1/')
fBin.write(chr(int(cand\[-2\])))
augImage.astype('float32').tofile(fBin)
# real time plot the task progress
if (line - line_start) % 10 == 0:
view\_bar(1, line - line\_start, line\_end - line\_start - 1)
elif line == line\_end - line\_start - 1:
view\_bar(1, line - line\_start, line\_end - line\_start - 1)
fBin.close()
# Global parameters
csvName = ‘/home/kong/400G/new_detect/LUNASTable.csv’
augName = ‘/home/kong/400G/new_detect/augTable.csv’
originalPath = ‘/home/kong/4T/nodule_project’
binOutputPath = ‘/home/kong/400G/new_detect’
csvLines = readCSV(csvName)[1:]
augLines = readCSV(augName)[1:]
ProcessorNum = pprocess.get_number_of_cores() - 1
def main():
series do
st = time.time()
write_csv(0, 0, len(csvLines))
print(‘\nseries costs {:.2f} seconds…’.format(time.time() - st))
parallel do
define parallel parameters
list_of_args = range(ProcessorNum)
group_num = len(csvLines) // ProcessorNum
cutPoint = np.empty([ProcessorNum, 2], dtype=int)
for row in range(ProcessorNum):
# start point
cutPoint\[row, 0\] = row * group_num
if row == ProcessorNum - 1:
# stop point
cutPoint\[row, 1\] = len(csvLines)
else:
# stop point
cutPoint\[row, 1\] = row * group\_num + group\_num - 1 + 1
starting parallel reading
st = time.time()
results = pprocess.Map()
parallel_function = results.manage(pprocess.MakeParallel(write_csv))
for args in list_of_args:
parallel_function(args, cutPoint[args, 0], cutPoint[args, 1])
print(‘\nStarting Parallel time {:.2f} seconds…’.format(time.time() - st))
st = time.time()
results[:]
parallel_results = results[:]
print(‘\nParallel costs {:.2f} seconds…’.format(time.time() - st))
if __name__ == ‘__main__‘:
main()
a=np.zeros((10,10),np.float)
a[0:8,1:6]=1
a[8,8]=1
print a
b=transform.rotate(a,60)
print b
fig=plt.figure()
ax = fig.add_subplot(2,1,1)
ax.imshow(a)
ax = fig.add_subplot(2,1,2)
ax.imshow(b)
plt.show()
print ‘finished’