<commit_before>import os
import unittest
import vtk
import qt
import ctk
import slicer
from slicer.ScriptedLoadableModule import *
import logging
import dicom
import numpy
from vtk.util import numpy_support

#
# AutoSegmentation
#
class AutoSegmentation(ScriptedLoadableModule):
  def __init__(self, parent):
    ScriptedLoadableModule.__init__(self, parent)
    # TODO make this more human readable by adding spaces
    self.parent.title = "AutoSegmentation"
    self.parent.categories = ["Examples"]
    self.parent.dependencies = []
    self.parent.contributors = ["<NAME> (TU Wien)"]
    self.parent.helpText = """
    This is an example of scripted loadable module bundled in an extension. It performs a simple thresholding on the input volume and optionally captures a screenshot.
    """
    self.parent.acknowledgementText = """
    This file was originally developed by <NAME>, Kitware Inc.
    and <NAME>, Isomics, Inc. and was partially funded by NIH grant 3P41RR013218-12S1.
"""  # replace with organization, grant and thanks.

#
# AutoSegmentationWidget
#
class AutoSegmentationWidget(ScriptedLoadableModuleWidget):
  """Uses ScriptedLoadableModuleWidget base class, available at:
  https://github.com/Slicer/Slicer/blob/master/Base/Python/slicer/ScriptedLoadableModule.py
  """
  def setup(self):
    # gui for testing and reloading, /todo remove at end
    ScriptedLoadableModuleWidget.setup(self)

    # Collapsible button
    self.collapsibleButton = ctk.ctkCollapsibleButton()
    self.collapsibleButton.text = "AutoSegmentation"
    self.layout.addWidget(self.collapsibleButton)

    # Layout for path input and process button
    self.formLayout = qt.QFormLayout(self.collapsibleButton)

    # Path input for dicom data to analyze
    self.inputPath = qt.QFileDialog()
    self.inputPath.setFileMode(qt.QFileDialog.Directory)

    #self.inputPathCtk = ctk.ctkDICOMAppWidget()
    
    self.formLayout.addWidget(self.inputPath)
    #self.formLayout.addWidget(self.inputPathCtk)

    # add vertical spacer
    self.layout.addStretch(1)

    # connect directory widget with function
    self.inputPath.connect('accepted()', self.createLogic)
    #self.inputPathCtk.connect('onFileIndexed(const QString &filePath)', self.onModelSelectedInput)

  def createLogic(self):
    pathToDICOM = self.inputPath.directory().absolutePath()
    self.logic = AutoSegmentationLogic(pathToDICOM, 0.75, 0.2, -0.2)

#
# AutoSegmentationLogic
#
class AutoSegmentationLogic(ScriptedLoadableModuleLogic):

  def __init__(self, pathToDICOM, minTreshold, curve1Minimum, curve3Maximum):
    self.pathToDICOM = pathToDICOM
    self.minTreshold = minTreshold
    self.curve1Minimum = curve1Minimum
    self.curve3Maximum = curve3Maximum
    logging.info("logic created")
    logging.info("parameters:")
    logging.info(pathToDICOM)
    logging.info(minTreshold)
    logging.info(curve1Minimum)
    logging.info(curve3Maximum)
    
    self.dicomDataNumpyArrays = self.readData()
    self.initialRiseArray = self.calcInitialRise()
    self.slopeArray = self.calcSlope()
    logging.info("initial rise and slope calculated")

    self.targetVoxels = self.getTargetedVoxels()
    #self.persistenceVoxels = self.getPersistanceVoxels()
    #self.plateauVoxels = self.getPlateauVoxels()
    #self.washoutVoxels = self.getWashoutVoxels()

    VTKtargetVoxels = numpy_support.numpy_to_vtk(self.targetVoxels.ravel(), deep=True, array_type=vtk.VTK_FLOAT)
    VTKtargetVoxels.contour


  def readData(self):
    filesDCM = []
    for dirName, subdirList, fileList in os.walk(self.pathToDICOM):
      for filename in fileList:
        if ".dcm" in filename.lower():
          filesDCM.append(os.path.join(dirName, filename))


    #load all dicom data 
    dicomData = []
    for files in filesDCM:
      dicomData.append(dicom.read_file(files)) 

    #sort dicom data by instance number
    dicomData.sort(key=lambda dicomData: dicomData.InstanceNumber)

    #separate data to their contrast volumes
    dicomDataContrastVolumes = []
    firstContrastVolume = []
    dicomDataContrastVolumes.append(firstContrastVolume)
    contrastVolumeIndexHelper = 0
    for i in range(0,len(dicomData)-1):
      if dicomData[i].SliceLocation < dicomData[i+1].SliceLocation:
        dicomDataContrastVolumes[contrastVolumeIndexHelper].append(dicomData[i])
      else:
        dicomDataContrastVolumes[contrastVolumeIndexHelper].append(dicomData[i])
        contrastVolumeIndexHelper = contrastVolumeIndexHelper+1
        nextContrastVolume = []
        dicomDataContrastVolumes.append(nextContrastVolume)
    #assign last element
    dicomDataContrastVolumes[contrastVolumeIndexHelper].append(dicomData[-1])

    logging.info("Dicom data seperated into %i sets.", len(dicomDataContrastVolumes))

    #convert sets into numpy Arrays
    dicomNumpyArrays = []
    for dicomSeries in dicomDataContrastVolumes:
      dicomNumpyArrays.append(self.createNumpyArray(dicomSeries))

    logging.info("3D Numpy Arrays created")

    return dicomNumpyArrays

  # Computes percentage increase from baseline (pre-contrast) at each voxel for each volume as numpy arrays.
  def calcInitialRise(self):
    # Initial Rise at each voxel (percentage increase from pre-contrast to first post-contrast)
    return ((self.dicomDataNumpyArrays[1]).__truediv__(self.dicomDataNumpyArrays[0]+1.0))-1.0

  def calcSlope(self):
    # Compute slope at each voxel from first to fourth volume to determine curve type
    return (self.dicomDataNumpyArrays[-1] - self.dicomDataNumpyArrays[1]).__truediv__(self.dicomDataNumpyArrays[1]+1.0)

  def getTargetedVoxels(self):
    targetVoxels = (self.initialRiseArray > self.minTreshold) & (self.dicomDataNumpyArrays[0] > 100)
    return targetVoxels.astype(bool)

  def getPersistanceVoxels(self):
    persistenceVoxels = numpy.where((self.slopeArray > self.curve1Minimum) & (self.targetVoxels))
    return persistenceVoxels.astype(bool)

  def getPlateauVoxels(self):
    plateuaVoxels = numpy.where( (self.slopeArray > self.curve3Maximum) & (self.slopeArray < self.curve1Minimum) & (self.targetVoxels) )
    return plateauVoxels.astype(bool)

  def getWashoutVoxels(self):
    washoutVoxels = numpy.where((self.slopeArray < self.curve3Maximum) & (self.targetVoxels))
    return washoutVoxels.astype(bool)



  # def arrayProcessing(self):
  #   # Create Boolean array, target_voxels, with target voxel indices highlighted as True 
  #   # Assigns color to SegmentCAD Label map index if corresponding slope condition is satisfied where target_voxel is True 
      
  #   target_voxels = (self.nodeArrayInitialRise > self.minTreshold) & (self.dicomNumpyArrays[0] > 100)
  
  #   # yellow (Plateau Slope)
  #   self.nodeArraySegmentCADLabel[numpy.where( (self.slopeArray1_4 > -0.2) & (self.slopeArray1_4 < 0.2) & (target_voxels) )] = 291
    
  #   # blue (slope of curve1 min = 0.2(default), Persistent Slope)
  #   self.nodeArraySegmentCADLabel[numpy.where((self.slopeArray1_4 > 0.2) & (target_voxels))] = 306
    
  #   # red (slope of curve3 max = -0.2(default), Washout Slope )
  #   self.nodeArraySegmentCADLabel[numpy.where((self.slopeArray1_4 < -0.2) & (target_voxels))] = 32

<commit_msg>numpy to stl<commit_after>import os
import unittest
import vtk
import qt
import ctk
import slicer
from slicer.ScriptedLoadableModule import *
import logging
import dicom
import numpy
from vtk.util import numpy_support

#
# AutoSegmentation
#
class AutoSegmentation(ScriptedLoadableModule):
  def __init__(self, parent):
    ScriptedLoadableModule.__init__(self, parent)
    # TODO make this more human readable by adding spaces
    self.parent.title = "AutoSegmentation"
    self.parent.categories = ["Examples"]
    self.parent.dependencies = []
    self.parent.contributors = ["<NAME> (TU Wien)"]
    self.parent.helpText = """
    This is an example of scripted loadable module bundled in an extension. It performs a simple thresholding on the input volume and optionally captures a screenshot.
    """
    self.parent.acknowledgementText = """
    This file was originally developed by <NAME>, Kitware Inc.
    and <NAME>, Isomics, Inc. and was partially funded by NIH grant 3P41RR013218-12S1.
"""  # replace with organization, grant and thanks.

#
# AutoSegmentationWidget
#
class AutoSegmentationWidget(ScriptedLoadableModuleWidget):
  """Uses ScriptedLoadableModuleWidget base class, available at:
  https://github.com/Slicer/Slicer/blob/master/Base/Python/slicer/ScriptedLoadableModule.py
  """
  def setup(self):
    # gui for testing and reloading, /todo remove at end
    ScriptedLoadableModuleWidget.setup(self)

    # Collapsible button
    self.collapsibleButton = ctk.ctkCollapsibleButton()
    self.collapsibleButton.text = "AutoSegmentation"
    self.layout.addWidget(self.collapsibleButton)

    # Layout for path input and process button
    self.formLayout = qt.QFormLayout(self.collapsibleButton)

    # Path input for dicom data to analyze
    self.inputPath = qt.QFileDialog()
    self.inputPath.setFileMode(qt.QFileDialog.Directory)

    #self.inputPathCtk = ctk.ctkDICOMAppWidget()
    
    self.formLayout.addWidget(self.inputPath)
    #self.formLayout.addWidget(self.inputPathCtk)

    # add vertical spacer
    self.layout.addStretch(1)

    # connect directory widget with function
    self.inputPath.connect('accepted()', self.createLogic)
    #self.inputPathCtk.connect('onFileIndexed(const QString &filePath)', self.onModelSelectedInput)

  def createLogic(self):
    pathToDICOM = self.inputPath.directory().absolutePath()
    self.logic = AutoSegmentationLogic(pathToDICOM, 0.75, 0.2, -0.2)

#
# AutoSegmentationLogic
#
class AutoSegmentationLogic(ScriptedLoadableModuleLogic):

  def __init__(self, pathToDICOM, minTreshold, curve1Minimum, curve3Maximum):
    self.pathToDICOM = pathToDICOM
    self.minTreshold = minTreshold
    self.curve1Minimum = curve1Minimum
    self.curve3Maximum = curve3Maximum
    logging.info("logic created")
    logging.info("parameters:")
    logging.info(pathToDICOM)
    logging.info(minTreshold)
    logging.info(curve1Minimum)
    logging.info(curve3Maximum)
    
    self.dicomDataNumpyArrays = self.readData()
    self.initialRiseArray = self.calcInitialRise()
    self.slopeArray = self.calcSlope()
    logging.info("initial rise and slope calculated")

    #boolean array with targeted voxels
    self.targetVoxels = self.getTargetedVoxels()
    
    self.persistenceVoxels = self.getPersistanceVoxels()

    #self.plateauVoxels = self.getPlateauVoxels()
    #self.washoutVoxels = self.getWashoutVoxels()

    #convert numpy to vtkImageData
    VTKTargetVoxelsImageImport =  vtk.vtkImageImport()
    array_string = self.persistenceVoxels.tostring()

    VTKTargetVoxelsImageImport.CopyImportVoidPointer(array_string, len(array_string))
    VTKTargetVoxelsImageImport.SetDataScalarTypeToUnsignedChar()
    VTKTargetVoxelsImageImport.SetNumberOfScalarComponents(1)

    w, d, h = self.persistenceVoxels.shape
    VTKTargetVoxelsImageImport.SetDataExtent(0,w-1,0,d-1,0,h-1)
    VTKTargetVoxelsImageImport.SetWholeExtent(0,w-1,0,d-1,0,h-1)

    VTKTargetVoxelsImageImport.SetDataSpacing(0.5,0.5,2)

    logging.info("imageImporter set up")

    # alphaChannelFunc = vtk.vtkPiecewiseFunction()
    # alphaChannelFunc.AddPoint(0, 0.0)
    # alphaChannelFunc.AddPoint(1, 1)

    # colorFunc = vtk.vtkColorTransferFunction()
    # colorFunc.AddRGBPoint(1, 1.0, 0.0, 0.0)

    # volumeProperty = vtk.vtkVolumeProperty()
    # volumeProperty.SetColor(colorFunc)
    # volumeProperty.SetScalarOpacity(alphaChannelFunc)

    # compositeFunction = vtk.vtkVolumeRayCastCompositeFunction()

    # volumeMapper = vtk.vtkVolumeRayCastMapper()
    # volumeMapper.SetVolumeRayCastFunction(compositeFunction)
    # volumeMapper.SetInputConnection(dataImporter.GetOutputPort())

    # volume = vtk.vtkVolume()
    # volume.SetMapper(volumeMapper)
    # volume.SetProperty(volumeProperty)


    threshold = vtk.vtkImageThreshold()
    threshold.SetInputConnection(VTKTargetVoxelsImageImport.GetOutputPort())
    threshold.ThresholdByLower(1)
    threshold.ReplaceInOn()
    threshold.SetInValue(0)
    threshold.SetOutValue(1)
    threshold.Update()
    
    dmc = vtk.vtkDiscreteMarchingCubes()
    dmc.SetInputConnection(threshold.GetOutputPort())
    dmc.GenerateValues(1,1,1)
    dmc.Update()

    logging.info("marching curbes applied")

    writer = vtk.vtkSTLWriter()
    writer.SetInputConnection(dmc.GetOutputPort())
    writer.SetFileTypeToBinary()
    writer.SetFileName("test.stl")
    writer.Write()

    logging.info("finished writing to file")

  def readData(self):
    filesDCM = []
    for dirName, subdirList, fileList in os.walk(self.pathToDICOM):
      for filename in fileList:
        if ".dcm" in filename.lower():
          filesDCM.append(os.path.join(dirName, filename))


    #load all dicom data 
    dicomData = []
    for files in filesDCM:
      dicomData.append(dicom.read_file(files)) 

    #sort dicom data by instance number
    dicomData.sort(key=lambda dicomData: dicomData.InstanceNumber)

    #separate data to their contrast volumes
    dicomDataContrastVolumes = []
    firstContrastVolume = []
    dicomDataContrastVolumes.append(firstContrastVolume)
    contrastVolumeIndexHelper = 0
    for i in range(0,len(dicomData)-1):
      if dicomData[i].SliceLocation < dicomData[i+1].SliceLocation:
        dicomDataContrastVolumes[contrastVolumeIndexHelper].append(dicomData[i])
      else:
        dicomDataContrastVolumes[contrastVolumeIndexHelper].append(dicomData[i])
        contrastVolumeIndexHelper = contrastVolumeIndexHelper+1
        nextContrastVolume = []
        dicomDataContrastVolumes.append(nextContrastVolume)
    #assign last element
    dicomDataContrastVolumes[contrastVolumeIndexHelper].append(dicomData[-1])

    logging.info("Dicom data seperated into %i sets.", len(dicomDataContrastVolumes))

    #convert sets into numpy Arrays
    dicomNumpyArrays = []
    for dicomSeries in dicomDataContrastVolumes:
      dicomNumpyArrays.append(self.createNumpyArray(dicomSeries))

    logging.info("3D Numpy Arrays created")

    return dicomNumpyArrays

  # Computes percentage increase from baseline (pre-contrast) at each voxel for each volume as numpy arrays.
  def calcInitialRise(self):
    # Initial Rise at each voxel (percentage increase from pre-contrast to first post-contrast)
    return ((self.dicomDataNumpyArrays[1]).__truediv__(self.dicomDataNumpyArrays[0]+1.0))-1.0

  def calcSlope(self):
    # Compute slope at each voxel from first to fourth volume to determine curve type
    return (self.dicomDataNumpyArrays[-1] - self.dicomDataNumpyArrays[1]).__truediv__(self.dicomDataNumpyArrays[1]+1.0)

  def getTargetedVoxels(self):
    targetVoxels = numpy.zeros((256,256,60))
    for x in range(0,255):
      for y in range(0,255):
        for z in range(0,59):
          if (self.initialRiseArray[x,y,z] > self.minTreshold):
            targetVoxels[x,y,z] = 100
    return targetVoxels

    #targetVoxels = (self.initialRiseArray > self.minTreshold) & (self.dicomDataNumpyArrays[0] > 100)
    #return targetVoxels

  def getPersistanceVoxels(self):
    persistenceVoxels = numpy.zeros((256,256,60))
    for x in range(0,255):
      for y in range(0,255):
        for z in range(0,59):
          if (self.slopeArray[x,y,z] > self.curve3Maximum) & (self.slopeArray[x,y,z] < self.curve1Minimum) & (self.targetVoxels[x,y,z] > 99):
            persistenceVoxels[x,y,z] = 100
    return persistenceVoxels

    #persistenceVoxels = numpy.where((self.slopeArray > self.curve1Minimum) & (self.targetVoxels))
    #return persistenceVoxels

  def getPlateauVoxels(self):
    plateuaVoxels = numpy.where( (self.slopeArray > self.curve3Maximum) & (self.slopeArray < self.curve1Minimum) & (self.targetVoxels) )
    return plateauVoxels

  def getWashoutVoxels(self):
    washoutVoxels = numpy.where((self.slopeArray < self.curve3Maximum) & (self.targetVoxels))
    return washoutVoxels



  # def arrayProcessing(self):
  #   # Create Boolean array, target_voxels, with target voxel indices highlighted as True 
  #   # Assigns color to SegmentCAD Label map index if corresponding slope condition is satisfied where target_voxel is True 
      
  #   target_voxels = (self.nodeArrayInitialRise > self.minTreshold) & (self.dicomNumpyArrays[0] > 100)
  
  #   # yellow (Plateau Slope)
  #   self.nodeArraySegmentCADLabel[numpy.where( (self.slopeArray1_4 > -0.2) & (self.slopeArray1_4 < 0.2) & (target_voxels) )] = 291
    
  #   # blue (slope of curve1 min = 0.2(default), Persistent Slope)
  #   self.nodeArraySegmentCADLabel[numpy.where((self.slopeArray1_4 > 0.2) & (target_voxels))] = 306
    
  #   # red (slope of curve3 max = -0.2(default), Washout Slope )
  #   self.nodeArraySegmentCADLabel[numpy.where((self.slopeArray1_4 < -0.2) & (target_voxels))] = 32


{'fraction_non_alphanumeric': 0.09870607787274453, 'fraction_numerical': 0.015432098765432098, 'mean_word_length': 4.227738132174992, 'pattern_counts': {'":': 0, '<': 20, '<?xml version=': 0, '>': 28, 'https://': 2, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '664275', 'n_tokens_mistral': 5790, 'n_tokens_neox': 5502, 'n_words': 1074}