-
Notifications
You must be signed in to change notification settings - Fork 1
/
Copy pathmylib_grid.py
executable file
·123 lines (104 loc) · 3.88 KB
/
mylib_grid.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
#!/usr/bin/env python
#
# Author: Patrick Brockmann
# Contact: Patrick.Brockmann@ipsl.jussieu.fr
# $Date: $
# $Name: $
# $Revision: $
# History:
# Modification:
#
import numpy
import vtk
# ********************************
def longitude_create(
arg_lon,
arg_projection,
arg_ratioxy=1):
if arg_projection == 'linear':
n_polypoints = 2
polygonpoints = vtk.vtkPoints()
polygonpoints.SetNumberOfPoints(n_polypoints)
polygonpoints.InsertPoint(0, arg_lon, -90*arg_ratioxy, 0)
polygonpoints.InsertPoint(1, arg_lon, 90*arg_ratioxy, 0)
else:
n_polypoints = 360
polygonpoints = vtk.vtkPoints()
polygonpoints.SetNumberOfPoints(n_polypoints)
deg2rad = numpy.pi/180.
for n in range(n_polypoints):
theta = 2*numpy.pi*n/(n_polypoints-1)
x = numpy.sin(arg_lon*deg2rad)*numpy.sin(theta)
y = numpy.cos(arg_lon*deg2rad)*numpy.sin(theta)
z = numpy.cos(theta)
polygonpoints.InsertPoint(n, x, y, z)
polygonlines = vtk.vtkCellArray()
polygonlines.InsertNextCell(n_polypoints)
for n in range(n_polypoints):
polygonlines.InsertCellPoint(n)
lonpolydata = vtk.vtkPolyData()
lonpolydata.SetPoints(polygonpoints)
lonpolydata.SetLines(polygonlines)
return lonpolydata
# ********************************
def latitude_create(
arg_lat,
arg_projection):
if arg_projection == 'linear':
n_polypoints = 2
polygonpoints = vtk.vtkPoints()
polygonpoints.SetNumberOfPoints(n_polypoints)
polygonpoints.InsertPoint(0, -180, arg_lat, 0)
polygonpoints.InsertPoint(1, 540, arg_lat, 0)
else:
n_polypoints = 360
polygonpoints = vtk.vtkPoints()
polygonpoints.SetNumberOfPoints(n_polypoints)
deg2rad = numpy.pi/180.
for n in range(n_polypoints):
theta = 2*numpy.pi*n/(n_polypoints-1)
x = numpy.cos(arg_lat*deg2rad)*numpy.cos(theta)
y = numpy.cos(arg_lat*deg2rad)*numpy.sin(theta)
z = numpy.sin(arg_lat*deg2rad)
polygonpoints.InsertPoint(n, x, y, z)
polygonlines = vtk.vtkCellArray()
polygonlines.InsertNextCell(n_polypoints)
for n in range(n_polypoints):
polygonlines.InsertCellPoint(n)
latpolydata = vtk.vtkPolyData()
latpolydata.SetPoints(polygonpoints)
latpolydata.SetLines(polygonlines)
return latpolydata
# ********************************
def grid_create(
arg_grid_delta,
arg_grid_color,
arg_grid_width,
arg_projection,
arg_ratioxy):
gridpolydata = vtk.vtkAppendPolyData()
for i in range(0, int(90*arg_ratioxy), int(arg_grid_delta*arg_ratioxy)):
latpolydata = latitude_create(i, arg_projection)
gridpolydata.AddInputData(latpolydata)
for i in range(0, int(-90*arg_ratioxy), int(-arg_grid_delta*arg_ratioxy)):
latpolydata = latitude_create(i, arg_projection)
gridpolydata.AddInputData(latpolydata)
if arg_projection == 'linear':
for i in range(-180, 540, arg_grid_delta):
lonpolydata = longitude_create(i, arg_projection, arg_ratioxy)
gridpolydata.AddInputData(lonpolydata)
else:
for i in range(0, 180, arg_grid_delta):
lonpolydata = longitude_create(i, arg_projection)
gridpolydata.AddInputData(lonpolydata)
gridpolygonmapper = vtk.vtkDataSetMapper()
gridpolygonmapper.SetInputConnection(gridpolydata.GetOutputPort())
gridpolygonactor = vtk.vtkActor()
gridpolygonactor.SetMapper(gridpolygonmapper)
gridpolygonactor.GetProperty().SetColor(arg_grid_color)
gridpolygonactor.GetProperty().SetLineWidth(arg_grid_width)
gridpolygonactor.GetProperty().SetAmbient(1)
gridpolygonactor.GetProperty().SetDiffuse(0)
gridpolygonactor.PickableOff()
return(gridpolygonactor)
# ********************************