importDicom.py 12 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359
  1. import dicom
  2. import vtkInterface as vi
  3. import os
  4. import vtk, qt, ctk, slicer
  5. from slicer.ScriptedLoadableModule import *
  6. import slicerNetwork
  7. import json
  8. import loadDicom
  9. import DICOMLib
  10. class importDicom(slicer.ScriptedLoadableModule.ScriptedLoadableModule):
  11. def __init__(self,parent):
  12. slicer.ScriptedLoadableModule.ScriptedLoadableModule.__init__(self, parent)
  13. self.className="importDICOM"
  14. self.parent.title="importDICOM"
  15. self.parent.categories = ["LabKey"]
  16. self.parent.dependencies = []
  17. self.parent.contributors = ["Andrej Studen (UL/FMF)"] # replace with "Firstname Lastname (Organization)"
  18. self.parent.helpText = """
  19. utilities for parsing dicom entries
  20. """
  21. self.parent.acknowledgementText = """
  22. Developed within the medical physics research programme of the Slovenian research agency.
  23. """ # replace with organization, grant and thanks.
  24. class importDicomWidget(slicer.ScriptedLoadableModule.ScriptedLoadableModuleWidget):
  25. """Uses ScriptedLoadableModuleWidget base class, available at:
  26. https://github.com/Slicer/Slicer/blob/master/Base/Python/slicer/ScriptedLoadableModule.py
  27. """
  28. def setup(self):
  29. slicer.ScriptedLoadableModule.ScriptedLoadableModuleWidget.setup(self)
  30. self.logic=importDicomLogic(self)
  31. self.network=slicerNetwork.labkeyURIHandler()
  32. connectionCollapsibleButton = ctk.ctkCollapsibleButton()
  33. connectionCollapsibleButton.text = "Connection"
  34. self.layout.addWidget(connectionCollapsibleButton)
  35. connectionFormLayout = qt.QFormLayout(connectionCollapsibleButton)
  36. self.loadConfigButton=qt.QPushButton("Load configuration")
  37. self.loadConfigButton.toolTip="Load configuration"
  38. self.loadConfigButton.connect('clicked(bool)',self.onLoadConfigButtonClicked)
  39. connectionFormLayout.addRow("Connection:",self.loadConfigButton)
  40. self.DICOMDirectory=qt.QLineEdit("Test/Temp/%40files/TEST/MLEM")
  41. connectionFormLayout.addRow("LabKey directory:",self.DICOMDirectory)
  42. #loadDICOMButton=qt.QPushButton("Load")
  43. #loadDICOMButton.toolTip="Load DICOM"
  44. #loadDICOMButton.clicked.connect(self.onLoadDICOMButtonClicked)
  45. #connectionFormLayout.addRow("DICOM:",loadDICOMButton)
  46. self.DICOMFilter=qt.QLineEdit('{"seriesNumber":"SeriesLabel"}')
  47. connectionFormLayout.addRow("Filter(JSON):",self.DICOMFilter)
  48. loadDICOMFilterButton=qt.QPushButton("Load with filter")
  49. loadDICOMFilterButton.toolTip="Load DICOM with filter"
  50. loadDICOMFilterButton.clicked.connect(self.onLoadDICOMFilterButtonClicked)
  51. connectionFormLayout.addRow("DICOM:",loadDICOMFilterButton)
  52. def onLoadConfigButtonClicked(self):
  53. filename=qt.QFileDialog.getOpenFileName(None,'Open configuration file (JSON)',
  54. os.path.join(os.path.expanduser('~'),'.labkey'), '*.json')
  55. self.network.parseConfig(filename)
  56. self.network.initRemote()
  57. self.loadConfigButton.setText(os.path.basename(filename))
  58. def onLoadDICOMFilterButtonClicked(self):
  59. filter=json.loads(self.DICOMFilter.text)
  60. #print("Filter is {}".format(filter))
  61. self.logic.loadVolumes(self.network,self.DICOMDirectory.text,filter)
  62. def onLoadDICOMButtonClicked(self):
  63. self.logic.load(self.network,self.DICOMDirectory.text)
  64. #equivalent of loadable + labkey interface
  65. class dicomSeries():
  66. def __init__(self):
  67. self.data = []
  68. self.idx = []
  69. self.z = []
  70. self.pixel_size = [0,0,0]
  71. self.lpsOrigin = [0,0,0]
  72. self.lpsOrigin[2]=1e30
  73. self.lpsOrientation=[0,0,0,0,0,0]
  74. self.local=False
  75. def getfile(self,net,relativePath):
  76. if self.local:
  77. return open(relativePath,'rb')
  78. return net.readFileToBuffer(relativePath)
  79. def addFile(self,f):
  80. try:
  81. self.files.append(f)
  82. except:
  83. self.files=[f]
  84. def setLabel(self,label):
  85. self.label=label
  86. def getLabel(self):
  87. try:
  88. return self.label
  89. except:
  90. return None
  91. def setMetadata(self,key,value):
  92. try:
  93. self.metadata[key]=value
  94. except:
  95. self.metadata={key:value}
  96. def getMetadata(self):
  97. try:
  98. return self.metadata
  99. except:
  100. return {}
  101. def load(self,net):
  102. for f in self.files:
  103. print '{}:'.format(f)
  104. fileBuffer=self.getfile(net,f)
  105. self.loadFile(fileBuffer)
  106. nz=len(self.idx)
  107. sh=self.data[-1].shape
  108. sh_list=list(sh)
  109. sh_list.append(nz)
  110. data_array=np.zeros(sh_list)
  111. for k in range(0,nz):
  112. kp=int(np.round((self.z[k]-self.center[2])/self.pixel_size[2]))
  113. data_array[:,:,kp]=np.transpose(self.data[k])
  114. try:
  115. nodeName='Series'+self.label
  116. except:
  117. print('Could not set series label')
  118. nodeName='UnknownSeries'
  119. newNode=slicer.vtkMRMLScalarVolumeNode()
  120. newNode.SetName(nodeName)
  121. ijkToRAS = vtk.vtkMatrix4x4()
  122. #think how to do this with image orientation
  123. rasOrientation=[-self.lpsOrientation[i] if (i%3 < 2) else self.lpsOrientation[i]
  124. for i in range(0,len(self.lpsOrientation))]
  125. rasOrigin=[-self.lpsOrigin[i] if (i%3<2) else self.lpsOrigin[i] for i in range(0,len(self.lpsOrigin))]
  126. for i in range(0,3):
  127. for j in range(0,3):
  128. ijkToRAS.SetElement(i,j,self.pixel_size[i]*rasOrientation[3*j+i])
  129. ijkToRAS.SetElement(i,3,rasOrigin[i])
  130. newNode.SetIJKToRASMatrix(ijkToRAS)
  131. v=vtk.vtkImageData()
  132. v.GetPointData().SetScalars(
  133. vtk.util.numpy_support.numpy_to_vtk(
  134. np.ravel(self.data,order='F'),deep=True, array_type=vtk.VTK_FLOAT))
  135. v.SetOrigin(0,0,0)
  136. v.SetSpacing(1,1,1)
  137. v.SetDimensions(self.data.shape)
  138. newNode.SetAndObserveImageData(v)
  139. slicer.mrmlScene.AddNode(newNode)
  140. volume={'node':newNode,'metadata':self.metadata}
  141. return volume
  142. def loadFile(self,fileBuffer):
  143. plan=dicom.read_file(fileBuffer)
  144. self.data.append(plan.pixel_array)
  145. self.idx.append(plan.InstanceNumber)
  146. self.z.append(plan.ImagePositionPatient[2])
  147. #pixelSize
  148. pixel_size=[plan.PixelSpacing[0],plan.PixelSpacing[1],
  149. plan.SliceThickness]
  150. for i in range(0,3):
  151. if self.pixel_size[i] == 0:
  152. self.pixel_size[i] = float(pixel_size[i])
  153. if abs(self.pixel_size[i]-pixel_size[i]) > 1e-3:
  154. print 'Pixel size mismatch {.2f}/{.2f}'.format(self.pixel_size[i],
  155. pixel_size[i])
  156. #origin
  157. for i in range(0,2):
  158. if self.lpsOrigin[i] == 0:
  159. self.lpsOrigin[i] = float(plan.ImagePositionPatient[i])
  160. if abs(self.lpsOrigin[i]-plan.ImagePositionPatient[i]) > 1e-3:
  161. print 'Image center mismatch {.2f}/{.2f}'.format(self.lpsOrigin[i],
  162. plan.ImagePositionPatient[i])
  163. #not average, but minimum (!) why??
  164. if plan.ImagePositionPatient[2]<self.lpsOrigin[2]:
  165. self.lpsOrigin[2]=plan.ImagePositionPatient[2]
  166. #orientation
  167. for i in range(0,6):
  168. if self.lpsOrientation[i] == 0:
  169. self.lpsOrientation[i] = float(plan.ImageOrientationPatient[i])
  170. if abs(self.lpsOrientation[i]-plan.ImageOrientationPatient[i]) > 1e-3:
  171. print 'Image orientation mismatch {0:.2f}/{1:.2f}'.format(self.lpsOrientation[i],
  172. plan.ImageOrientationPatient[i])
  173. return True
  174. class importDicomLogic(slicer.ScriptedLoadableModule.ScriptedLoadableModuleLogic):
  175. def __init__(self,parent):
  176. slicer.ScriptedLoadableModule.ScriptedLoadableModuleLogic.__init__(self, parent)
  177. self.local=False
  178. self.tag={
  179. 'studyInstanceUid':0x0020000d,
  180. 'seriesInstanceUid':0x0020000e,
  181. 'patientId':0x00100020,
  182. 'patientName':0x00100010,
  183. 'sequenceName':0x00180024,
  184. 'seriesNumber':0x00200011,
  185. 'percentPhaseFieldOfView':0x00180094,
  186. 'modality': 0x00080060,
  187. 'patientSex': 0x00100040,
  188. 'patientBirthDate': 0x00100030,
  189. 'patientComments': 0x00104000,
  190. 'studyDescription': 0x00081030,
  191. 'studyDate': 0x00080020,
  192. 'studyId': 0x00200010,
  193. 'studyTime': 0x00080030,
  194. 'frameOfReferenceInstanceUid':0x00200052}
  195. def setLocal(self, basePath):
  196. self.local=True
  197. self.basePath=basePath
  198. def loadVolumes(self,net,directory,filter):
  199. #mimic examineForImport
  200. seriesList=self.examineForImport(net,directory,filter)
  201. print("Got {} series").format(len(seriesList))
  202. volumes=[]
  203. for s in seriesList:
  204. try:
  205. volumes.append(s.load(net))
  206. #often fails, e.g. JPEGLossles
  207. except:
  208. loadable=DICOMLib.DICOMLoadable()
  209. loadable.name='Series'+str(s.getLabel())
  210. print("Loading for {} number of files (pre-load) {}").format(loadable.name,len(s.files))
  211. loadable.files=[net.DownloadFileToCache(f) for f in s.files]
  212. print("Loading for {} number of files (pre-sort) {}").format(loadable.name,len(loadable.files))
  213. loadable.files,distances,loadable.warning=DICOMLib.DICOMUtils.getSortedImageFiles(loadable.files,1e-3)
  214. print("Loading for {} number of files {}").format(loadable.name,len(loadable.files))
  215. try:
  216. volumeNode=self.volumePlugin.load(loadable)
  217. except:
  218. self.volumePlugin=slicer.modules.dicomPlugins['DICOMScalarVolumePlugin']()
  219. volumeNode=self.volumePlugin.load(loadable)
  220. volume={'node':volumeNode,'metadata':s.getMetadata()}
  221. volumes.append(volume)
  222. return volumes
  223. def listdir(self,net,relativeDirectory):
  224. if self.local:
  225. dirs=os.listdir(os.path.join(self.basePath,relativeDirectory))
  226. return [os.path.join(relativeDirectory,dir) for dir in dirs]
  227. return net.listRelativeDir(relativeDirectory)
  228. def getfile(self,net,relativePath):
  229. if self.local:
  230. return open(os.path.join(self.basePath,relativePath),'rb')
  231. return net.readFileToBuffer(relativePath)
  232. def examineForImport(self,net,directory,filter):
  233. #split by series
  234. seriesList=[]
  235. files=self.listdir(net,directory)
  236. if len(files)==0:
  237. print("No input found in {}".format(directory))
  238. return seriesList
  239. for f in files:
  240. fileBuffer=self.getfile(net,f)
  241. #validate
  242. try:
  243. plan = dicom.read_file(fileBuffer)
  244. except:
  245. print ("{}: Not a dicom file")
  246. continue
  247. #determine validity first
  248. fileValid=True
  249. for key in filter:
  250. if filter[key]==None:
  251. continue
  252. if filter[key]=='SeriesLabel':
  253. seriesTag=self.tag[key]
  254. continue
  255. v=plan[self.tag[key]].value
  256. if not v==filter[key]:
  257. print('Filter mismatch {}{:x}: {}/{}').format(key,self.tag[key],v,filter[key])
  258. fileValid=False
  259. if not fileValid:
  260. continue
  261. #determine serieslabel second
  262. seriesLabel=plan[seriesTag].value
  263. try:
  264. if series.getLabel()==seriesLabel:
  265. series.addFile(f)
  266. continue
  267. except NameError:
  268. pass
  269. #add new series
  270. seriesList.append(dicomSeries())
  271. series=seriesList[-1]
  272. series.local=self.local
  273. #set series parameters
  274. series.addFile(f)
  275. series.setLabel(seriesLabel)
  276. for key in filter:
  277. if not filter[key]==None:
  278. continue
  279. v=plan[self.tag[key]].value
  280. series.setMetadata(key,v)
  281. return seriesList