diff --git a/README.md b/README.md
index 6b66051..a9ba348 100644
--- a/README.md
+++ b/README.md
@@ -119,4 +119,9 @@ For instructions on how to run TOPAS, please see the user guide of OpenTOPAS (op
 > We recommend to run several of the PHSP provided by Varian and combine the results for better statistics. To run the same simulation with different PHSP files, in the _Main.txt" TOPAS file, replace the name of the current path to the PHSP file in the command named s:So/phsp/PhaseSpaceFileName_ by the new PHSP filename.
 
 > [!WARNING]
-> The calculation of dose distributions in TOPAS relies on a calibration factor, dependent on the beam energy, that relates the dose per primary history in the simulation to the dose per MU (i.e., calibration factor of the system). This calibration factor may need to be modified if the calibration method employed at your institution is different from the implemented by default in TPS2TOPAS (see [https://doi.org/10.1016/j.ejmp.2024.104485](https://doi.org/10.1016/j.ejmp.2024.104485)) or the phsp files used are different from the recommended ones. Note that the calibration factors included in this version of TPS2TOPAS are for 6MV and 10MV (**not FFF**) beams. 
+> The calculation of dose distributions in TOPAS relies on a calibration factor, dependent on the beam energy, that relates the dose per primary history in the simulation to the dose per MU (i.e., calibration factor of the system). This calibration factor may need to be modified if the calibration method employed at your institution is different from the implemented by default in TPS2TOPAS (see [https://doi.org/10.1016/j.ejmp.2024.104485](https://doi.org/10.1016/j.ejmp.2024.104485)) or the phsp files used are different from the recommended ones. Note that the calibration factors included in this version of TPS2TOPAS are for 6MV and 10MV (**not FFF**) beams.
+> 
+
+
+## Updates by Justin Stanton-Sheetz
+Updated to add PhaseSpace scoring and also WSL compatability for those storing files in different drives than running TOPAS
diff --git a/TPS2TOPAS.py b/TPS2TOPAS.py
index 4560393..7ad662a 100644
--- a/TPS2TOPAS.py
+++ b/TPS2TOPAS.py
@@ -121,8 +121,16 @@ def main():
 
     # Create the project directory
     os.system('mkdir %s' %DATA["project_name"])
-    os.system('mkdir %s/output' %DATA["project_name"])
-    os.system('cp HUtoMaterialSchneider.txt %s' %DATA["project_name"])
+    if DATA["OS"] == "linux/macos":
+        os.system('mkdir %s/output' %DATA["project_name"]) 
+        os.system('cp HUtoMaterialSchneider.txt %s' %DATA["project_name"]) 
+    elif DATA["OS"] == "windows":
+        os.system('mkdir %s\\output' %DATA["project_name"]) 
+        os.system('copy HUtoMaterialSchneider.txt %s' %DATA["project_name"]) #Changed cp to copy for windows
+    else:
+        print('OS selection error. Linux/MacOS assumed')
+        os.system('mkdir %s/output' %DATA["project_name"]) 
+        os.system('cp HUtoMaterialSchneider.txt %s' %DATA["project_name"]) 
 
     # Retrieve data from files exported from TPS
     CT_DATA = RetrieveCTData(DATA)
@@ -149,4 +157,4 @@ def main():
 ############################################################################################################################################
 
 if __name__ == "__main__":
-    main()
\ No newline at end of file
+    main()
diff --git a/TPS2TOPAS_Track_Changes.md b/TPS2TOPAS_Track_Changes.md
new file mode 100644
index 0000000..3353f26
--- /dev/null
+++ b/TPS2TOPAS_Track_Changes.md
@@ -0,0 +1,58 @@
+track changes to TPS2TOPAS
+
+Changes implemented by Justin Stanton-Sheetz 1/14/2026
+
+
+
+
+
+\####################### Changes/additions ###############################
+
+**inputfile\_template:**
+
+ 	Added OS to the template
+
+**input\_handling:**
+
+ 	Added "PhaseSpace" to the options for Scoring Quantity dropdown
+
+ 	Added OS option to inputFile mode
+
+ 	Adjusted numbers in inputfile mode to fit OS choice in
+
+ 	Added "ASCII" to output format
+
+ 	Added if statement to score "PhaseSpace" 
+
+ 	Added a selection button for OS. 2 choices Linux/iOS or Windows. This will allow for some specific changes built for the wsl environment to be 			able to access files stored in windows drives
+
+ 	Added ASCII and PhaseSpace to the checks for errors and added input file length check dependent on if PhaseSpace is chosen (13 inputs instead of 		12 so you can choose ROIs to score)
+
+ 	For inputfile mode: Added reading, converting, and checking of line 13. This checks that the ROIs listed are within the structure set and TOPAS 		can read them (lines 122-170)
+
+ 	For GUI Mode: If scoring quantity is phasespace, user will be prompted with a second gui that allows for ROI selection from the structure set.
+
+ 		Added ChooseROIs function to prompt said gui
+
+ 
+
+**TPS2TOPAS:**
+
+ 	Changed cp to copy for windows (line 125) with if statement
+
+**write\_PCF:**
+
+ 	changed theta view angle for phase space to 0 degrees. Honestly I think this should be the case for all, but I don't want to mess up what is 			already here if there is a good reason for 89.9 degrees.
+
+ 	added if statement for file path in TOPAS file. If OS is selected to be windows it will change the file path from a drive (like C:/... or H:/... 		to something WSL can read like "/mnt/c/..."
+
+ 	Added if statement to adjust the scoring for PhaseSpace if selected. This will make scoring for all separate structures selected
+
+ 	Added If statement to visualization file so that ROIs are colored on DICOM. They all are set to yellow. User can change the topas file if desired.
+
+
+
+**plan\_data:**
+
+ 	Added a warning if a beam has 0 MU
+
diff --git a/input_handling.py b/input_handling.py
index 66eb244..a162fea 100644
--- a/input_handling.py
+++ b/input_handling.py
@@ -46,12 +46,14 @@ def EvaluateCorrectInputData(DATA):
         DATA["scoring_quantity"] = "DoseToWater"
     elif DATA["scoring_quantity"] == 'dosetomedium':
         DATA["scoring_quantity"] = "DoseToMedium"
+    elif DATA["scoring_quantity"] == 'phasespace':
+        DATA["scoring_quantity"] = "PhaseSpace"
     else:
         print("--- WARNING: No valid scoring quantity selected (%s) - DoseToMedium used by default" %DATA["scoring_quantity"])
         DATA["scoring_quantity"] = "DoseToMedium"
         warnings +=1
     
-    if DATA["output_format"] != "binary" and DATA["output_format"] != "csv" and DATA["output_format"] != "dicom" and DATA["output_format"] != "root" and DATA["output_format"] != "xml":
+    if DATA["output_format"] != "binary" and DATA["output_format"] != "csv" and DATA["output_format"] != "dicom" and DATA["output_format"] != "root" and DATA["output_format"] != "ascii" and DATA["output_format"] != "xml":
         print("--- WARNING: No valid Output Format selected (%s) - binary format used by default" %DATA["output_format"])
         DATA["output_format"] = "binary"
         warnings +=1
@@ -73,6 +75,8 @@ def EvaluateCorrectInputData(DATA):
     print("     (9) Scoring quantity:       %s" % DATA["scoring_quantity"])
     print("    (10) Output file name:       %s" % DATA["output_file"])
     print("    (11) Output format:          %s" % DATA["output_format"])
+    if DATA["scoring_quantity"] == "PhaseSpace":
+        print("    (12) ROIs to score:          %s" % DATA["ROIs_to_score"])
 
     return DATA
 
@@ -85,27 +89,86 @@ def EvaluateCorrectInputData(DATA):
 def InputDataInputFileMode(inputFile):
     inputInfo = open(inputFile,'r').read().split('\n')
     inputInfo = list(filter(None, inputInfo)) # filter empty lines
-    if len(inputInfo) != 11:
+    if len(inputInfo) != 13 and inputInfo[9].lower() == 'phasespace':
         print("######")
-        print('ERROR! %s arguments found in input file. %s arguments required' %(len(inputInfo),11))
+        print('ERROR! %s arguments found in input file. %s arguments required' %(len(inputInfo),13))
+        print("######")
+        printHelp()
+        exit(0)
+    elif len(inputInfo) != 12 and inputInfo[9].lower() != 'phasespace':
+        print("######")
+        print('ERROR! %s arguments found in input file. %s arguments required' %(len(inputInfo),12))
         print("######")
         printHelp()
         exit(0)
     
     DATA = {}
+    ROI_List = {}
     project_name_temp = inputInfo[0]
     DATA["project_name"]        = project_name_temp.replace(" ", "_")
-    DATA["dicom_dirname"]       = inputInfo[1]
-    DATA["RS_filename"]         = inputInfo[2] 
-    DATA["RD_filename"]         = inputInfo[3]
-    DATA["RP_filename"]         = inputInfo[4]
-    DATA["phsp_filename"]       = inputInfo[5]
-    DATA["MLC_model"]           = inputInfo[6].lower()
-    DATA["multipleUse"]         = inputInfo[7]
-    DATA["scoring_quantity"]    = inputInfo[8].lower()
-    output_file_temp = inputInfo[9]
+    DATA["OS"]                  = inputInfo[1]
+    DATA["dicom_dirname"]       = inputInfo[2]
+    DATA["RS_filename"]         = inputInfo[3] 
+    DATA["RD_filename"]         = inputInfo[4]
+    DATA["RP_filename"]         = inputInfo[5]
+    DATA["phsp_filename"]       = inputInfo[6]
+    DATA["MLC_model"]           = inputInfo[7].lower()
+    DATA["multipleUse"]         = inputInfo[8]
+    DATA["scoring_quantity"]    = inputInfo[9].lower() 
+    output_file_temp = inputInfo[10]
     DATA["output_file"]         = output_file_temp.replace(" ", "_")
-    DATA["output_format"]    = inputInfo[10].lower()
+    DATA["output_format"]    = inputInfo[11].lower()
+
+
+    if DATA["scoring_quantity"] == 'phasespace':
+        DATA["ROIs_to_score"]     = inputInfo[12].split(',')
+        for roi in DATA["ROIs_to_score"]:
+            newroi = roi.strip()
+            for letter in range(len(newroi)):
+                if newroi[letter] == ' ':
+                    newroi = newroi.replace(' ','_')
+                if '.' in newroi:
+                    newroi = newroi.replace('.','_')
+                if '+' in newroi:
+                    newroi = newroi.replace('+','_')
+            DATA["ROIs_to_score"][DATA["ROIs_to_score"].index(roi)] = newroi
+    
+        ds_RS = pydicom.dcmread(DATA["RS_filename"])
+        RoiGeometries = ds_RS.RTROIObservationsSequence
+        ROI_List = {}
+        for roi in RoiGeometries:
+            for region in ds_RS.StructureSetROISequence:
+                if region.ROINumber == roi.ReferencedROINumber:
+                    roiName = region.ROIName
+                    for letter in range(len(roiName)):
+                        if roiName[letter] == ' ':
+                            
+                            roiName = roiName.replace(' ','_')
+                        if '.' in roiName:
+                            roiName = roiName.replace('.','_')
+                        if '+' in roiName:
+                            roiName = roiName.replace('+','_')
+                    if len(roiName) > 16:
+                        roiName = roiName[:16]
+                        print("--- WARNING: The ROI %s has a name longer than 16 characters. The name will be truncated to %s. Please rename in structure set" %(roiName,roiName[:16]))
+                    ROI_List[region.ROINumber] = roiName
+
+    
+        for roi in DATA["ROIs_to_score"]:
+            score = 0 
+            for region in ds_RS.StructureSetROISequence:
+                if roi == ROI_List[region.ROINumber]:
+                    score = score + 1
+            if score == 0:
+                print("######")
+                print('ERROR! The ROI %s specified in the input file do not match the ROIs found in the DICOM-structure file'%roi)
+                print("######")
+                exit(0)
+            elif score > 1:
+                print("######")
+                print('ERROR! The ROI %s specified in the input file match multiple ROIs found in the DICOM-structure file'%roi)
+                print("######")
+                exit(0)
 
     return EvaluateCorrectInputData(DATA)
 
@@ -170,6 +233,21 @@ def browse_phsp_button():
     ################
     ################
     i += 1
+    label13=StringVar()
+    label13.set("Operating System")
+    lbl13 = Label(master=root, textvariable=label13)
+    lbl13.grid(row=i, column=1)
+    OperatingSys = StringVar()
+    OperatingSys.set('select')
+    options = ["Linux/MacOS","Windows"]
+    OSentry = OptionMenu(root , OperatingSys , *options)
+    OSentry.grid(row=i, column=2)
+
+    lbl13 = Label()
+    lbl13.grid(row=i+1, column=1)
+    ################
+    ################
+    i += 1
     label1=StringVar()
     label1.set('DICOM directory')
     lbl1 = Label(master=root, textvariable=label1) 
@@ -275,7 +353,7 @@ def browse_phsp_button():
     lbl8.grid(row=i, column=1)
     Scoring = StringVar()
     Scoring.set('select')
-    options = ["DoseToWater", "DoseToMedium"]
+    options = ["DoseToWater", "DoseToMedium", "PhaseSpace"]
     scoringentry = OptionMenu(root , Scoring , *options)
     scoringentry.grid(row=i, column=2)
 
@@ -304,7 +382,7 @@ def browse_phsp_button():
     lbl10.grid(row=i, column=1)
     outputFormat = StringVar()
     outputFormat.set('select')
-    options = ["binary", "DICOM", "csv", "root", "xml"]
+    options = ["binary", "DICOM", "csv", "root", "xml", "ASCII"]
     Outpytentry = OptionMenu(root , outputFormat , *options)
     Outpytentry.grid(row=i, column=2)
 
@@ -330,6 +408,7 @@ def browse_phsp_button():
     DATA = {}
     project_name_temp = projectName.get()
     DATA["project_name"]        = project_name_temp.replace(" ", "_")
+    DATA["OS"]                  = OperatingSys.get().lower()
     DATA["dicom_dirname"]       = dicomDirectoryName.get()
     DATA["RS_filename"]         = dicomStructureFileName.get()
     DATA["RD_filename"]         = dicomDoseFileName.get()
@@ -341,7 +420,65 @@ def browse_phsp_button():
     output_file_temp = outFileName.get()
     DATA["output_file"]         = output_file_temp.replace(" ", "_")
     DATA["output_format"]    = outputFormat.get().lower()
+    
+    
+    if DATA["scoring_quantity"] == "phasespace":
+        ds_RS = pydicom.dcmread(DATA["RS_filename"])
+        RoiGeometries = ds_RS.RTROIObservationsSequence
+        ROI_List = {}
+        for roi in RoiGeometries:
+            for region in ds_RS.StructureSetROISequence:
+                if region.ROINumber == roi.ReferencedROINumber:
+                    roiName = region.ROIName
+                    for letter in range(len(roiName)):
+                        if roiName[letter] == ' ':
+                            
+                            roiName = roiName.replace(' ','_')
+                        if '.' in roiName:
+                            roiName = roiName.replace('.','_')
+                        if '+' in roiName:
+                            roiName = roiName.replace('+','_')
+                    if len(roiName) > 16:
+                        roiName = roiName[:16]
+                    ROI_List[region.ROINumber] = roiName
+
+        DATA["ROIs_to_score"] = []
+        DATA["ROIs_to_score"] = ChooseROIs(ROI_List)
 
     return EvaluateCorrectInputData(DATA)
+############################################################################################################################################
+############################################################################################################################################
+############################################Choose ROIs to score from DICOM-structure file##################################################
+############################################################################################################################################
+############################################################################################################################################
+
+def ChooseROIs(ROI_List):
+    selected_ROIs = []
+    ROI_List = dict(sorted(ROI_List.items()))
+    def _quit_ROI():
+        nonlocal selected_ROIs
+        selected_ROIs = [ROI_List[i+1] for i in listbox.curselection()]
+        roi_root.quit()
+        roi_root.destroy() 
+        
+    global roi_root; roi_root = Tk()
+    roi_root.title("Select ROIs to score")
+    roi_root.geometry('180x200')
+
+    listbox = Listbox(roi_root, width=40,height=10, selectmode=MULTIPLE)
+    idx = 0
+    for key in ROI_List:
+        
+        listbox.insert(idx, ROI_List[key])
+        idx+=1
+    # print(ROI_List)
+    # print(selected_ROIs)
+    listbox.pack()
+
+    exitButton = Button(text="Confirm selection", command=_quit_ROI)
+    exitButton.pack()
+
+    mainloop()
 
+    return selected_ROIs
 
diff --git a/inputfile_Example.txt b/inputfile_Example.txt
new file mode 100644
index 0000000..0e95104
--- /dev/null
+++ b/inputfile_Example.txt
@@ -0,0 +1,13 @@
+Test
+windows
+C:/DICOM_dir_file_path/
+C:/DICOM_dir_file_path/RTst_file.dcm
+C:/DICOM_dir_file_path/RT_dose_file.dcm
+C:/DICOM_dir_file_path/RT_plan_file.dcm
+C:/phase_space_file_path/phasespace_file.dcm
+genericHD
+1
+PhaseSpace
+output
+ASCII
+GTV04_L_TEMP, GTV05_R_CER, GTV03_L_FR
diff --git a/inputfile_template.txt b/inputfile_template.txt
index e7fce77..d319b41 100644
--- a/inputfile_template.txt
+++ b/inputfile_template.txt
@@ -1,4 +1,5 @@
 Project name
+Operating System: "linux/macos" or "windows"
 path to DICOM directory
 path to DICOM-structure file
 path to DICOM-dose file
@@ -8,4 +9,5 @@ MLC model: "generic", "genericHD", "Varian" or "VarianHD"
 Geometrical particle splitting factor (integer number)
 Scoring quantity: "DoseToWater" or "DoseToMedium
 Output file name
-Output file format: "binary", "DICOM", "csv", "root" or "xml"
\ No newline at end of file
+Output file format: "binary", "DICOM", "csv", "root", "ASCII", or "xml"
+ROI_Names (You can score multiple here, separated by commas)
diff --git a/plan_data.py b/plan_data.py
index 85f9520..42a4321 100644
--- a/plan_data.py
+++ b/plan_data.py
@@ -295,7 +295,9 @@ def RetrievePlanData(DATA):
                         TsMLCX2[key].append(MLCX2[key])
                     
                     dt += CONSTANTS.dt 
-    
+        else:
+            print("--- WARNING: Beam %s has 0 MU and will be skipped." %beam.BeamName) 
+            
     phspChunk = int(CONSTANTS.primaryHistories / len(TsWeights))
     for i in range(len(TsWeights)):
         TsPrimaries.append(phspChunk)
diff --git a/write_PCF.py b/write_PCF.py
index b2ba919..a24e7d7 100644
--- a/write_PCF.py
+++ b/write_PCF.py
@@ -1,4 +1,5 @@
 import pathlib
+import wslPath
 
 ############################################################################################################################################
 ############################################################################################################################################
@@ -131,8 +132,23 @@ def WritePlanParameterFile(DATA,CT_DATA,PLAN_DATA):
     parFile.write('includeFile = Isocenter_TimeFeatures.txt\n')
     parFile.write('includeFile = CalibrationFactor_TimeFeatures.txt\n')
     parFile.write('\n')
-    parFile.write('s:RS/DicomDirectory = \"%s\"\n' % DATA["dicom_dirname"]) 
-    parFile.write('s:RS/DicomDoseFileName = \"%s\"\n' % DATA["RD_filename"]) 
+    if DATA["OS"] == "linux/macos":
+        parFile.write('s:RS/DicomDirectory = \"%s\"\n' % DATA["dicom_dirname"]) 
+        parFile.write('s:RS/DicomDoseFileName = \"%s\"\n' % DATA["RD_filename"]) 
+    elif DATA["OS"] == "windows":
+        windows_dicom_dir = DATA["dicom_dirname"]
+        windows_dicom_dir = windows_dicom_dir.replace('/', '\\')
+        dicom_dir = wslPath.to_posix(windows_dicom_dir)
+        windows_dose_file = DATA["RD_filename"]
+        windows_dose_file = windows_dose_file.replace('/','\\')
+        dose_file = wslPath.to_posix(windows_dose_file)
+        parFile.write('s:RS/DicomDirectory = \"%s\"\n' % dicom_dir) 
+        parFile.write('s:RS/DicomDoseFileName = \"%s\"\n' % dose_file) 
+        
+    else:
+        print('OS selection error. Dicom and dose file names left blank')
+        parFile.write('s:RS/DicomDirectory = \"%s\"\n') 
+        parFile.write('s:RS/DicomDoseFileName = \"%s\"\n') 
     parFile.write('d:Tf/VirtualSimulationTimeEnd = %1.2f s\n' % (PLAN_DATA["finalTime"]-1))
     parFile.write('i:Tf/VirtualSimulationNumberOfSequentialTimes = %d\n' % int(PLAN_DATA["finalTime"]-1))
     parFile.write('\n') 
@@ -416,7 +432,10 @@ def WriteMainWithVisualizationFile(DATA):
     parFile.write('s:Gr/ViewA/Type        = "OpenGL"\n')
     parFile.write('i:Gr/ViewA/WindowSizeX = 900\n')
     parFile.write('i:Gr/ViewA/WindowSizeY = 900\n')
-    parFile.write('d:Gr/ViewA/Theta          = 89.9 deg\n')
+    if DATA['scoring_quantity'] == 'phasespace':
+        parFile.write('d:Gr/ViewA/Theta          = 0 deg\n')
+    else:
+        parFile.write('d:Gr/ViewA/Theta          = 89.9 deg\n')
     parFile.write('d:Gr/ViewA/Phi            = -90 deg\n')
     parFile.write('b:Gr/ViewA/IncludeAxes    = "true"\n')
     parFile.write('d:Gr/ViewA/AxesSize       = 100 cm \n')
@@ -427,6 +446,29 @@ def WriteMainWithVisualizationFile(DATA):
     parFile.write('#iv:Ge/Patient/ShowSpecificSlicesX = 1 -1\n')
     parFile.write('#iv:Ge/Patient/ShowSpecificSlicesY = 1 -1\n')
     parFile.write('iv:Ge/Patient/ShowSpecificSlicesZ = 1 -1\n')
+    if DATA['scoring_quantity'] == 'phasespace':
+        parFile.write('# Show structures you want to score\n')
+        if DATA["OS"] == "linux/macos":
+            parFile.write('s:Ge/Patient/DicomRTStructFile = "%s"\n') %DATA["RS_filename"]
+        elif DATA["OS"] == "windows":
+            windows_RS_file = DATA["RS_filename"]
+            windows_RS_file = windows_RS_file.replace('/', '\\')
+            RS_file = wslPath.to_posix(windows_RS_file) 
+            parFile.write('s:Ge/Patient/DicomRTStructFile = "%s"\n'%RS_file) 
+        else:
+            parFile.write('s:Ge/Patient/DicomRTStructFile = "%s"\n') %DATA["RS_filename"]
+        roi_list = str(DATA["ROIs_to_score"])
+        roi_list = roi_list.replace('[','')
+        roi_list = roi_list.replace(']','')
+        parFile.write('sv:Ge/Patient/ColorByRTStructNames = %d %s \n'%(len(DATA["ROIs_to_score"]), roi_list))
+        color_list = ''
+        for color in range(0,len(DATA["ROIs_to_score"]),1):
+            if color != len(DATA["ROIs_to_score"])-1:
+                color_list += '"yellow", '
+            else:
+                color_list += '"yellow"'
+        parFile.write('sv:Ge/Patient/ColorByRTStructColors = %d %s\n'%(len(DATA["ROIs_to_score"]),color_list))
+        parFile.write('\n')
     parFile.write('\n')
     parFile.write('# comment the follow line in Geometry.txt to avoid displaying the RT Dose Grid\n')
     parFile.write('#s:Ge/Patient/CloneRTDoseGridFrom = RS/DicomDoseFileName \n')
@@ -539,21 +581,44 @@ def WriteMainFile(DATA,PLAN_DATA):
     parFile.write('\n')
     parFile.write('s:So/phsp/Type = "PhaseSpace"\n')
     parFile.write('s:So/phsp/Component = "source"\n')
-    parFile.write('s:So/phsp/PhaseSpaceFileName = "%s"\n' % DATA["phsp_filename"])
+    if DATA["OS"] == "linux/macos":
+        parFile.write('s:So/phsp/PhaseSpaceFileName = "%s"\n' % DATA["phsp_filename"])
+    elif DATA["OS"] == "windows":
+        windows_PHSP_file = DATA["phsp_filename"]
+        windows_PHSP_file = windows_PHSP_file.replace('/', '\\')
+        PHSP_file = wslPath.to_posix(windows_PHSP_file) 
+        parFile.write('s:So/phsp/PhaseSpaceFileName = "%s"\n'%PHSP_file) 
+    else:
+        parFile.write('s:So/phsp/PhaseSpaceFileName = "%s"\n' % DATA["phsp_filename"])
     parFile.write('b:So/phsp/PhaseSpacePrecheck = "F"\n')
     parFile.write('i:So/phsp/PhaseSpaceMultipleUse = 0\n')
     parFile.write('i:So/phsp/NumberOfHistoriesInRun = Tf/PrimariesPerSegment/Value #10 * Ts/ShowHistoryCountAtInterval\n') 
     parFile.write('b:So/phsp/LimitedAssumeFirstParticleIsNewHistory = "True" #if IAEA phsp\n') 
     parFile.write('b:So/phsp/LimitedAssumePhotonIsNewHistory = "True" #if IAEA phsp\n') 
     parFile.write('\n')
-    parFile.write('s:Sc/Dose/Quantity = "%s" \n' % (DATA["scoring_quantity"]))
-    parFile.write('u:Sc/Dose/OutputWeightingFactor = So/ScalingFactor\n')
-    parFile.write('b:Sc/Dose/PreCalculateStoppingPowerRatios = "True"\n')
-    parFile.write('s:Sc/Dose/Component = "Patient/RTDoseGrid" \n')
-    parFile.write('s:Sc/Dose/IfOutputfileAlreadyExists = "Overwrite"\n')
-    parFile.write('s:Sc/Dose/OutputType = "%s" \n' % DATA["output_format"])
-    parFile.write('s:Sc/Dose/OutputFile = "./output/%s" \n' % (DATA["output_file"]))
-    parFile.write('\n')
+    if DATA["scoring_quantity"] == "PhaseSpace":
+        roinum = 0
+        parFile.write('# To score the phase space in a surface on a DICOM you need OpenTOPAS v4.2.2 or later\n')
+        for roi in DATA["ROIs_to_score"]:
+            roinum += 1
+            parFile.write('s:Sc/PhaseSpace%d/Quantity = "%s" \n' %(roinum,DATA["scoring_quantity"]))
+            parFile.write('b:Sc/PhaseSpace%d/OutputToConsole = "False" \n'%roinum)
+            parFile.write('s:Sc/PhaseSpace%d/Surface = "Patient/AnySurface" \n'%roinum) 
+            parFile.write('sv:Sc/PhaseSpace%d/OnlyIncludeIfInRTStructure = 1 "%s" \n' %(roinum,roi))
+            parFile.write('b:Sc/PhaseSpace%d/KillAfterPhaseSpace = "False"\n'%roinum)
+            parFile.write('s:Sc/PhaseSpace%d/IfOutputfileAlreadyExists = "Overwrite"\n'%roinum)
+            parFile.write('s:Sc/PhaseSpace%d/OutputType = "%s" \n' %(roinum, DATA["output_format"]))
+            parFile.write('s:Sc/PhaseSpace%d/OutputFile = "./output/%s_%s" \n' %(roinum, DATA["output_file"], roi))
+            parFile.write('\n')
+    else:
+        parFile.write('s:Sc/Dose/Quantity = "%s" \n' % (DATA["scoring_quantity"]))
+        parFile.write('u:Sc/Dose/OutputWeightingFactor = So/ScalingFactor\n')
+        parFile.write('b:Sc/Dose/PreCalculateStoppingPowerRatios = "True"\n')
+        parFile.write('s:Sc/Dose/Component = "Patient/RTDoseGrid" \n')
+        parFile.write('s:Sc/Dose/IfOutputfileAlreadyExists = "Overwrite"\n')
+        parFile.write('s:Sc/Dose/OutputType = "%s" \n' % DATA["output_format"])
+        parFile.write('s:Sc/Dose/OutputFile = "./output/%s" \n' % (DATA["output_file"]))
+        parFile.write('\n')
     parFile.write('d:Tf/TimeLineEnd = Tf/VirtualSimulationTimeEnd s\n')
     parFile.write('i:Tf/NumberOfSequentialTimes = Tf/VirtualSimulationNumberOfSequentialTimes\n')
     parFile.write('i:Ts/Seed = 1\n')