diff --git a/AmpScan/Measurements3.pdf b/AmpScan/Measurements3.pdf new file mode 100644 index 0000000000000000000000000000000000000000..664dd34f5d7b114e3f6130dfbefeaf7923df253a Binary files /dev/null and b/AmpScan/Measurements3.pdf differ diff --git a/AmpScan/Output_Template.pdf b/AmpScan/Output_Template.pdf new file mode 100644 index 0000000000000000000000000000000000000000..6cc1bc37a03df22bb883699f48d60d9b30213391 Binary files /dev/null and b/AmpScan/Output_Template.pdf differ diff --git a/AmpScan/ampVis.py b/AmpScan/ampVis.py index b33ab416f7e56f45c566a7d906a3ecbf5e3ca86b..9d8584294de3570d188baf76212843c009dcfa4b 100644 --- a/AmpScan/ampVis.py +++ b/AmpScan/ampVis.py @@ -36,6 +36,8 @@ class vtkRenWin(vtk.vtkRenderWindow): # self.rens[0].SetBackground(1.0,1.0,1.0) #self.rens[0].SetActiveCamera(self.cams[0]) # self.axes.append(vtk.vtkCubeAxesActor()) + self.markers = [] + self.labels = [] def renderActors(self, actors, viewport=0, zoom=1.0): r""" @@ -267,6 +269,88 @@ class vtkRenWin(vtk.vtkRenderWindow): writer.SetInputConnection(w2if.GetOutputPort()) writer.Write() + def Pick_point(self, loc): + """ + This receives coordinates in the GUI where user has picked, converts it + to mesh coordinates using vtkCellPicker, and places a sphere at picked + location as a visual aid. Also places a label at the point in the + render window. + TO-DO: Add functionality so user can type the label in, rather than + have it read 'Mid Patella' every time + """ + + x, y = loc + renderer = self.rens[0] + picker = vtk.vtkCellPicker() + picker.SetTolerance(0.01) + picker.Pick(x, y, 0, renderer) + points = picker.GetPickedPositions() + numPoints = points.GetNumberOfPoints() + #if no points selected, exits function + if numPoints<1: return + # Calls function to create a sphere at selected point + pnt = points.GetPoint(0) + self.mark(pnt[0], pnt[1], pnt[2]) + # Creating label at selected point + label = vtk.vtkStringArray() + label.SetName('label') + label.InsertNextValue(" Mid Patella") + lPoints = vtk.vtkPolyData() + lPoints.SetPoints(points) + lPoints.GetPointData().AddArray(label) + + hier = vtk.vtkPointSetToLabelHierarchy() + hier.SetInputData(lPoints) + hier.SetLabelArrayName('label') + hier.GetTextProperty().SetColor(0,0,0) + hier.GetTextProperty().SetFontSize(30) + + lMapper = vtk.vtkLabelPlacementMapper() + lMapper.SetInputConnection(hier.GetOutputPort()) + lMapper.SetBackgroundColor(0.3,0.3,0.3) + lMapper.SetBackgroundOpacity(0.8) + lMapper.SetMargin(10) + + lActor = vtk.vtkActor2D() + lActor.SetMapper(lMapper) + self.labels.append(lActor) # keep track of all label actors + + self.rens[0].AddActor(lActor) + self.Render() + return pnt + + + def mark(self, x,y,z): + """ + mark the picked point with a sphere + """ + sphere = vtk.vtkSphereSource() + sphere.SetRadius(3) + res = 20 + sphere.SetThetaResolution(res) + sphere.SetPhiResolution(res) + sphere.SetCenter(x,y,z) + mapper = vtk.vtkPolyDataMapper() + mapper.SetInputConnection(sphere.GetOutputPort()) + + self.marker = vtk.vtkActor() + self.marker.SetMapper(mapper) + self.rens[0].AddActor(self.marker) + self.marker.GetProperty().SetColor( (1,0,0) ) + self.markers.append(self.marker) #keep track of all marker actors + self.Render() + + def delMarker(self): + """ + removes the sphere marker and label from the renderer + """ + for i in self.markers: + self.rens[0].RemoveActor(i) + for i in self.labels: + self.rens[0].RemoveActor(i) + self.markers = [] + self.labels = [] + class qtVtkWindow(QVTKRenderWindowInteractor): r""" @@ -280,6 +364,7 @@ class qtVtkWindow(QVTKRenderWindowInteractor): self.iren = self._RenderWindow.GetInteractor() self.iren.Initialize() + class visMixin(object): r""" Set of visualisation methods that are contained within the AmpActor @@ -289,7 +374,7 @@ class visMixin(object): def genIm(self, size=[512, 512], views=[[0, -1, 0]], background=[1.0, 1.0, 1.0], projection=True, shading=True, mag=10, out='im', fh='test.tiff', - zoom=1.0, az = 0, crop=False, cam=None): + zoom=1.0, az = 0, el=0,crop=False, cam=None): r""" Creates a temporary off screen vtkRenWin which is then either returned as a numpy array or saved as a .png file @@ -335,6 +420,7 @@ class visMixin(object): win.setProjection(projection) win.SetSize(size[0], size[1]) win.Modified() + win.OffScreenRenderingOn() for i, view in enumerate(views): # win.addAxes([self.actor,], color=[0.0, 0.0, 0.0], viewport=i) @@ -386,7 +472,7 @@ class visMixin(object): win.setView() win.renderActors([self.actor,]) win.Render() - win.rens[0].GetActiveCamera().Azimuth(180) + win.rens[0].GetActiveCamera().Azimuth(0) win.rens[0].GetActiveCamera().SetParallelProjection(True) win.Render() return win @@ -605,6 +691,3 @@ class ampActor(vtk.vtkActor): self.GetProperty().LightingOn() if shading is False: self.GetProperty().LightingOff() - - - \ No newline at end of file diff --git a/AmpScan/analyse.py b/AmpScan/analyse.py index 81005f0c34b1ce00cbbffac0bb8971d26172a323..daa76dd58c69537db98031a7b61996f19e35ce42 100644 --- a/AmpScan/analyse.py +++ b/AmpScan/analyse.py @@ -7,19 +7,23 @@ Copyright: Joshua Steer 2018, Joshua.Steer@soton.ac.uk import numpy as np import matplotlib.pyplot as plt +import matplotlib.colors as clr +import matplotlib.colorbar as clb from mpl_toolkits.mplot3d import Axes3D from collections import defaultdict +from .output import getPDF +from math import floor #from .cython_ext import planeEdgeIntersect_cy, logEuPath_cy class analyseMixin(object): - r""" + """ Analysis methods to act upon a single AmpObject and generate a mpl figure """ - + def plot_slices(self, axis=2, slWidth=10): r""" Generate a mpl figure with information about the AmpObject @@ -279,3 +283,192 @@ class analyseMixin(object): (edges[:, axis+3] - edges[:, axis])) return intersectPoints + + + def MeasurementsOut(self, pos): + """ + Calculates perimeter of limb/cast at intervals from mid-patella to the + end of stump + Takes position of mid-patella (x,y,z) coordinates as input + Also creates images of limb views and graphs of CSA/Widths, which are + used in the PDF. + Calls the function responsible for adding the information to the PDF + template. + TODO: Split this into functions for each part i.e. Volume measure, CSA, + widths + """ + print(pos) + maxZ = [] + for i in [0,1,2]: + maxZ.append((self.vert[:, i]).max() - (self.vert[:, i]).min()) + #slice in longest axis of scan + self.axis = maxZ.index(max(maxZ)) + maxZ = max(maxZ) + zval = pos[self.axis] + # Get 6 equally spaced pts between mid-patella and stump end + slices = np.linspace(zval, (self.vert[:,self.axis]).min()+0.1, 6) + # uses create_slices + polys = self.create_slices(slices, axis=self.axis) + # calc perimeter of slices + perimeter = np.zeros([len(polys)]) + for i,poly in enumerate(polys): + nverts = np.arange(len(poly)-1) + dists = [] + for x in nverts: + xc = (poly[x][0] - poly[x+1][0])**2 + yc = (poly[x][1] - poly[x+1][1])**2 + zc = (poly[x][2] - poly[x+1][2])**2 + dist = np.sqrt(xc+yc+zc) + dists.append(dist) + perimeter[i] = sum(dists) / 10 + # distance between slice and mid-patella + lngth = (slices - zval) / 10 + #print(lngth, perimeter) + #generate png files of anterior and lateral views + self.genIm(out='fh',fh='lat.png',az=-90) + self.genIm(mag=1,out='fh',fh='ant.png') + #calculations at %length intervals of 10% + L = maxZ - ((self.vert[:,self.axis]).max()-zval)-10 + pL = np.linspace(0,1.2,13) + slices2 = [] + for i in pL: + slices2.append((self.vert[:,self.axis]).min()+10+(i*L)) + polys2 = self.create_slices(slices2,self.axis) + PolyArea = np.zeros([len(polys2)]) + MLWidth = np.zeros([len(polys2)]) + APWidth = np.zeros([len(polys2)]) + for i, poly in enumerate(polys2): + # Compute area of slice + area = 0.5*np.abs(np.dot(poly[:,0], np.roll(poly[:,1], 1)) - + np.dot(poly[:,1], np.roll(poly[:,0], 1))) + PolyArea[i] = area/100 + APW = poly[:,0].max() - poly[:,0].min() + APWidth[i] = APW/10 + MLW = poly[:,1].max() - poly[:,1].min() + MLWidth[i] = MLW/10 + print(PolyArea, MLWidth, APWidth) + fig = plt.figure() + fig.set_size_inches(7.5, 4.5) + ax = fig.add_subplot(221) + ax.plot(pL*100, PolyArea) + ax.set_xlabel("% length") + ax.set_ylabel("Area (cm^2)") + ax2 = fig.add_subplot(222) + ax2.plot(pL*100, MLWidth, 'ro',label='Medial-Lateral') + ax2.plot(pL*100, APWidth, 'b.',label='Anterior-Posterior') + ax2.set_xlabel("% length") + ax2.set_ylabel("width (cm)") + ax2.legend() + fig.savefig("figure.png") + getPDF(lngth, perimeter, PolyArea, APWidth,MLWidth) #PDF Creation function (in output.py) + #divided by 10 to convert to cms, assumes stl files are in mm + #TO-DO: Some sort of metric conversion function? + + + def CMapOut(self, colors): + """ + Colour Map with 4 views (copied Josh's code) + """ + titles = ['Anterior', 'Medial', 'Proximal', 'Lateral'] + fig,axes = plt.subplots(ncols=5) + cmap = clr.ListedColormap(colors, name='Amp') + norm = clr.Normalize(vmin=-10,vmax=10) + cb1 = clb.ColorbarBase(axes[-1], cmap=cmap,norm=norm) + cb1.set_label('Shape deviation / mm') + for i, ax in enumerate(axes[:-1]): + im = self.genIm(size=[3200, 8000],crop=True, az = i*90) + ax.imshow(im) + ax.set_title(titles[i]) + ax.set_axis_off() + #plt.colorbar(CMap) + fig.set_size_inches([12.5, 4]) + plt.savefig("Limb Views.png", dpi=600) + + + def volumeMeasure(self, zval, axis=2, slWidth=0.1): + """ + volume estimation using slices + """ + ind = [0,1,2] + ind.remove(axis) + slices = np.arange((self.vert[:,axis]).min()+slWidth, zval, slWidth) + polys = self.create_slices(slices, axis) + PArea = np.zeros(len(polys)) + for i, poly in enumerate(polys): + # Compute area of slice + area = 0.5*np.abs(np.dot(poly[:,ind[0]], np.roll(poly[:,ind[1]], 1)) - + np.dot(poly[:,ind[1]], np.roll(poly[:,ind[0]], 1))) + PArea[i] = area/100 + return sum(PArea*slWidth/10)+(PArea[-1]*(zval-slices[-1])) + + + def CSAMeasure(self,zval, axis=2, interval = 0.05): + """ + Measure CSA at 5% increments (0-95%) from selected point to stump + """ + ind = [0,1,2] + ind.remove(axis) + percents = np.arange(0.0,0.955,interval) + distump = (zval-(self.vert[:,axis]).min()) + slices = [] + [slices.append(zval - (distump*i)) for i in percents] + PArea = np.zeros(len(slices)) + for i,j in enumerate(slices): + try: + polys = self.create_slices([j], axis) + poly = polys[0] + area = 0.5*np.abs(np.dot(poly[:,ind[0]], np.roll(poly[:,ind[1]], 1)) - + np.dot(poly[:,ind[1]], np.roll(poly[:,ind[0]], 1))) + PArea[i] = area/100 + except IndexError: + PArea[i] = None + return PArea, distump + + def widthsMeasure(self,zval,axis=2,interval=0.05): + """ + Measure Coronal and Sagittal widths at intervals along limb (0-95%) + """ + percents = np.arange(0.0,0.955,interval) + distump = (zval-(self.vert[:,axis]).min()) + slices = [] + [slices.append(zval - (distump*i)) for i in percents] + MLWidth = np.zeros([len(slices)]) + APWidth = np.zeros([len(slices)]) + for i,j in enumerate(slices): + try: + polys = self.create_slices([j], axis) + poly = polys[0] + APW = poly[:,0].max() - poly[:,0].min() + APWidth[i] = APW + MLW = poly[:,1].max() - poly[:,1].min() + MLWidth[i] = MLW + except IndexError: + MLWidth[i] = None + APWidth[i] = None + return APWidth, MLWidth + + def perimeterMeasure(self,zval,axis=2,interval=0.05): + """ + Measure Coronal and Sagittal widths at intervals along limb (0-95%) + """ + percents = np.arange(0.0,0.955,interval) + distump = (zval-(self.vert[:,axis]).min()) + slices = [] + [slices.append(zval - (distump*i)) for i in percents] + perimeter = np.zeros([len(slices)]) + for i,j in enumerate(slices): + try: + polys = self.create_slices([j], axis) + poly = polys[0] + nverts = np.arange(len(poly)-1) + dists = [] + for x in nverts: + xc = (poly[x][0] - poly[x+1][0])**2 + yc = (poly[x][1] - poly[x+1][1])**2 + zc = (poly[x][2] - poly[x+1][2])**2 + dist = np.sqrt(xc+yc+zc) + dists.append(dist) + perimeter[i] = sum(dists) / 10 + except IndexError: + perimeter[i] = None + return perimeter, distump \ No newline at end of file diff --git a/AmpScan/output.py b/AmpScan/output.py new file mode 100644 index 0000000000000000000000000000000000000000..19cc0f8e3869eb639b0f138e7c70c799e80a1e82 --- /dev/null +++ b/AmpScan/output.py @@ -0,0 +1,49 @@ +from PyPDF2 import PdfFileReader, PdfFileWriter, PdfFileMerger +from reportlab.pdfgen import canvas +import io + +def getPDF(lngths, perimeters, CSA,APW,MLW): + """ + creates a PDF file containing information about the limb in correct + locations on the page + then merges the PDF file with the existing template to create the output + file + """ + packet = io.BytesIO() + c = canvas.Canvas(packet) + for i in range(1, len(lngths)-1): + stringl = "{}".format(abs(round(lngths[i],1))) + stringp = "{}".format(abs(round(perimeters[i],1))) + c.drawString(360+((i-1)*27), 474-((i-1)*41.5), stringl) + c.drawString(88, 524.5- ((i-1)*74.5), stringp) + stringmaxZ = "{}".format(abs(round(lngths[len(lngths)-1],1))) + c.drawString(514, 419, stringmaxZ) + c.setFont("Courier-Bold", 12) + c.drawString(65, 575, "Perimeter / cm") + c.drawString(400, 520, "Distance / cm") + c.showPage() + c.drawImage("ant.png", 38,225, 256,256) + c.drawImage("lat.png", 300,225,256,256) + c.drawImage("figure.png", -2.5,-50, 334,200) + for i in range(1,len(CSA),2): + sCSA = "{}".format(round(CSA[i],1)) + sAPW = "{}".format(round(APW[i],1)) + sMLW = "{}".format(round(MLW[i],1)) + c.drawString(403, 145-((i-1)*11.5), sCSA) + c.drawString(465, 145-((i-1)*11.5), sAPW) + c.drawString(520, 145-((i-1)*11.5), sMLW) + c.save() + packet.seek(0) + newpdf = PdfFileReader(packet) + template = PdfFileReader(open("Measurements3.pdf", "rb")) + t2 = PdfFileReader(open("Output Template.pdf", "rb")) + output = PdfFileWriter() + page = t2.getPage(0) + page.mergePage(newpdf.getPage(1)) + page2 = template.getPage(0) + page2.mergePage(newpdf.getPage(0)) + output.addPage(page) + output.addPage(page2) + outputStream = open("Output.pdf", "wb") + output.write(outputStream) + outputStream.close \ No newline at end of file diff --git a/GUIs/AmpScanGUI.py b/GUIs/AmpScanGUI.py index e74ef0cd4a492ddc51b1197cc57acc305ae3355c..555d4a50eda1cde64c61ef7a8db3191604cb295e 100644 --- a/GUIs/AmpScanGUI.py +++ b/GUIs/AmpScanGUI.py @@ -4,7 +4,7 @@ import vtk from AmpScan import AmpObject from AmpScan.registration import registration from AmpScan.align import align -from AmpScan.ampVis import qtVtkWindow +from AmpScan.ampVis import qtVtkWindow, vtkRenWin from PyQt5.QtCore import QPoint, QSize, Qt, QTimer, QRect, pyqtSignal from PyQt5.QtGui import (QColor, QFontMetrics, QImage, QPainter, QIcon, QOpenGLVersionProfile) @@ -12,7 +12,7 @@ from PyQt5.QtWidgets import (QAction, QApplication, QGridLayout, QHBoxLayout, QMainWindow, QMessageBox, QComboBox, QButtonGroup, QOpenGLWidget, QFileDialog, QLabel, QPushButton, QSlider, QWidget, QTableWidget, QTableWidgetItem, - QAbstractButton, QErrorMessage) + QAbstractButton, QCheckBox, QErrorMessage) class AmpScanGUI(QMainWindow): @@ -88,6 +88,12 @@ class AmpScanGUI(QMainWindow): return moving = str(self.alCont.moving.currentText()) self.files[moving].save(fname[0]) + try: + f = open(fname[0]+'.txt','w+') + f.write('{}'.format(self.pnt)) + except AttributeError: + print('A point has not been selected') + def display(self): render = [] @@ -120,6 +126,37 @@ class AmpScanGUI(QMainWindow): self.alCont.ztraButton.buttonClicked[QAbstractButton].connect(self.transz) else: show_message("Must be at least 1 object loaded to run align") + + def Point_Pick(self): + """ + Waits for a point click to occur before calling further functions + TODO: Create 'Picker controls'? Similar to Alignment controls, but where + user can enter the name of the point they select - this can allow + multiple landmark locations to be stored and marked? + """ + self.vtkWidget.iren.AddObserver('RightButtonPressEvent', self.pick_loc) + self.renWin.Render() + + def pick_loc(self, event, x): + """ + calcs the location of click in GUI (x,y) + calls function in ampVis.py which converts from GUI coordinates to + mesh coordinates and marks the point + """ + #print(event, x) + self.vtkWidget.iren.RemoveObservers('RightButtonPressEvent') + loc = event.GetEventPosition() + self.pnt = vtkRenWin.Pick_point(self.renWin, loc) + #vtkRenWin.mark(self.renWin,self.pnt[0],self.pnt[1],self.pnt[2]) + print(self.pnt) + + def removePick(self): + """ + delete all marked points and labels + TODO: be able to delete individual points? + """ + self.pnt = None + vtkRenWin.delMarker(self.renWin) def rotatex(self, button): moving = str(self.alCont.moving.currentText()) @@ -206,7 +243,6 @@ class AmpScanGUI(QMainWindow): def runRegistration(self): if self.objectsReady(2): - # Needs to be at least 2 files to run registration c1 = [31.0, 73.0, 125.0] c3 = [170.0, 75.0, 65.0] c2 = [212.0, 221.0, 225.0] @@ -221,7 +257,8 @@ class AmpScanGUI(QMainWindow): self.fileManager.setTable(target, [0,0,1], 0.5, 0) reg = registration(self.files[baseline], self.files[target], steps = 5, smooth=1).reg - reg.addActor(CMap = self.CMap02P) + #reg.addActor(CMap = self.CMap02P) + reg.addActor(CMap = self.CMapN2P) regName = target + '_reg' self.files[regName] = reg self.filesDrop.append(regName) @@ -230,7 +267,12 @@ class AmpScanGUI(QMainWindow): self.alCont.getNames() if hasattr(self, 'regCont'): self.regCont.getNames() - + #im = [] + if self.regCont.tick.isChecked() is True: + reg.actor.setScalarRange([-10,10]) + reg.actor.setShading(False) + reg.CMapOut(colors=self.CMapN2P) + reg.plotResults(name="distributionofshapevariance.png") print('Run the Registration code between %s and %s' % (baseline, target)) else: show_message("Must be at least 2 objects loaded to run registration") @@ -300,6 +342,14 @@ class AmpScanGUI(QMainWindow): self.renWin.renderActors(self.AmpObj.actors, ['socket', 'antS']) self.renWin.setScalarBar(self.AmpObj.actors['antS']) + def measure(self): + #if no point selected condition move to analyse.py + if self.pnt is None: + print("Please select a reference point first.") + else: + [name, _, color, opacity, display] = self.fileManager.getRow(0) + self.files[name].MeasurementsOut(self.pnt) + def createActions(self): """ Numpy style docstring. @@ -323,6 +373,12 @@ class AmpScanGUI(QMainWindow): triggered=self.register) self.analyse = QAction(QIcon('open.png'), 'Analyse', self, triggered=self.analyse) + self.pick = QAction(QIcon('open.png'), 'Pick', self, + triggered=self.Point_Pick) + self.removePick = QAction(QIcon('open.png'), 'Clear all picked points', self, + triggered = self.removePick) + self.Measure = QAction(QIcon('open.png'), 'Generate Measurements', self, + triggered = self.measure) self.openObjectManager = QAction(QIcon('open.png'), 'Show Object Manager', self, triggered=self.openAmpObjectManager) @@ -346,6 +402,11 @@ class AmpScanGUI(QMainWindow): self.analyseMenu.addAction(self.analyse) self.kineticMenu = self.menuBar().addMenu("&Kinetic Measurements") self.kineticMenu.addAction(self.openPress) + self.PointMenu = self.menuBar().addMenu("&Pick Point") + self.PointMenu.addAction(self.pick) + self.PointMenu.addAction(self.removePick) + self.measureMenu = self.menuBar().addMenu("Measure") + self.measureMenu.addAction(self.Measure) self.viewMenu = self.menuBar().addMenu("&View") self.viewMenu.addAction(self.openObjectManager) @@ -355,7 +416,7 @@ class AmpScanGUI(QMainWindow): def objectsReady(self, num): return len(self.files) >= num - + class fileManager(QMainWindow): """ Controls to manage the displayed @@ -365,7 +426,7 @@ class fileManager(QMainWindow): Perhaps an example implementation: >>> from AmpScan.AmpScanGUI import AmpScanGUI - +* """ def __init__(self, parent = None): @@ -476,7 +537,17 @@ class AlignControls(QMainWindow): self.moving.clear() self.moving.addItems(self.names) - + + + def getNames(self): + """ + """ + self.baseline.clear() + self.baseline.addItems(self.names) + self.target.clear() + self.target.addItems(self.names) + + class RegistrationControls(QMainWindow): """ Pop up for controls to align the @@ -496,13 +567,15 @@ class RegistrationControls(QMainWindow): self.baseline = QComboBox() self.target = QComboBox() self.reg = QPushButton("Run Registration") + self.tick = QCheckBox("Generate Output File for Comparison?") self.setCentralWidget(self.main) self.layout = QGridLayout() self.layout.addWidget(QLabel('Baseline'), 0, 0) self.layout.addWidget(QLabel('Target'), 1, 0) self.layout.addWidget(self.baseline, 0, 1) self.layout.addWidget(self.target, 1, 1) - self.layout.addWidget(self.reg, 2, 0, 1, -1) + self.layout.addWidget(self.tick, 2,1) + self.layout.addWidget(self.reg, 3, 0, 1, -1) self.main.setLayout(self.layout) self.setWindowTitle("Registration Manager") self.getNames()