cardiacSPECT.py 24 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728
  1. import os
  2. import sys
  3. import unittest
  4. import vtk, qt, ctk, slicer
  5. from slicer.ScriptedLoadableModule import *
  6. import logging
  7. import parseDicom
  8. import vtkInterface as vi
  9. import fileIO
  10. import slicer
  11. import numpy as np
  12. import slicerNetwork
  13. #
  14. # cardiacSPECT
  15. #
  16. class cardiacSPECT(ScriptedLoadableModule):
  17. """Uses ScriptedLoadableModule base class, available at:
  18. https://github.com/Slicer/Slicer/blob/master/Base/Python/slicer/ScriptedLoadableModule.py
  19. """
  20. def __init__(self, parent):
  21. ScriptedLoadableModule.__init__(self, parent)
  22. parent.title = "Cardiac SPECT"
  23. parent.categories = ["Examples"]
  24. parent.dependencies = []
  25. parent.contributors = ["Andrej Studen (FMF/JSI)"] # replace with "Firstname Lastname (Org)"
  26. parent.helpText = """
  27. Load dynamic cardiac SPECT data to Slicer
  28. """
  29. parent.acknowledgementText = """
  30. This module was developed within the frame of the ARRS sponsored medical
  31. physics research programe to investigate quantitative measurements of cardiac
  32. function using sestamibi-like tracers
  33. """ # replace with organization, grant and thanks.
  34. self.parent.helpText += self.getDefaultModuleDocumentationLink()
  35. self.parent = parent
  36. #
  37. # cardiacSPECTWidget
  38. #
  39. class cardiacSPECTWidget(ScriptedLoadableModuleWidget):
  40. """Uses ScriptedLoadableModuleWidget base class, available at:
  41. https://github.com/Slicer/Slicer/blob/master/Base/Python/slicer/ScriptedLoadableModule.py
  42. """
  43. def setup(self):
  44. ScriptedLoadableModuleWidget.setup(self)
  45. self.selectRemote=fileIO.remoteFileSelector()
  46. try:
  47. self.network=slicer.modules.labkeySlicerPythonExtensionWidget.network
  48. except:
  49. self.network=slicerNetwork.labkeyURIHandler()
  50. self.logic=cardiacSPECTLogic(self)
  51. self.logic.setURIHandler(self.network)
  52. self.selectRemote.setMaster(self)
  53. # Instantiate and connect widgets ...
  54. dataButton = ctk.ctkCollapsibleButton()
  55. dataButton.text = "Data"
  56. self.layout.addWidget(dataButton)
  57. # Layout within the sample collapsible button
  58. dataFormLayout = qt.QFormLayout(dataButton)
  59. #pathGuess="file://"+os.environ['HOME']+"/SPECT"
  60. pathGuess="labkey://" + "dinamic_spect/%40files/Dinamika%20test2/SPECT_Dinamika_Rekonstruirano"
  61. self.dataPath=qt.QLineEdit(pathGuess)
  62. dataFormLayout.addRow("Data location",self.dataPath)
  63. self.remotePath=qt.QLineEdit();
  64. dataFormLayout.addRow('Remote Path', self.remotePath)
  65. self.remotePath.textChanged.connect(self.onRemotePathTextChanged)
  66. browseButton = qt.QPushButton("Browse local")
  67. browseButton.toolTip="Set file location"
  68. dataFormLayout.addRow("Select local",browseButton)
  69. browseButton.connect('clicked(bool)',self.onBrowseButtonClicked)
  70. browseRemoteButton = qt.QPushButton("Browse remote")
  71. browseRemoteButton.toolTip="Set remote location"
  72. dataFormLayout.addRow("Select remote",browseRemoteButton)
  73. browseRemoteButton.connect('clicked(bool)',self.onRemoteBrowseButtonClicked)
  74. dataLoadButton = qt.QPushButton("Load")
  75. dataLoadButton.toolTip="Load data from DICOM"
  76. dataFormLayout.addRow("Data",dataLoadButton)
  77. dataLoadButton.connect('clicked(bool)',self.onDataLoadButtonClicked)
  78. self.dataLoadButton = dataLoadButton
  79. self.patientId=qt.QLineEdit();
  80. dataFormLayout.addRow('Patient ID', self.patientId)
  81. patientLoadButton = qt.QPushButton("Load")
  82. patientLoadButton.toolTip="Load data from DICOM"
  83. dataFormLayout.addRow("Load Patient",patientLoadButton)
  84. patientLoadButton.clicked.connect(self.onPatientLoadButtonClicked)
  85. saveVolumeButton = qt.QPushButton("Save Volume to NRRD")
  86. saveVolumeButton.toolTip="Save to NRRD"
  87. dataFormLayout.addRow("Save to NRRD",saveVolumeButton)
  88. saveVolumeButton.clicked.connect(self.onSaveVolumeButtonClicked)
  89. transformNodeButton = qt.QPushButton("Transform Node")
  90. transformNodeButton.toolTip="Transform node with patient based transform"
  91. dataFormLayout.addRow("Transform Nodes",transformNodeButton)
  92. transformNodeButton.clicked.connect(self.onTransformNodeButtonClicked)
  93. # Add vertical spacer
  94. self.layout.addStretch(1)
  95. #addFrameButton=qt.QPushButton("Add Frame")
  96. #addFrameButton.toolTip="Add frame to VTK"
  97. #dataFormLayout.addWidget(addFrameButton)
  98. #addFrameButton.connect('clicked(bool)',self.onAddFrameButtonClicked)
  99. #addCTButton=qt.QPushButton("Add CT")
  100. #addCTButton.toolTip="Add CT to VTK"
  101. #dataFormLayout.addWidget(addCTButton)
  102. #addCTButton.connect('clicked(bool)',self.onAddCTButtonClicked)
  103. #
  104. # Parameters Area
  105. #
  106. parametersCollapsibleButton = ctk.ctkCollapsibleButton()
  107. parametersCollapsibleButton.text = "Parameters"
  108. self.layout.addWidget(parametersCollapsibleButton)
  109. # Layout within the dummy collapsible button
  110. parametersFormLayout = qt.QFormLayout(parametersCollapsibleButton)
  111. #
  112. # check box to trigger taking screen shots for later use in tutorials
  113. #
  114. hbox1=qt.QHBoxLayout()
  115. frameLabel = qt.QLabel()
  116. frameLabel.setText("Select frame")
  117. hbox1.addWidget(frameLabel)
  118. self.time_frame_select=qt.QSlider(qt.Qt.Horizontal)
  119. self.time_frame_select.valueChanged.connect(self.onTimeFrameSelect)
  120. #self.time_frame_select.connect('valueChanged()', self.onTimeFrameSelect)
  121. self.time_frame_select.setMinimum(0)
  122. self.time_frame_select.setMaximum(0)
  123. self.time_frame_select.setValue(0)
  124. self.time_frame_select.setTickPosition(qt.QSlider.TicksBelow)
  125. self.time_frame_select.setTickInterval(5)
  126. self.time_frame_select.toolTip = "Select the time frame"
  127. hbox1.addWidget(self.time_frame_select)
  128. parametersFormLayout.addRow(hbox1)
  129. hbox2 = qt.QHBoxLayout()
  130. meanROILabel = qt.QLabel()
  131. meanROILabel.setText("MeanROI")
  132. hbox2.addWidget(meanROILabel)
  133. self.meanROIVolume = qt.QLineEdit()
  134. self.meanROIVolume.setText("testVolume15")
  135. hbox2.addWidget(self.meanROIVolume)
  136. self.meanROISegment = qt.QLineEdit()
  137. self.meanROISegment.setText("Segment_1")
  138. hbox2.addWidget(self.meanROISegment)
  139. computeMeanROI = qt.QPushButton("Compute mean ROI")
  140. computeMeanROI.connect('clicked(bool)',self.onComputeMeanROIClicked)
  141. hbox2.addWidget(computeMeanROI)
  142. self.meanROIResult = qt.QLineEdit()
  143. self.meanROIResult.setText("0")
  144. hbox2.addWidget(self.meanROIResult)
  145. parametersFormLayout.addRow(hbox2)
  146. #row 3
  147. hbox3 = qt.QHBoxLayout()
  148. drawTimePlot=qt.QPushButton("Draw ROI time plot")
  149. drawTimePlot.connect('clicked(bool)',self.onDrawTimePlotClicked)
  150. hbox3.addWidget(drawTimePlot)
  151. parametersFormLayout.addRow(hbox3)
  152. #dataFormLayout.addWidget(hbox)
  153. #row 4
  154. hbox4 = qt.QHBoxLayout()
  155. countSegments=qt.QPushButton("Count segmentation segments")
  156. countSegments.connect('clicked(bool)',self.onCountSegmentsClicked)
  157. hbox4.addWidget(countSegments)
  158. self.countSegmentsDisplay=qt.QLineEdit()
  159. self.countSegmentsDisplay.setText("0")
  160. hbox4.addWidget(self.countSegmentsDisplay)
  161. parametersFormLayout.addRow(hbox4)
  162. #
  163. # Apply Button
  164. #
  165. self.applyButton = qt.QPushButton("Apply")
  166. self.applyButton.toolTip = "Run the algorithm."
  167. self.applyButton.enabled = False
  168. parametersFormLayout.addRow(self.applyButton)
  169. # connections
  170. self.applyButton.connect('clicked(bool)', self.onApplyButton)
  171. #self.inputSelector.connect("currentNodeChanged(vtkMRMLNode*)", self.onSelect)
  172. #self.outputSelector.connect("currentNodeChanged(vtkMRMLNode*)", self.onSelect)
  173. # Add vertical spacer
  174. self.layout.addStretch(1)
  175. self.resetPosition=1
  176. def cleanup(self):
  177. pass
  178. def onApplyButton(self):
  179. pass
  180. #logic = cardiacSPECTLogic()
  181. #imageThreshold = self.imageThresholdSliderWidget.value
  182. def onBrowseButtonClicked(self):
  183. startDir=self.dataPath.text
  184. inputDir=qt.QFileDialog.getExistingDirectory(None,
  185. 'Select DICOM directory',startDir)
  186. self.dataPath.setText("file://"+inputDir)
  187. def onRemoteBrowseButtonClicked(self):
  188. self.selectRemote.show()
  189. def onDataLoadButtonClicked(self):
  190. self.logic.loadData(self)
  191. def onRemotePathTextChanged(self,str):
  192. self.dataPath.setText('labkey://'+str)
  193. def onTimeFrameSelect(self):
  194. it=self.time_frame_select.value
  195. selectionNode = slicer.app.applicationLogic().GetSelectionNode()
  196. print("Propagating CT volume")
  197. nodeName=self.patientId.text+'CT'
  198. node=slicer.mrmlScene.GetFirstNodeByName(nodeName)
  199. selectionNode.SetReferenceActiveVolumeID(node.GetID())
  200. if self.resetPosition==1:
  201. self.resetPosition=0
  202. slicer.app.applicationLogic().PropagateVolumeSelection(1)
  203. else:
  204. slicer.app.applicationLogic().PropagateVolumeSelection(0)
  205. print("Propagating SPECT volume")
  206. nodeName=self.patientId.text+'Volume'+str(it)
  207. node=slicer.mrmlScene.GetFirstNodeByName(nodeName)
  208. selectionNode.SetSecondaryVolumeID(node.GetID())
  209. slicer.app.applicationLogic().PropagateForegroundVolumeSelection(0)
  210. node.GetDisplayNode().SetAndObserveColorNodeID('vtkMRMLColorTableNodeRed')
  211. lm = slicer.app.layoutManager()
  212. sID=['Red','Yellow','Green']
  213. for s in sID:
  214. sliceLogic = lm.sliceWidget(s).sliceLogic()
  215. compositeNode = sliceLogic.GetSliceCompositeNode()
  216. compositeNode.SetForegroundOpacity(0.5)
  217. #make sure the viewer is matched to the volume
  218. print("Done")
  219. #to access sliceLogic (slice control) use
  220. #lcol=slicer.app.layoutManager().mrmlSliceLogics() (vtkCollection)
  221. #vtkMRMLSliceLogic are named by colors (Red,Green,Blue)
  222. def onComputeMeanROIClicked(self):
  223. s=self.logic.meanROI(self.meanROIVolume.text,self.meanROISegment.text)
  224. self.meanROIResult.setText(str(s))
  225. def onDrawTimePlotClicked(self):
  226. n=self.time_frame_select.maximum+1
  227. ft=self.logic.frame_time
  228. #find number of segments
  229. ns = self.logic.countSegments()
  230. #add the chart node
  231. cn = slicer.mrmlScene.AddNode(slicer.vtkMRMLChartNode())
  232. for j in range(0,ns):
  233. #add node for data
  234. dn = slicer.mrmlScene.AddNode(slicer.vtkMRMLDoubleArrayNode())
  235. dn.SetName(self.patientId.text+'_'+self.logic.getSegmentName(j))
  236. a = dn.GetArray()
  237. a.SetNumberOfTuples(n)
  238. dt=0;
  239. t0=0;
  240. for i in range(0,n):
  241. vol=self.patientId.text+"Volume"+str(i)
  242. fx=ft[i]
  243. fy=self.logic.meanROI(vol,j)
  244. dt=2*ft[i]-t0
  245. t0+=dt
  246. a.SetComponent(i, 0, fx)
  247. a.SetComponent(i, 1, fy/dt)
  248. a.SetComponent(i, 2, 0)
  249. print("[{0} at {1:.2f}:{2:.2f}]".format(vol,fx,fy))
  250. #fish the number of the segment
  251. cn.AddArray(self.logic.getSegmentName(j), dn.GetID())
  252. cn.SetProperty('default', 'title', 'ROI time plot')
  253. cn.SetProperty('default', 'xAxisLabel', 'time [ms]')
  254. cn.SetProperty('default', 'yAxisLabel', 'Activity (arb)')
  255. #update the chart node
  256. cvns = slicer.mrmlScene.GetNodesByClass('vtkMRMLChartViewNode')
  257. if cvns.GetNumberOfItems() == 0:
  258. cvn = slicer.mrmlScene.AddNode(slicer.vtkMRMLChartViewNode())
  259. else:
  260. cvn = cvns.GetItemAsObject(0)
  261. cvn.SetChartNodeID(cn.GetID())
  262. def onCountSegmentsClicked(self):
  263. self.countSegmentsDisplay.setText(self.logic.countSegments())
  264. def onPatientLoadButtonClicked(self):
  265. self.logic.loadPatient(self,self.patientId.text)
  266. def onSaveVolumeButtonClicked(self):
  267. self.logic.storeVolumeNodes(self.patientId.text,
  268. self.time_frame_select.minimum,self.time_frame_select.maximum)
  269. def onTransformNodeButtonClicked(self):
  270. self.logic.applyTransform(self.patientId.text,
  271. self.time_frame_select.minimum,self.time_frame_select.maximum)
  272. #def onAddFrameButtonClicked(self):
  273. # it=int(self.time_frame_select.text)
  274. # self.logic.addFrame(it)
  275. # def onAddCTButtonClicked(self):
  276. # self.logic.addCT()
  277. #
  278. #
  279. # cardiacSPECTLogic
  280. #
  281. class cardiacSPECTLogic(ScriptedLoadableModuleLogic):
  282. """This class should implement all the actual
  283. computation done by your module. The interface
  284. should be such that other python code can import
  285. this class and make use of the functionality without
  286. requiring an instance of the Widget.
  287. Uses ScriptedLoadableModuleLogic base class, available at:
  288. https://github.com/Slicer/Slicer/blob/master/Base/Python/slicer/ScriptedLoadableModule.py
  289. """
  290. def __init__(self,parent):
  291. ScriptedLoadableModuleLogic.__init__(self, parent)
  292. self.pd=parseDicom.parseDicomLogic(self)
  293. def setURIHandler(self,net):
  294. self.net=net
  295. self.pd.setURIHandler(net)
  296. def loadData(self,widget):
  297. inputDir=str(widget.dataPath.text)
  298. self.pd.readMasterDirectory(inputDir)
  299. self.frame_data, self.frame_time, self.frame_origin, \
  300. self.frame_pixel_size, self.frame_orientation=self.pd.readNMDirectory(inputDir)
  301. self.ct_data,self.ct_origin,self.ct_pixel_size, \
  302. self.ct_orientation=self.pd.readCTDirectory(inputDir)
  303. self.ct_orientation=vi.completeOrientation(self.ct_orientation)
  304. self.frame_orientation=vi.completeOrientation(self.frame_orientation)
  305. self.addCT('test')
  306. self.addFrames('test')
  307. widget.time_frame_select.setMaximum(self.frame_data.shape[3]-1)
  308. #additional message via qt
  309. qt.QMessageBox.information(
  310. slicer.util.mainWindow(),
  311. 'Slicer Python','Data loaded')
  312. def loadPatient(self,widget,patientId):
  313. print("Loading {}").format(patientId)
  314. ds=self.net.loadDataset("dinamic_spect/Patients","Imaging")
  315. for r in ds['rows']:
  316. print r
  317. if r['aliasID']==patientId:
  318. visit=r
  319. print visit
  320. dicoms=self.net.loadDataset("Test/Transfer","Imaging")
  321. for r in dicoms['rows']:
  322. if not r['PatientId']==visit['aliasID']:
  323. continue
  324. if abs(r['SequenceNum']-float(visit['nmMaster']))<0.1:
  325. path=r['_labkeyurl_Series']
  326. path=path[:path.rfind('/')]
  327. masterPath="labkey://"+path
  328. if abs(r['SequenceNum']-float(visit['nmCorrData']))<0.1:
  329. path=r['_labkeyurl_Series']
  330. path=path[:path.rfind('/')]
  331. nmPath="labkey://"+path
  332. if abs(r['SequenceNum']-float(visit['ctData']))<0.1:
  333. path=r['_labkeyurl_Series']
  334. path=path[:path.rfind('/')]
  335. ctPath="labkey://"+path
  336. self.pd.readMasterDirectory(masterPath)
  337. self.frame_data, self.frame_time, self.frame_origin, \
  338. self.frame_pixel_size, self.frame_orientation=self.pd.readNMDirectory(nmPath)
  339. self.ct_data,self.ct_origin,self.ct_pixel_size, \
  340. self.ct_orientation=self.pd.readCTDirectory(ctPath)
  341. self.ct_orientation=vi.completeOrientation(self.ct_orientation)
  342. self.frame_orientation=vi.completeOrientation(self.frame_orientation)
  343. self.addCT(patientId)
  344. self.addFrames(patientId)
  345. widget.time_frame_select.setMaximum(self.frame_data.shape[3]-1)
  346. def addNode(self,nodeName,v, lpsOrigin, pixel_size, lpsOrientation, dataType):
  347. # if dataType=0 it is CT data, which gets propagated to background an is
  348. #used to fit the view field dimensions
  349. # if dataType=1, it is SPECT data, which gets propagated to foreground
  350. #and is not fit; keeping window set from CT
  351. #nodeName='testVolume'+str(it)
  352. newNode=slicer.vtkMRMLScalarVolumeNode()
  353. newNode.SetName(nodeName)
  354. #pixel_size=[0,0,0]
  355. #pixel_size=v.GetSpacing()
  356. #print(pixel_size)
  357. #origin=[0,0,0]
  358. #origin=v.GetOrigin()
  359. v.SetOrigin([0,0,0])
  360. v.SetSpacing([1,1,1])
  361. ijkToRAS = vtk.vtkMatrix4x4()
  362. #think how to do this with image orientation
  363. rasOrientation=[-lpsOrientation[i] if (i%3 < 2) else lpsOrientation[i]
  364. for i in range(0,len(lpsOrientation))]
  365. rasOrigin=[-lpsOrigin[i] if (i%3<2) else lpsOrigin[i] for i in range(0,len(lpsOrigin))]
  366. for i in range(0,3):
  367. for j in range(0,3):
  368. ijkToRAS.SetElement(i,j,pixel_size[i]*rasOrientation[3*j+i])
  369. ijkToRAS.SetElement(i,3,rasOrigin[i])
  370. newNode.SetIJKToRASMatrix(ijkToRAS)
  371. newNode.SetAndObserveImageData(v)
  372. slicer.mrmlScene.AddNode(newNode)
  373. def addFrames(self,patientId):
  374. #convert data from numpy.array to vtkImageData
  375. #use time point it
  376. print "NFrames: {}".format(self.frame_data.shape[3])
  377. for it in range(0,self.frame_data.shape[3]):
  378. frame_data=self.frame_data[:,:,:,it];
  379. nodeName=patientId+'Volume'+str(it)
  380. self.addNode(nodeName,
  381. vi.numpyToVTK(frame_data,frame_data.shape),
  382. self.frame_origin,
  383. self.frame_pixel_size,
  384. self.frame_orientation,1)
  385. def addCT(self,patientId):
  386. nodeName=patientId+'CT'
  387. self.addNode(nodeName,
  388. #vi.numpyToVTK3D(self.ct_data,
  389. # self.ct_origin,self.ct_pixel_size),
  390. vi.numpyToVTK(self.ct_data,self.ct_data.shape),
  391. self.ct_origin,self.ct_pixel_size,
  392. self.ct_orientation,0)
  393. def rFromI(i,volumeNode):
  394. ijkToRas = vtk.vtkMatrix4x4()
  395. volumeNode.GetIJKToRASMatrix(ijkToRas)
  396. vImage=volumeNode.GetImageData()
  397. i1=list(vImage.GetPoint(i))
  398. i1=i1.append(1)
  399. #ras are global coordinates (in mm)
  400. position_ras=ijkToRas.MultiplyPoint(i1)
  401. return position_ras[0:3]
  402. def IfromR(pos,volumeNode):
  403. fM=vtk.vtkMatrix4x4()
  404. volumeNode.GetRASToIJKMatrix(fM)
  405. fM.MultiplyPoint(pos)
  406. vImage=volumeNode.GetImageData()
  407. #nearest neighbor
  408. return vImage.FindPoint(pos[0:3])
  409. def meanROI(self, volName1, i):
  410. s=0
  411. #get the segmentation mask
  412. fNode=slicer.mrmlScene.GetNodesByClass("vtkMRMLSegmentationNode").GetItemAsObject(0)
  413. print "Found segmentation node: {}".format(fNode.GetName())
  414. segNode=slicer.vtkMRMLSegmentationNode.SafeDownCast(fNode)
  415. #no python bindings for vtkSegmentation
  416. #if segNode.GetSegmentation().GetNumberOfSegments()==0 :
  417. # print("No segments available")
  418. # return 0
  419. #edit here to change for more segments
  420. segment=segNode.GetSegmentation().GetNthSegmentID(int(i))
  421. mask = segNode.GetBinaryLabelmapRepresentation(segment)
  422. if mask==None:
  423. print("Segment {} not found".format(segment))
  424. return s
  425. print "Got mask for segment {}".format(segment)
  426. #get mask at (x,y,z)
  427. #mask.GetPointData().GetScalars().GetTuple1(mask.FindPoint([x,y,z]))
  428. #get the image data
  429. dataNode=slicer.mrmlScene.GetFirstNodeByName(volName1)
  430. dataImage=dataNode.GetImageData()
  431. # use IJK2RAS to get global coordinates
  432. ijkToRas = vtk.vtkMatrix4x4()
  433. dataNode.GetIJKToRASMatrix(ijkToRas)
  434. #iterate over volume pixelData
  435. n=dataImage.GetNumberOfPoints()
  436. extent=mask.GetExtent()
  437. fM=vtk.vtkMatrix4x4()
  438. mask.GetWorldToImageMatrix(fM)
  439. ns=0
  440. for i in range(0,n):
  441. #get global coordinates of point i
  442. [ix,iy,iz]=dataImage.GetPoint(i)
  443. position_ijk=[ix, iy, iz, 1]
  444. #ras are global coordinates (in mm)
  445. position_ras=ijkToRas.MultiplyPoint(position_ijk)
  446. fpos=[int(np.round(x)) for x in fM.MultiplyPoint(position_ras)]
  447. outOfRange=False
  448. for k in range(0,3):
  449. if fpos[k]<extent[2*k] or fpos[k]>extent[2*k+1]:
  450. outOfRange=True
  451. break;
  452. if outOfRange:
  453. continue
  454. #find point in mask with the same global coordinates
  455. maskValue=mask.GetPointData().GetScalars().GetTuple1(mask.ComputePointId(fpos[0:3]))
  456. if maskValue == 0:
  457. continue
  458. #use maskValue to project ROI data
  459. s+=maskValue*dataImage.GetPointData().GetScalars().GetTuple1(i)
  460. ns+=1
  461. return s/ns
  462. def countSegments(self):
  463. segNodeList=slicer.mrmlScene.GetNodesByClass("vtkMRMLSegmentationNode")
  464. if segNodeList.GetNumberOfItems()==0:
  465. return 0
  466. fNode=segNodeList.GetItemAsObject(0)
  467. segNode=slicer.vtkMRMLSegmentationNode.SafeDownCast(fNode)
  468. if fNode==None:
  469. return 0
  470. return segNode.GetSegmentation().GetNumberOfSegments()
  471. def getSegmentName(self,i):
  472. segNodeList=slicer.mrmlScene.GetNodesByClass("vtkMRMLSegmentationNode")
  473. if segNodeList.GetNumberOfItems()==0:
  474. return "NONE"
  475. fNode=segNodeList.GetItemAsObject(0)
  476. segNode=slicer.vtkMRMLSegmentationNode.SafeDownCast(fNode)
  477. if fNode==None:
  478. return "NONE"
  479. return segNode.GetSegmentation().GetSegment(segNode.GetSegmentation().GetNthSegmentID(i)).GetName()
  480. def storeNodeRemote(self,relativePath,nodeName):
  481. labkeyPath=self.pd.net.GetLabkeyWebdavUrl()+'/'+relativePath
  482. print ("Remote: {}").format(labkeyPath)
  483. #checks if exists
  484. self.pd.net.mkdir(labkeyPath)
  485. localPath=self.pd.net.GetLocalCacheDirectory()+'/'+relativePath
  486. node=slicer.mrmlScene.GetFirstNodeByName(nodeName)
  487. if node==None:
  488. print("Node {} not found").format(nodeName)
  489. return
  490. suffix=".nrrd"
  491. if node.__class__.__name__=="vtkMRMLDoubleArrayNode":
  492. suffix=".mcsv"
  493. if node.__class__.__name__=="vtkMRMLTransformNode":
  494. suffix=".h5"
  495. fileName=nodeName+suffix
  496. file=localPath+'/'+fileName
  497. slicer.util.saveNode(node,file)
  498. print("Stored to: {}").format(file)
  499. f=open(file)
  500. remoteFile=labkeyPath+'/'+fileName
  501. self.pd.net.put(remoteFile,f.read())
  502. def storeVolumeNodes(self,patientId,n1,n2):
  503. #n1=self.time_frame.minimum;
  504. #n2=self.time_frame.maximum
  505. project="dinamic_spect/Patients"
  506. relativePath=project+'/@files/'+patientId
  507. print("Store CT")
  508. nodeName=patientId+'CT'
  509. self.storeNodeRemote(relativePath,nodeName)
  510. print("Storing NM from {} to {}").format(n1,n2)
  511. n=n2-n1+1
  512. for i in range(n):
  513. it=i+n1
  514. nodeName=patientId+'Volume'+str(it)
  515. self.storeNodeRemote(relativePath,nodeName)
  516. segNodeName="Heart"
  517. self.storeNodeRemote(relativePath,segNodeName)
  518. doubleArrayNodeName=patientId+'_Ventricle'
  519. self.storeNodeRemote(relativePath,doubleArrayNodeName)
  520. transformNodeName=patientId+"_DF"
  521. self.storeNodeRemote(relativePath,transformNodeName)
  522. def applyTransform(self, patientId,n1,n2):
  523. transformNodeName=patientId+"_DF"
  524. transformNode=slicer.util.getFirstNodeByName(transformNodeName)
  525. if transformNode==None:
  526. print("Node {} not found").format(transformNodeName)
  527. return
  528. try:
  529. if self.transformApplied:
  530. print("Transform already applied")
  531. return
  532. except AttributeError:
  533. print("Not defined yet")
  534. self.transformApplied=True
  535. n=n2-n1+1
  536. for i in range(n):
  537. it=i+n1
  538. nodeName=patientId+'Volume'+str(it)
  539. node=slicer.util.getFirstNodeByName(nodeName)
  540. if node==None:
  541. continue
  542. node.ApplyTransform(transformNode.GetTransformToParent())
  543. nodeName=patientId+'CT'
  544. node=slicer.util.getFirstNodeByName(nodeName)
  545. if not node==None:
  546. node.ApplyTransform(transformNode.GetTransformToParent())
  547. class cardiacSPECTTest(ScriptedLoadableModuleTest):
  548. """
  549. This is the test case for your scripted module.
  550. Uses ScriptedLoadableModuleTest base class, available at:
  551. https://github.com/Slicer/Slicer/blob/master/Base/Python/slicer/ScriptedLoadableModule.py
  552. """
  553. def setUp(self):
  554. """ Do whatever is needed to reset the state - typically a scene clear will be enough.
  555. """
  556. slicer.mrmlScene.Clear(0)
  557. def runTest(self):
  558. """Run as few or as many tests as needed here.
  559. """
  560. self.setUp()
  561. self.test_cardiacSPECT1()
  562. def test_cardiacSPECT1(self):
  563. """ Ideally you should have several levels of tests. At the lowest level
  564. tests should exercise the functionality of the logic with different inputs
  565. (both valid and invalid). At higher levels your tests should emulate the
  566. way the user would interact with your code and confirm that it still works
  567. the way you intended.
  568. One of the most important features of the tests is that it should alert other
  569. developers when their changes will have an impact on the behavior of your
  570. module. For example, if a developer removes a feature that you depend on,
  571. your test should break so they know that the feature is needed.
  572. """
  573. self.delayDisplay("Starting the test")
  574. #
  575. # first, get some data
  576. #
  577. self.delayDisplay('Test passed!')