source: lambert/lamberttest.py @ f4ab898fcdc0208d315c5ec7ca8ae38ab1160e27

Revision f4ab898fcdc0208d315c5ec7ca8ae38ab1160e27, 13.7 KB checked in by Administrateur <Administrateur@…>, 6 years ago (diff)

Moved remotely

git-svn-id:  http://kawa.selfip.org/svn/projects@5 532e76a9-45ae-c241-9c6d-8309ac6440cd

  • Property mode set to 100644
Line 
1import math
2from graphics import PyGameGLWindow
3import pygame
4from OpenGL.GL import *
5from OpenGL.GLU import *
6import lambert
7import pprint
8import sys
9import Image
10import PngImagePlugin
11
12
13LATITUDE = 0
14LONGITUDE = 1
15MESH = 2
16POS = 0
17NORMAL = 1
18COORDS = 2
19PHI = 0
20LAMBDA = 1
21
22def length(vector):
23  squarredLength = 0.0
24  for dimension in vector:
25    squarredLength += dimension * dimension
26  return math.sqrt(squarredLength)
27
28def normalize(vector):
29  len = length(vector)
30  result = []
31  for dimension in vector:
32    result.append(dimension / len)
33  return result
34
35def vecsub(v1, v2):
36  result = []
37  for d1, d2 in zip(v1, v2):
38    result.append(d1 - d2)
39  return result
40
41def vecadd(v1, v2):
42  result = []
43  for d1, d2 in zip(v1, v2):
44    result.append(d1 + d2)
45  return result
46
47def extends(vectors):
48  minimum = array(vectors[0])
49  maximum = array(vectors[0])
50  i = 1
51  while i < len(vectors):
52    for j, (minComp, maxComp, comp) in enumerate(zip(minimum, maximum, vectors[i])):
53      if comp < minComp:
54        minimum[j] = comp
55      elif comp > maxComp:
56        maximum[j] = comp
57    i += 1
58  return minimum, maximum
59
60
61class Texture:
62  def __init__(self, file):
63    self.load(file)
64
65  def load(self, file):
66    self.texture = Image.open(file)
67    ix, iy = self.texture.size
68    image = self.texture.tostring("raw", "RGBX", 0, -1)
69    self.id = glGenTextures(1)
70    self.bind()
71    glPixelStorei(GL_UNPACK_ALIGNMENT, 1)
72    glTexImage2D(GL_TEXTURE_2D, 0, GL_RGBA, ix, iy, 0, GL_RGBA, GL_UNSIGNED_BYTE, image)
73    glTexParameteri(GL_TEXTURE_2D, GL_TEXTURE_MIN_FILTER, GL_LINEAR)
74    glTexParameteri(GL_TEXTURE_2D, GL_TEXTURE_MAG_FILTER, GL_LINEAR)
75
76  def bind(self):
77    glBindTexture(GL_TEXTURE_2D, self.id)
78   
79
80class Window(PyGameGLWindow):
81  MB1 = 1
82  MB2 = 2
83  MB3 = 4
84 
85  def __init__(self, app, width, height, name):
86    self.app = app
87    PyGameGLWindow.__init__(self, width, height, name)
88    self.mouseButtons = 0
89
90  def setup(self, width, height):
91    # General setup
92    glShadeModel(GL_FLAT)
93    glClearColor(0.0, 0.0, 0.0, 0.0)
94    glClearDepth(1.0)
95    glEnable(GL_DEPTH_TEST)
96    glDepthFunc(GL_LEQUAL)
97    glEnable(GL_LIGHTING)
98    glEnable(GL_BLEND)
99    # Enable polygon offset for filled primitives
100    glEnable(GL_POLYGON_OFFSET_FILL)
101    glPolygonOffset(1, 1)
102    # Rendering hints
103    glEnable(GL_LINE_SMOOTH)
104    glBlendFunc(GL_SRC_ALPHA, GL_ONE_MINUS_SRC_ALPHA)
105    glHint(GL_LINE_SMOOTH_HINT, GL_NICEST)
106    glHint(GL_PERSPECTIVE_CORRECTION_HINT, GL_NICEST)
107
108  def onKeyPressed(self, event):
109    if event.key == pygame.K_ESCAPE:
110      sys.exit(0)
111
112  def onMouseDown(self, event):
113    self.mouseButtons |= 1 << (event.button - 1)
114
115  def onMouseUp(self, event):
116    self.mouseButtons &= ~(1 << (event.button - 1))
117
118  def onMouseMotion(self, event):
119    if self.mouseButtons & self.MB1:
120      self.app.mouseRot = (self.app.mouseRot[0] + self.app.mouseRotIncrement[0] * event.rel[0], 
121                           self.app.mouseRot[1] + self.app.mouseRotIncrement[1] * event.rel[1])
122    if self.mouseButtons & self.MB3:
123      self.app.mouseZoom += self.app.mouseZoomIncrement * event.rel[1]
124    if self.mouseButtons & self.MB2:
125      self.app.mousePan = (self.app.mousePan[0] + self.app.mousePanIncrement[0] * event.rel[0], 
126                           self.app.mousePan[1] + self.app.mousePanIncrement[1] * event.rel[1])
127
128
129class Sphere:
130  POINTS = GL_POINTS
131  LINE_LOOP = GL_LINE_LOOP
132  LINES = GL_LINES
133  LINE_STRIP = GL_LINE_STRIP
134
135  def points(self, radius, latitudeRes, longitudeRes):
136    # Generate latitudes and longitudes
137    latitudes = []
138    longitudes = [[None for y in range(longitudeRes + 2)] for x in range(latitudeRes)]
139    for j in range(longitudeRes):
140      latitude = [None for x in range(latitudeRes)]
141      floatJ = (j + 1) / float(longitudeRes + 1)
142      phi = (floatJ * math.pi) - (math.pi / 2.0)     
143      #y = (floatJ * radius * 2.0) - radius
144      y = radius * math.sin(phi)
145      latitudeRadius = math.cos(math.asin(y))
146      for i in range(latitudeRes / 2):
147        floatI = i / float(latitudeRes)
148        lambd1 = (floatI * math.pi * 2.0) - math.pi
149        lambd2 = math.fmod(lambd1 + math.pi, math.pi)
150        x = latitudeRadius * math.sin(lambd1)
151        z = latitudeRadius * math.cos(lambd1)
152        point1 = ((x, y, z), normalize((x, y, z)), (phi, lambd1))
153        point2 = ((-x, y, -z), normalize((-x, y, -z)), (phi, lambd2))
154        top1 = ((0.0, radius, 0.0), (0.0, 1.0, 0.0), (math.pi / 2.0, lambd1))
155        bottom1 = ((0.0, -radius, 0.0), (0.0, -1.0, 0.0), (-math.pi / 2.0, lambd1))
156        top2 = ((0.0, radius, 0.0), (0.0, 1.0, 0.0), (math.pi / 2.0, lambd2))
157        bottom2 = ((0.0, -radius, 0.0), (0.0, -1.0, 0.0), (-math.pi / 2.0, lambd2))
158        latitude[i] = point1
159        latitude[i + (latitudeRes / 2)] = point2       
160        longitudes[i][0] = bottom1
161        longitudes[i][-1] = top1
162        longitudes[i][j + 1] = point1
163        longitudes[i + (latitudeRes / 2)][0] = bottom2
164        longitudes[i + (latitudeRes / 2)][-1] = top2
165        longitudes[i + (latitudeRes / 2)][j + 1] = point2
166      latitudes.append(latitude)
167     
168    # Generate mesh
169    mesh = {}
170    # Generate strips
171    mesh["strips"] = []
172    prevLatitude = None
173    for latitude in latitudes:
174      if prevLatitude is not None:
175        strip = []
176        for point1, point2 in zip(prevLatitude, latitude):
177          strip.append(point1)
178          strip.append(point2)
179        strip.append((prevLatitude[0][POS], prevLatitude[0][NORMAL], (prevLatitude[0][COORDS][PHI], math.pi)))
180        strip.append((latitude[0][POS], latitude[0][NORMAL], (latitude[0][COORDS][PHI], math.pi)))
181        mesh["strips"].append(strip)
182      prevLatitude = latitude
183    # Generate fans
184    mesh["fans"] = []
185    fan = [((0.0, -radius, 0.0), (0.0, -1.0, 0.0), (-math.pi / 2.0, 0.0))]
186    for point in latitudes[0]:
187      fan.append(point)
188    fan.append((latitudes[0][0][POS], latitudes[0][0][NORMAL], (latitudes[0][0][COORDS][PHI], math.pi)))
189    mesh["fans"].append(fan)
190    fan = [((0.0, radius, 0.0), (0.0, 1.0, 0.0), (math.pi / 2.0, 0.0))]
191    for point in latitudes[-1]:
192      fan.append(point)
193    fan.append((latitudes[-1][0][POS], latitudes[-1][0][NORMAL], (latitudes[-1][0][COORDS][PHI], math.pi)))
194    mesh["fans"].append(fan)
195     
196    return latitudes, longitudes, mesh
197
198  def renderArray(self, points, primitive):
199    glBegin(primitive)
200    for point in points:
201      glNormal3fv(point[1])
202      glVertex3fv(point[0])
203    glEnd()
204
205  def render2DArray(self, points, primitive):
206    glBegin(primitive)
207    for point in points:
208      glVertex2fv(point)
209    glEnd()
210
211  def render(self, points, primitive=GL_LINE_LOOP):
212    for latitude in points[LATITUDE]:
213      self.renderArray(latitude, self.LINE_LOOP)
214    for longitude in points[LONGITUDE]:
215      self.renderArray(longitude, self.LINE_STRIP)
216
217  def renderMesh(self, points):
218    for strip in points[MESH]["strips"]:
219      glBegin(GL_QUAD_STRIP)
220      for point in strip:
221        glTexCoord2f((point[COORDS][LAMBDA] + math.pi) / (math.pi * 2.0), (point[COORDS][PHI] + math.pi / 2.0) / math.pi)
222        glNormal3fv(point[NORMAL])
223        glVertex3fv(point[POS])
224      glEnd()
225    for fan in points[MESH]["fans"]:
226      glBegin(GL_TRIANGLE_FAN)
227      for point in fan:
228        glTexCoord2f((point[COORDS][LAMBDA] + math.pi) / (math.pi * 2.0), (point[COORDS][PHI] + math.pi / 2.0) / math.pi)
229        glNormal3fv(point[NORMAL])
230        glVertex3fv(point[POS])
231      glEnd()
232
233  def renderClipped(self, points):
234    for latitude in points[LATITUDE]:
235      self.renderArray(latitude, self.LINE_STRIP)
236    for longitude in points[LONGITUDE]:
237      self.renderArray(longitude, self.LINE_STRIP)
238
239  def render2D(self, points, primitive=GL_LINE_STRIP):
240    for latitude in points[LATITUDE]: 
241      self.render2DArray(latitude, primitive)
242    for longitude in points[LONGITUDE]:
243      self.render2DArray(longitude, primitive)
244
245  def clip(self, points, minCoords, maxCoords):
246    # Get points lying between minPos and maxPos
247    clippedLatitudes = []
248    for latitude in points[LATITUDE]:
249      if latitude[0][COORDS][PHI] >= minCoords[PHI] and latitude[0][COORDS][PHI] <= maxCoords[PHI]:
250        clippedLatitude = []
251        for point in latitude:
252          if point[COORDS][LAMBDA] >= minCoords[LAMBDA] and point[COORDS][LAMBDA] <= maxCoords[LAMBDA]:
253            clippedLatitude.append(point)
254        clippedLatitudes.append(clippedLatitude)
255     
256    clippedLongitudes = []
257    for longitude in points[LONGITUDE]:
258      if longitude[0][COORDS][LAMBDA] >= minCoords[LAMBDA] and longitude[0][COORDS][LAMBDA] <= maxCoords[LAMBDA]:
259        clippedLongitude = []
260        for point in longitude:
261          if point[COORDS][PHI] >= minCoords[PHI] and point[COORDS][PHI] <= maxCoords[PHI]:
262            clippedLongitude.append(point)
263        clippedLongitudes.append(clippedLongitude)
264
265    return clippedLatitudes, clippedLongitudes, points[MESH]
266       
267  def project(self, points, projection, normalize=True, shift=True):
268    projectedLatitudes = []
269    projectedLongitudes = []
270    min = [float(2**30), float(2**30)]
271    max = [-float(2**30), -float(2**30)]
272    for latitude in points[LATITUDE]:
273      projectedLatitude = []
274      for point in latitude:
275        ppoint = projection.project(point[COORDS])
276        for i in range(2):
277          if ppoint[i] < min[i]:
278            min[i] = ppoint[i]
279          if ppoint[i] > max[i]:
280            max[i] = ppoint[i]
281        projectedLatitude.append(ppoint)
282      projectedLatitudes.append(projectedLatitude)
283    for longitude in points[LONGITUDE]:
284      projectedLongitude = []
285      for point in longitude:
286        projectedLongitude.append(projection.project(point[COORDS]))
287      projectedLongitudes.append(projectedLongitude)
288
289    scale = length(vecsub(max, min))
290    if scale == 0.0:
291      scale = 1.0
292    print max, min
293    if normalize:
294      for j in range(len(projectedLatitudes)):
295        for i in range(len(projectedLatitudes[j])):
296          if shift:
297            projectedLatitudes[j][i] = (projectedLatitudes[j][i][0] - min[0], projectedLatitudes[j][i][1] - min[1])
298          projectedLatitudes[j][i] = (projectedLatitudes[j][i][0] / scale, projectedLatitudes[j][i][1] / scale)
299      for j in range(len(projectedLongitudes)):
300        for i in range(len(projectedLongitudes[j])):
301          if shift:
302            projectedLongitudes[j][i] = (projectedLongitudes[j][i][0] - min[0], projectedLongitudes[j][i][1] - min[1])
303          projectedLongitudes[j][i] = (projectedLongitudes[j][i][0] / scale, projectedLongitudes[j][i][1] / scale)
304
305    pp = pprint.PrettyPrinter()
306    #pp.pprint(projectedLongitudes)
307    print len(points[LATITUDE]), len(points[LONGITUDE]), len(projectedLatitudes), len(projectedLongitudes)
308    return projectedLatitudes, projectedLongitudes
309       
310
311class App:
312  def __init__(self):
313    self.mouseRot = (0.0, 0.0)
314    self.mouseRotIncrement = (0.1, 0.1)
315    self.mouseZoom = -3.0
316    self.mouseZoomIncrement = 0.01
317    self.mousePan = (0.0, 0.0)
318    self.mousePanIncrement = (0.005, -0.005)
319    self.window = Window(self, 640, 480, "Lambert")
320    self.run()
321
322  def run(self):
323    glLightfv(GL_LIGHT0, GL_AMBIENT, (0.5, 0.5, 0.5, 1.0))              # Setup The Ambient Light
324    glLightfv(GL_LIGHT0, GL_DIFFUSE, (1.0, 1.0, 1.0, 1.0))              # Setup The Diffuse Light
325    glLightfv(GL_LIGHT0, GL_POSITION, (10.0, 20.0, 10.0, 1.0))    # Position The Light
326    glEnable(GL_LIGHT0)                                                               # Enable Light One
327
328    worldTexture = Texture("EarthMap.png")
329   
330    quadratic = gluNewQuadric()
331    gluQuadricNormals(quadratic, GLU_SMOOTH)
332    gluQuadricTexture(quadratic, GL_TRUE)
333
334    franceExtends = {"min": (math.radians(42), math.radians(-5)), "max": (math.radians(52), math.radians(10))}
335    largeExtends = {"min": (math.radians(10), math.radians(-50)), "max": (math.radians(52), math.radians(10))}
336    sphere = Sphere()
337    points = sphere.points(1.0, 360, 180)
338    points = sphere.clip(points, franceExtends["min"], franceExtends["max"])
339    #points = sphere.clip(points, largeExtends["min"], largeExtends["max"])
340    #points = sphere.clip(points, (0.0, 0.0), (math.pi / 4.0, math.pi / 2.0))
341    projectedPoints = sphere.project(points, lambert.STD_PROJECTIONS["Lambert 93"], normalize=True, shift=True)
342
343    cycleLength = (len(points[LATITUDE]), len(points[LONGITUDE]))   
344
345    xrot = yrot = zrot = 0.0
346
347    while True:
348      self.window.processEvents()
349      self.window.clear()     
350     
351      glMatrixMode(GL_MODELVIEW)
352      glLoadIdentity()
353      glTranslate(self.mousePan[0], self.mousePan[1], self.mouseZoom)
354      glRotatef(self.mouseRot[0], 0.0, 1.0, 0.0)
355      glRotatef(self.mouseRot[1], 1.0, 0.0, 0.0)
356
357      # gluSphere(quadratic, 1.0, 32, 32)
358      glEnable(GL_LIGHTING)
359      glEnable(GL_TEXTURE_2D)
360      glPolygonMode(GL_BACK, GL_FILL)
361      sphere.renderMesh(points)
362     
363      glDisable(GL_LIGHTING)
364      glDisable(GL_TEXTURE_2D)
365      glPolygonMode(GL_BACK, GL_LINE)
366      glDepthMask(GL_FALSE)
367      #sphere.renderMesh(points)
368      sphere.renderClipped(points)
369      glDepthMask(GL_TRUE)
370
371      glTranslate(1.5, 0.0, 0.0)
372      sphere.render2D(projectedPoints, Sphere.LINE_STRIP)
373      """sphere.renderClipped(lambert.STD_PROJECTIONS["Lambert II extended"],
374          points,
375          (math.radians(42), math.radians(-5)),
376          (math.radians(52), math.radians(10)))     """
377      """sphere.renderClipped(lambert.STD_PROJECTIONS["Lambert II extended"],
378          points,
379          (0.0, 0.0),
380          (math.pi / 2.0, math.pi / 2.0))"""
381     
382      self.window.swap()
383
384if __name__ == "__main__":
385  app = App()
Note: See TracBrowser for help on using the repository browser.