preprocess.py 6.6 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239
  1. import os
  2. import json
  3. import re
  4. import subprocess
  5. import nibabel
  6. import shutil
  7. import sys
  8. shome=os.path.expanduser('~nixUser')
  9. fhome=os.path.expanduser('~')
  10. with open(os.path.join(fhome,".labkey","setup.json"),"w") as f:
  11. setup=json.load(f)
  12. sys.path.insert(1,setup["paths"]["labkeyInterface"])
  13. import labkeyInterface
  14. import labkeyDatabaseBrowser
  15. import labkeyFileBrowser
  16. sys.path.insert(1,setup["paths"]["orthancInterface"])
  17. import orthancInterface
  18. import orthancFileBrowser
  19. #sys.path.insert(1,shome+'/software/src/IPNUMM/dicomUtils')
  20. #import loadDicom
  21. fconfig=os.path.join(fhome,'.labkey','network.json')
  22. #matlab=os.path.join("/","data","software","install","matlab","bin","matlab")
  23. matlab=setup["paths"]["matlab"]
  24. #os.path.join(fhome,"software","install","matlab","bin","matlab")
  25. generalCodes=setup["paths"]["generalCodes"]
  26. #ios.path.join(fhome,"software","src","generalCodes")
  27. niftiTools=setup["paths"]["niftiTools"]
  28. #niftiTools=os.path.join(fhome,"software","src","NifTiScripts")
  29. net=labkeyInterface.labkeyInterface()
  30. net.init(fconfig)
  31. db=labkeyDatabaseBrowser.labkeyDB(net)
  32. fb=labkeyFileBrowser.labkeyFileBrowser(net)
  33. onet=orthancInterface.orthancInterface()
  34. onet.init(fconfig)
  35. ofb=orthancFileBrowser.orthancFileBrowser(onet)
  36. hi=0
  37. project='iPNUMMretro/Study'
  38. dataset='Imaging'
  39. h
  40. tempBase=os.path.join(fhome,'temp')
  41. #all images from database
  42. ds=db.selectRows(project,'study',dataset,[])
  43. imageSelector={"CT":"CT","PETWB":"PET"};
  44. imageResampledField={"CT":"ctResampled","PETWB":"petResampled"}
  45. #projectNIfTIBase=os.path.join(labkeyBase,'files',project,'@files/nifti')
  46. #use webdav to transfer file (even though it is localhost)
  47. def getPatientLabel(row):
  48. return row['PatientId'].replace('/','_')
  49. def getVisitLabel(row):
  50. return 'VISIT_'+str(int(row['SequenceNum']))
  51. def getStudyLabel(row):
  52. return getPatientLabel(row)+'-'+getVisitLabel(row)
  53. def runPreprocess_DM(matlab,generalCodes,niftiTools,studyDir):
  54. #run after all directories have been assembled
  55. script="addpath('"+generalCodes+"');"
  56. script+="addpath('"+niftiTools+"');"
  57. script+="preprocess_DM('"+studyDir+"',0,0);"
  58. script+="test;"
  59. script+="quit();"
  60. #outText=subprocess.check_output(["/bin/echo",script])
  61. try:
  62. outText=subprocess.check_output([matlab,"-nojvm","-r",script])
  63. except subprocess.CalledProcessError as e:
  64. print("Failed with:\n{}".format(e.output.decode('utf-8')))
  65. return False
  66. print(outText.decode('utf-8'))
  67. return True
  68. def getDicom(ofb,row,zipDir,rawDir,im,imageSelector):
  69. #Load the dicom zip file and unzips it. If zip file is already at the expected path, it skips the loading step
  70. #Return True for valid outcome and False for troubles in row formating or unzip failures
  71. seriesId=row[im];
  72. if seriesId=="0":
  73. return False
  74. print("{}: {}".format(im,seriesId))
  75. fname=os.path.join(zipDir,getStudyLabel(row)+'_'+imageSelector[im]+".zip");
  76. #copy data from orthanc
  77. if os.path.isfile(fname):
  78. print("Data already loaded. Skipping")
  79. else:
  80. print("Loading data from orthanc")
  81. ofb.getZip('series',seriesId,fname)
  82. #unzip the zipped dicom series
  83. unzipDir=os.path.join(rawDir,imageSelector[im])
  84. if os.path.isdir(unzipDir):
  85. print("Data already unzipped")
  86. return True
  87. try:
  88. os.mkdir(unzipDir)
  89. except FileExistsError:
  90. shutil.rmtree(unzipDir)
  91. try:
  92. outTxt=subprocess.check_output(["unzip","-d",unzipDir,"-xj",fname])
  93. except subprocess.CalledProcessError:
  94. print("unzip failed for {}".format(fname))
  95. return False
  96. return True
  97. def updateRow(project,dataset,row,imageResampledField,gzFileNames):
  98. row['patientCode']=getPatientLabel(row)
  99. row['visitCode']=getVisitLabel(row)
  100. for im in imageResampledField:
  101. row[imageResampledField[im]]=gzFileNames[im]
  102. db.modifyRows('update',project,'study',dataset,[row])
  103. i=0
  104. for row in ds["rows"]:
  105. #interesting files are processedDir/studyName_CT_notCropped_2mmVoxel.nii
  106. #asn processedDir/studyName_PET_notCropped_2mmVoxel.nii
  107. volumeFileNames={im:\
  108. getStudyLabel(row)+'_'+imageSelector[im]+
  109. '_notCropped_2mmVoxel.nii'\
  110. for im in imageSelector}
  111. gzFileNames={im:f+".gz" \
  112. for (im,f) in volumeFileNames.items()}
  113. #build/check remote directory structure
  114. remoteDir=fb.buildPathURL(project,['preprocessedImages',getPatientLabel(row),getVisitLabel(row)])
  115. gzRemoteFiles={im:remoteDir+'/'+f\
  116. for (im,f) in gzFileNames.items()}
  117. remoteFilePresent=[fb.entryExists(f)\
  118. for f in gzRemoteFiles.values()]
  119. for f in gzRemoteFiles.values():
  120. print("[{}]: [{}]".format(f,fb.entryExists(f)))
  121. if all(remoteFilePresent):
  122. print("Entry for row done.")
  123. updateRow(project,dataset,row,imageResampledField,\
  124. gzFileNames)
  125. continue
  126. #setup the directory structure for preprocess_DM
  127. studyDir=os.path.join(tempBase,getStudyLabel(row))
  128. if not os.path.isdir(studyDir):
  129. os.mkdir(studyDir)
  130. rawDir=os.path.join(studyDir,'Raw')
  131. if not os.path.isdir(rawDir):
  132. os.mkdir(rawDir)
  133. zipDir=os.path.join(studyDir,'Zip')
  134. if not os.path.isdir(zipDir):
  135. os.mkdir(zipDir)
  136. processedDir=os.path.join(studyDir,'Processed')
  137. if not os.path.isdir(processedDir):
  138. os.mkdir(processedDir)
  139. #specify local file names with path
  140. volumeFiles={im:os.path.join(processedDir,f)\
  141. for (im.f) in volumeFileNames.items()}
  142. gzFiles={im:f+".gz"\
  143. for (im,f) in volumeFiles.items()}
  144. filesPresent=[os.path.isfile(f) for f in gzFiles.values()]
  145. if not all(filesPresent):
  146. for im in imageSelector:
  147. #checks if raw files are already loaded
  148. getDicom(ofb,row,zipDir,rawDir,im,imageSelector)
  149. #preprocess and zip
  150. ok=runPreprocess_DM(matlab,generalCodes,niftiTools,studyDir)
  151. if not ok:
  152. shutil.rmtree(studyDir)
  153. continue
  154. for f in volumeFiles.values():
  155. print("Running gzip {}".format(f))
  156. outText=subprocess.check_output(["/bin/gzip",f])
  157. print(outText.decode('utf-8'))
  158. #upload local files to remote
  159. for im in gzFiles:
  160. #for local,remote in zip(gzFiles,gzRemoteFiles):
  161. local=gzFiles[im]
  162. remote=gzRemoteFiles[im]
  163. print("Uploading {}".format(local))
  164. fb.writeFileToFile(local,remote)
  165. #update row and let it know where the processed files are
  166. updateRow(project,dataset,row,imageResampledField,gzFileNames)
  167. #cleanup
  168. shutil.rmtree(studyDir)
  169. if i==-1:
  170. break
  171. i=i+1
  172. print("Done")