Coverage for steam_pysigma\plotters\PlotterRoxie.py: 100%

1# STEAM PySigma is a python wrapper of STEAM-SIGMA written in Java. 

2# Copyright (C) 2023, CERN, Switzerland. All rights reserved. 

3# 

4# This program is free software: you can redistribute it and/or modify 

5# it under the terms of the GNU General Public License as published by 

6# the Free Software Foundation, version 3 of the License. 

7# 

8# This program is distributed in the hope that it will be useful, 

9# but WITHOUT ANY WARRANTY; without even the implied warranty of 

10# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 

11# GNU General Public License for more details. 

12# 

13# You should have received a copy of the GNU General Public License 

14# along with this program. If not, see <https://www.gnu.org/licenses/>. 

15 

16import copy 

17import math 

18import numpy as np 

19import matplotlib.pyplot as plt 

20import matplotlib.lines as lines 

21import matplotlib.patches as patches 

22 

23from steam_pysigma.data import DataRoxieParser as pd 

24from steam_pysigma.utils.Utils import displayWaitAndClose 

25 

26 

27def arcCenterFromThreePoints(a, b, c): 

28 ab = [a.x - b.x, a.y - b.y] 

29 ac = [a.x - c.x, a.y - c.y] 

30 sac = [a.x * a.x - c.x * c.x, a.y * a.y - c.y * c.y] 

31 sba = [b.x * b.x - a.x * a.x, b.y * b.y - a.y * a.y] 

32 yy = (sac[0] * ab[0] + sac[1] * ab[0] + sba[0] * ac[0] + sba[1] * ac[0]) / \ 

33 (2 * ((c.y - a.y) * ab[0] - (b.y - a.y) * ac[0])) 

34 xx = (sac[0] * ab[1] + sac[1] * ab[1] + sba[0] * ac[1] + sba[1] * ac[1]) / \ 

35 (2 * ((c.x - a.x) * ab[1] - (b.x - a.x) * ac[1])) 

36 return [-xx, -yy] 

37 

38 

39def plotIronGeometry(iron, selectedFont): 

40 plt.figure(figsize=(7.5, 7.5)) 

41 ax = plt.axes() 

42 

43 max_x = 0 

44 max_y = 0 

45 

46 for point_name, point in iron.key_points.items(): 

47 max_x = max(point.x, max_x) 

48 max_y = max(point.y, max_y) 

49 

50 for line_name, line in iron.hyper_lines.items(): 

51 if line.type == 'line': 

52 ax.add_line(lines.Line2D([iron.key_points[line.kp1].x, iron.key_points[line.kp2].x], 

53 [iron.key_points[line.kp1].y, iron.key_points[line.kp2].y], 

54 color='black', linewidth=1)) 

55 

56 elif line.type == 'arc': 

57 pt1 = iron.key_points[line.kp1] 

58 pt2 = iron.key_points[line.kp2] 

59 center = arcCenterFromThreePoints(pt1, iron.key_points[line.kp3], pt2) 

60 radius = (np.sqrt(np.square(pt1.x - center[0]) + np.square(pt1.y - center[1])) + 

61 np.sqrt(np.square(pt2.x - center[0]) + np.square(pt2.y - center[1]))) / 2 

62 if pt1.x < pt2.x and pt1.x < center[0] and pt1.y < pt2.y and pt1.y < center[1]: 

63 th1 = - np.arctan2(pt1.y - center[1], pt1.x - center[0]) * 180 / np.pi 

64 else: 

65 th1 = np.arctan2(pt1.y - center[1], pt1.x - center[0]) * 180 / np.pi 

66 th2 = np.arctan2(pt2.y - center[1], pt2.x - center[0]) * 180 / np.pi 

67 ax.add_patch(patches.Arc((center[0], center[1]), width=2 * radius, height=2 * radius, angle=0, 

68 theta1=min(th1, th2), theta2=max(th1, th2), color='blue', linewidth=1)) 

69 

70 # elif line.type == 'ellipticArc': 

71 # pt1 = iron.key_points[line.kp1] 

72 # pt2 = iron.key_points[line.kp2] 

73 # r1 = (pt1.x - pt2.x) / (2 * line.arg1) 

74 # r2 = (pt1.y - pt2.y) / (2 * line.arg2) 

75 # a1 = np.arctan2(- r1 / r2) # (t1 + t2) / 2 

76 # a2 = np.arcsin(np.sqrt(r1**2 + r2**2)) # (t1 - t2) / 2 

77 # center = [pt1.x - line.arg1 * np.cos(a1 + a2), pt1.y - line.arg2 * np.sin(a1 + a2)] 

78 # 

79 # if pt1.x < pt2.x and pt1.x < center[0] and pt1.y < pt2.y and pt1.y < center[1]: 

80 # th1 = - np.arctan2(pt1.y - center[1], pt1.x - center[0]) * 180 / np.pi 

81 # else: 

82 # th1 = np.arctan2(pt1.y - center[1], pt1.x - center[0]) * 180 / np.pi 

83 # th2 = np.arctan2(pt2.y - center[1], pt2.x - center[0]) * 180 / np.pi 

84 # ax.add_patch(patches.Arc((center[0], center[1]), width=2 * line.arg1, height=2 * line.arg2, angle=0, 

85 # theta1=min(th1, th2), theta2=max(th1, th2), color='purple', linewidth=1)) 

86 

87 elif line.type == 'circle': 

88 pt1 = iron.key_points[line.kp1] 

89 pt2 = iron.key_points[line.kp2] 

90 center = [(pt1.x + pt2.x) / 2, (pt1.y + pt2.y) / 2] 

91 radius = (np.sqrt(np.square(pt1.x - center[0]) + np.square(pt1.y - center[1])) + 

92 np.sqrt(np.square(pt2.x - center[0]) + np.square(pt2.y - center[1]))) / 2 

93 ax.add_patch(patches.Circle((center[0], center[1]), 

94 radius=radius, fill=False, edgecolor='green', linewidth=1)) 

95 

96 ax.set_xlim(0, 1.1 * max_x) 

97 ax.set_ylim(0, max_y + 0.1 * max_x) 

98 plt.xlabel('x [m]', **selectedFont) 

99 plt.ylabel('y [m]', **selectedFont) 

100 plt.title('Iron Yoke', **selectedFont) 

101 plt.set_cmap('jet') 

102 plt.rcParams.update({'font.size': 12}) 

103 

104 

105# def plotCoilGeometry(roxie_data, ax): 

106# max_x = 0 

107# max_y = 0 

108# for coil_nr, coil in roxie_data.coil.coils.items(): 

109# for pole_nr, pole in coil.poles.items(): 

110# for layer_nr, layer in pole.layers.items(): 

111# for winding_key, winding in layer.windings.items(): 

112# for block_key, block in winding.blocks.items(): 

113# for halfTurn_nr, halfTurn in block.half_turns.items(): 

114# iL = halfTurn.corners.insulated.iL 

115# iR = halfTurn.corners.insulated.iR 

116# oL = halfTurn.corners.insulated.oL 

117# oR = halfTurn.corners.insulated.oR 

118# max_x = max(oL.x, oR.x, max_x) 

119# max_y = max(oL.y, oR.y, max_y) 

120# 

121# ax.add_line(lines.Line2D([iL.x, iR.x], [iL.y, iR.y], color='red')) 

122# ax.add_line(lines.Line2D([oL.x, oR.x], [oL.y, oR.y], color='red')) 

123# ax.add_line(lines.Line2D([oR.x, iR.x], [oR.y, iR.y], color='red')) 

124# ax.add_line(lines.Line2D([iL.x, oL.x], [iL.y, oL.y], color='red')) 

125# # cc = roxie_data.coil.coils[1].bore_center 

126# # ax.set_xlim((len(roxie_data.coil.coils) == 1) * 2 * cc.x - (1.1 * (max_x - cc.x) + cc.x), 

127# # 1.1 * (max_x - cc.x) + cc.x) 

128# # ax.set_ylim((len(roxie_data.coil.coils) == 1) * 2 * cc.y - (max_y + 0.1 * (max_x - cc.x)), 

129# # max_y + 0.1 * (max_x - cc.x)) 

130 

131 

132# Plot conductors and their numbers 

133def plotEdges(xPos, yPos, xBarePos, yBarePos, iPos, selectedFont): 

134 plt.figure(figsize=(10, 10)) 

135 # Plot edges 

136 for c, (cXPos, cYPos) in enumerate(zip(xPos, yPos)): 

137 pt1, pt2, pt3, pt4 = (cXPos[0], cYPos[0]), (cXPos[1], cYPos[1]), (cXPos[2], cYPos[2]), (cXPos[3], cYPos[3]) 

138 if iPos[c] > 0: 

139 line = plt.Polygon([pt1, pt2, pt3, pt4], closed=True, fill=True, facecolor='r', edgecolor='k', alpha=.25) 

140 else: 

141 line = plt.Polygon([pt1, pt2, pt3, pt4], closed=True, fill=True, facecolor='b', edgecolor='k', alpha=.25) 

142 plt.gca().add_artist(line) 

143 

144 # Plot conductor numbers 

145 x_ave_cond, y_ave_cond = sum(cXPos) / len(cXPos), sum(cYPos) / len(cYPos) 

146 plt.text(x_ave_cond, y_ave_cond, '{}'.format(c + 1)) 

147 

148 # Plot edges of bare conductors 

149 for c, (cXBarePos, cYBarePos) in enumerate(zip(xBarePos, yBarePos)): 

150 pt1, pt2, pt3, pt4 = (cXBarePos[0], cYBarePos[0]), (cXBarePos[1], cYBarePos[1]), \ 

151 (cXBarePos[2], cYBarePos[2]), (cXBarePos[3], cYBarePos[3]) 

152 if iPos[c] > 0: 

153 line = plt.Polygon([pt1, pt2, pt3, pt4], closed=True, fill=False, facecolor='r', edgecolor='k', alpha=.25) 

154 else: 

155 line = plt.Polygon([pt1, pt2, pt3, pt4], closed=True, fill=False, facecolor='b', edgecolor='k', alpha=.25) 

156 plt.gca().add_artist(line) 

157 

158 plt.xlabel('x [m]', **selectedFont) 

159 plt.ylabel('y [m]', **selectedFont) 

160 plt.title('Conductors and their numbers', **selectedFont) 

161 plt.set_cmap('jet') 

162 plt.rcParams.update({'font.size': 12}) 

163 plt.axis('equal') 

164 plt.grid() 

165 

166 

167def plot_all(roxie_data: pd.RoxieData): 

168 """ 

169 Plot all default plots 

170 """ 

171 selectedFont = {'fontname': 'DejaVu Sans', 'size': 14} 

172 

173 if roxie_data.iron: 

174 plotIronGeometry(roxie_data.iron, selectedFont) 

175 # plotCoilGeometry(roxie_data, ax) 

176 

177 xPos = [] 

178 yPos = [] 

179 xBarePos = [] 

180 yBarePos = [] 

181 iPos = [] 

182 for eo in roxie_data.coil.physical_order: 

183 winding = roxie_data.coil.coils[eo.coil].poles[eo.pole].layers[eo.layer].windings[eo.winding] 

184 block = winding.blocks[eo.block] 

185 for halfTurn_nr, halfTurn in block.half_turns.items(): 

186 insu = halfTurn.corners.insulated 

187 bare = halfTurn.corners.bare 

188 

189 xPos.append([insu.iH.x, insu.oH.x, insu.oL.x, insu.iL.x]) 

190 yPos.append([insu.iH.y, insu.oH.y, insu.oL.y, insu.iL.y]) 

191 xBarePos.append([bare.iH.x, bare.oH.x, bare.oL.x, bare.iL.x]) 

192 yBarePos.append([bare.iH.y, bare.oH.y, bare.oL.y, bare.iL.y]) 

193 iPos.append(block.current_sign) 

194 plotEdges(xPos, yPos, xBarePos, yBarePos, iPos, selectedFont) 

195 

196 displayWaitAndClose(waitTimeBeforeMessage=.1, waitTimeAfterMessage=10)