From a67a8401933ce41a867fed5558beebf03d6112b5 Mon Sep 17 00:00:00 2001
From: Eduardo Eloy <m47215@alunos.uevora.pt>
Date: Wed, 21 Sep 2022 10:32:04 +0100
Subject: [PATCH] python script to plot 3d pareto frontiers

---
 MODFIRE-Prototype/plot3d.py | 110 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 1 file changed, 110 insertions(+), 0 deletions(-)
 create mode 100644 MODFIRE-Prototype/plot3d.py

diff --git a/MODFIRE-Prototype/plot3d.py b/MODFIRE-Prototype/plot3d.py
new file mode 100644
index 0000000..36a5393
--- /dev/null
+++ b/MODFIRE-Prototype/plot3d.py
@@ -0,0 +1,110 @@
+import numpy as np
+import matplotlib.pyplot as plt
+
+
+#val = input("1:WoodYield/SoilLoss \n 2:WoodYield/RiskPercentile \n 3:SoilLoss/RiskPercentile \n 4:All three
+
+f=open("paretoWSR.csv","r")
+#lines=sorted(list(map(str, f.readlines())))
+
+f2 =open("nonParetoWSR.csv","r")
+#lines2 = sorted(list(map(str, f2.readlines())))
+
+lines=f.readlines()
+lines2 = f2.readlines()
+
+
+paretoPoints = []
+dominatedPoints = []
+x = []
+y = []
+z = []
+lst = []
+scatterx = []
+scattery = []
+
+for i in lines:
+	result = i.split(',')
+	x.append(int(result[0]))
+	y.append(int(result[1]))
+	paretoPoints.append([int(result[0]),int(result[1]),int(result[2])])
+f.close()
+
+for i in lines2:
+	result = i.split(',')
+	x.append(int(result[0]))
+	y.append(int(result[1]))
+	z.append(int(result[2]))
+	if([int(result[0]),int(result[1]),int(result[2])] not in paretoPoints):	
+		dominatedPoints.append([int(result[0]),int(result[1]),int(result[2])])
+f2.close()
+
+#print(dominatedPoints)
+
+def simple_cull(inputPoints, dominates):
+    paretoPoints = set()
+    candidateRowNr = 0
+    dominatedPoints = set()
+    while True:
+        candidateRow = inputPoints[candidateRowNr]
+        inputPoints.remove(candidateRow)
+        rowNr = 0
+        nonDominated = True
+        while len(inputPoints) != 0 and rowNr < len(inputPoints):
+            row = inputPoints[rowNr]
+            if dominates(candidateRow, row):
+                # If it is worse on all features remove the row from the array
+                inputPoints.remove(row)
+                dominatedPoints.add(tuple(row))
+            elif dominates(row, candidateRow):
+                nonDominated = False
+                dominatedPoints.add(tuple(candidateRow))
+                rowNr += 1
+            else:
+                rowNr += 1
+
+        if nonDominated:
+            # add the non-dominated point to the Pareto frontier
+            paretoPoints.add(tuple(candidateRow))
+
+        if len(inputPoints) == 0:
+            break
+    return paretoPoints, dominatedPoints
+
+def dominates(row, candidateRow):
+    return sum([row[x] >= candidateRow[x] for x in range(len(row))]) == len(row)  
+
+import random
+inputPoints = [[random.randint(70,100) for i in range(3)] for j in range(500)]
+#paretoPoints, dominatedPoints = simple_cull(inputPoints, dominates)
+
+#print()
+
+#print(paretoPoints)
+
+import matplotlib.pyplot as plt
+from mpl_toolkits.mplot3d import Axes3D
+
+fig = plt.figure()
+ax = fig.add_subplot(111, projection='3d')
+dp = np.array(dominatedPoints)
+pp = np.array(paretoPoints)
+#print(dp)
+ax.scatter(dp[:,0],dp[:,1],dp[:,2])
+ax.scatter(pp[:,0],pp[:,1],pp[:,2],color='red')
+
+import matplotlib.tri as mtri
+triang = mtri.Triangulation(pp[:,0],pp[:,1])
+ax.plot_trisurf(triang,pp[:,2],color='red')
+
+ax.set_xlabel("Wood Yield",linespacing=5)
+ax.set_ylabel("Soil Loss",linespacing=5)
+ax.set_zlabel("Fire Risk Protection",linespacing=5)
+ax.xaxis.labelpad=7
+ax.yaxis.labelpad=7
+ax.zaxis.labelpad=7
+#plt.xlabel("Wood Yield")
+#plt.ylabel("Soil Loss")
+#plt.clabel("Fire Risk Protection")
+
+plt.show()
\ No newline at end of file
--
libgit2 0.21.2