Search code examples
pythonautomationsubdirectory

Run Python script in each subfolder automatically

I have this file structure:

  1. Main Folder: W001-W100.

    • Subfolder: W001

      • Around 500 DICOM files (I am working with CT data).

    • Subfolder: W002

      • Around 500 DICOM files (I am working with CT data).

    • Subfolder: ...

    • Subfolder: W100

      • Around 500 DICOM files (I am working with CT data).

enter image description here

I am using Slicer software and Python to automate the process of getting surface meshes automatically.

This is the Python script I wrote that needs to be run in each subfolder.

As you can see, there are two lines with the location of the subfolder:

  • Second line, specifying the location of the subfolder
  • Almost the last line, indicating where to save the .ply file.

As I have around 1000 subfolders, I don't want to run the script changing those lines. In other words, I wan to automate the process.

I am much more familiar with R than with Python. So:

  • could please indicate me how can I automate the process in detail by using Python?
  • could you please edit the code to adapt to my needs?
  • could you please show me a step-by-step process?
# Load DICOM files
dicomDataDir = "C:/Users/mario.modesto/Desktop/DICOM/W001"  # input folder with DICOM files
import os 
baboon_skull  = os.path.basename(dicomDataDir)
loadedNodeIDs = []  # this list will contain the list of all loaded node IDs

from DICOMLib import DICOMUtils
with DICOMUtils.TemporaryDICOMDatabase() as db:
  DICOMUtils.importDicom(dicomDataDir, db)
  patientUIDs = db.patients()
  for patientUID in patientUIDs:
    loadedNodeIDs.extend(DICOMUtils.loadPatientByUID(patientUID))

# Display volume rendering
logic = slicer.modules.volumerendering.logic()
volumeNode = slicer.mrmlScene.GetNodeByID('vtkMRMLScalarVolumeNode1')
displayNode = logic.CreateVolumeRenderingDisplayNode()
displayNode.UnRegister(logic)
slicer.mrmlScene.AddNode(displayNode)
volumeNode.AddAndObserveDisplayNodeID(displayNode.GetID())
logic.UpdateDisplayNodeFromVolumeNode(displayNode, volumeNode)

# find the files NodeID
volumeNode = getNode('2: Facial Bones  0.75  H70h')

#create a blank Markup ROI
roiNode = slicer.mrmlScene.AddNewNodeByClass("vtkMRMLMarkupsROINode")

#set the new markup ROI to the dimensions of the volume
cropVolumeParameters = slicer.mrmlScene.AddNewNodeByClass("vtkMRMLCropVolumeParametersNode")
cropVolumeParameters.SetInputVolumeNodeID(volumeNode.GetID())
cropVolumeParameters.SetROINodeID(roiNode.GetID())
slicer.modules.cropvolume.logic().SnapROIToVoxelGrid(cropVolumeParameters)  # optional (rotates the ROI to match the volume axis directions)
slicer.modules.cropvolume.logic().FitROIToInputVolume(cropVolumeParameters)
slicer.mrmlScene.RemoveNode(cropVolumeParameters)

#set the cropping parameters
cropVolumeLogic = slicer.modules.cropvolume.logic()
cropVolumeParameterNode = slicer.vtkMRMLCropVolumeParametersNode()
cropVolumeParameterNode.SetIsotropicResampling(True)

#set the output resolution to 2 millimeters. units in slicer is always in mm. 
cropVolumeParameterNode.SetSpacingScalingConst(2)
cropVolumeParameterNode.SetROINodeID(roiNode.GetID())
cropVolumeParameterNode.SetInputVolumeNodeID(volumeNode.GetID())

#do the cropping
cropVolumeLogic.Apply(cropVolumeParameterNode)

#obtain the nodeID of the cropped volume 
croppedVolume = slicer.mrmlScene.GetNodeByID(cropVolumeParameterNode.GetOutputVolumeNodeID())

# Segmentation
segmentationNode = slicer.mrmlScene.AddNewNodeByClass("vtkMRMLSegmentationNode")
segmentationNode.CreateDefaultDisplayNodes() # only needed for display
segmentationNode.SetReferenceImageGeometryParameterFromVolumeNode(croppedVolume)
addedSegmentID = segmentationNode.GetSegmentation().AddEmptySegment("skull")

# Create segment editor to get access to effects
segmentEditorWidget = slicer.qMRMLSegmentEditorWidget()
segmentEditorWidget.setMRMLScene(slicer.mrmlScene)
segmentEditorNode = slicer.mrmlScene.AddNewNodeByClass("vtkMRMLSegmentEditorNode")
segmentEditorWidget.setMRMLSegmentEditorNode(segmentEditorNode)
segmentEditorWidget.setSegmentationNode(segmentationNode)
segmentEditorWidget.setMasterVolumeNode(croppedVolume)

# Thresholding
segmentEditorWidget.setActiveEffectByName("Threshold")
effect = segmentEditorWidget.activeEffect()
effect.setParameter("MinimumThreshold","115")
effect.setParameter("MaximumThreshold","3071")
effect.self().onApply()

# Clean up
segmentEditorWidget = None
slicer.mrmlScene.RemoveNode(segmentEditorNode)

# Make segmentation results visible in 3D
segmentationNode.CreateClosedSurfaceRepresentation()

# Creatint surface mesh and saving

surfaceMesh = segmentationNode.GetClosedSurfaceInternalRepresentation(addedSegmentID)
writer = vtk.vtkPLYWriter()
writer.SetInputData(surfaceMesh)
writer.SetFileName("C:/Users/mario.modesto/Desktop/DICOM/"+baboon_skull+"_surfaceMesh.ply")
writer.Update()

I would really appreciate your help.

Solution

  • First, for readability reasons, I really would suggest to move the import statements to the top. (Except you have any on demand imports, which is not the case in your posted script)

    Second, screening of the directory and applying the function to each folder should do the job:

    import os
from DICOMLib import DICOMUtils

yourpath = r"<enter your path>/W001-W100"

#walk through DICOM directory
for dir in os.scandir(yourpath):
    # Load DICOM files
    dicomDataDir = dir.path  # path to input folder with DICOM files
    baboon_skull  = dir.name
    loadedNodeIDs = []  # this list will contain the list of all loaded node IDs

    with DICOMUtils.TemporaryDICOMDatabase() as db:
        DICOMUtils.importDicom(dicomDataDir, db)
        patientUIDs = db.patients()
        for patientUID in patientUIDs:
            loadedNodeIDs.extend(DICOMUtils.loadPatientByUID(patientUID))


    < rest of your code >

    writer.SetFileName(fr"{dicomDataDir}_surfaceMesh.ply")
    writer.Update()