GRASS GIS 8 Programmer's Manual  8.4.0dev(2024)-81301443e3
gsd_views.c
Go to the documentation of this file.
1 /*!
2  \file lib/ogsf/gsd_views.c
3 
4  \brief OGSF library - manipulating views (lower level functions)
5 
6  GRASS OpenGL gsurf OGSF Library
7 
8  (C) 1999-2008 by the GRASS Development Team
9 
10  This program is free software under the
11  GNU General Public License (>=v2).
12  Read the file COPYING that comes with GRASS
13  for details.
14 
15  \author Bill Brown USACERL (January 1993)
16  \author Doxygenized by Martin Landa <landa.martin gmail.com> (May 2008)
17  */
18 
19 #include <grass/config.h>
20 
21 #if defined(OPENGL_X11) || defined(OPENGL_WINDOWS)
22 #include <GL/gl.h>
23 #include <GL/glu.h>
24 #elif defined(OPENGL_AQUA)
25 #include <OpenGL/gl.h>
26 #include <OpenGL/glu.h>
27 #endif
28 
29 #include <grass/ogsf.h>
30 #include "math.h"
31 
32 /*!
33  \brief ADD
34 
35  \param vect
36  \param sx, sy screen coordinates
37 
38  \return 1
39  */
40 int gsd_get_los(float (*vect)[3], short sx, short sy)
41 {
42  double fx, fy, fz, tx, ty, tz;
43  GLdouble modelMatrix[16], projMatrix[16];
44  GLint viewport[4];
45 
46  GS_ready_draw();
47  glPushMatrix();
48 
49  gsd_do_scale(1);
50  glGetDoublev(GL_MODELVIEW_MATRIX, modelMatrix);
51  glGetDoublev(GL_PROJECTION_MATRIX, projMatrix);
52  glGetIntegerv(GL_VIEWPORT, viewport);
53  glPopMatrix();
54 
55  /* OGLXXX XXX I think this is backwards gluProject(XXX); */
56  /* WAS: mapw(Vobj, sx, sy, &fx, &fy, &fz, &tx, &ty, &tz); */
57  gluUnProject((GLdouble)sx, (GLdouble)sy, 0.0, modelMatrix, projMatrix,
58  viewport, &fx, &fy, &fz);
59  gluUnProject((GLdouble)sx, (GLdouble)sy, 1.0, modelMatrix, projMatrix,
60  viewport, &tx, &ty, &tz);
61  vect[FROM][X] = fx;
62  vect[FROM][Y] = fy;
63  vect[FROM][Z] = fz;
64  vect[TO][X] = tx;
65  vect[TO][Y] = ty;
66  vect[TO][Z] = tz;
67 
68  /* DEBUG - should just be a dot */
69  /* OGLXXX frontbuffer: other possibilities include GSD_BOTH */
71  glPushMatrix();
72  gsd_do_scale(1);
73  gsd_linewidth(3);
74  gsd_color_func(0x8888FF);
75 
76  /* OGLXXX for multiple, independent line segments: use GL_LINES */
77  glBegin(GL_LINE_STRIP);
78  glVertex3fv(vect[FROM]);
79  glVertex3fv(vect[TO]);
80  glEnd();
81 
82  gsd_linewidth(1);
83  glPopMatrix();
84 
85  /* OGLXXX frontbuffer: other possibilities include GSD_BOTH */
87 
88  return (1);
89 }
90 
91 #if 0
92 /*!
93  \brief Set view
94 
95  Establishes viewing & projection matrices
96 
97  \param gv view (geoview)
98  \param dp display (geodisplay)
99  */
100 void gsd_set_view(geoview * gv, geodisplay * gd)
101 {
102  double up[3];
103  GLint mm;
104 
105  /* will expand when need to check for in focus, ortho, etc. */
106 
107  gsd_check_focus(gv);
108  gsd_get_zup(gv, up);
109 
110  gd->aspect = GS_get_aspect();
111 
112  glGetIntegerv(GL_MATRIX_MODE, &mm);
113  glMatrixMode(GL_PROJECTION);
114  glLoadIdentity();
115  gluPerspective((double).1 * (gv->fov), (double)gd->aspect,
116  (double)gd->nearclip, (double)gd->farclip);
117 
118  glMatrixMode(mm);
119 
120  glLoadIdentity();
121 
122  /* update twist parm */
123  glRotatef((float)(gv->twist / 10.), 0.0, 0.0, 1.0);
124 
125  /* OGLXXX lookat: replace UPx with vector */
126  gluLookAt((double)gv->from_to[FROM][X], (double)gv->from_to[FROM][Y],
127  (double)gv->from_to[FROM][Z], (double)gv->from_to[TO][X],
128  (double)gv->from_to[TO][Y], (double)gv->from_to[TO][Z],
129  (double)up[X], (double)up[Y], (double)up[Z]);
130 
131  /* have to redefine clipping planes when view changes */
132 
134 
135  return;
136 }
137 #endif
138 /*!
139  \brief Set view
140 
141  Establishes viewing & projection matrices
142 
143  \param gv view (geoview)
144  \param dp display (geodisplay)
145  */
147 {
148  double up[3];
149  float pos[3];
150  int i;
151  GLdouble modelMatrix[16];
152  GLint mm;
153 
154  /* will expand when need to check for in focus, ortho, etc. */
155 
156  gsd_check_focus(gv);
157  gsd_get_zup(gv, up);
158 
159  gd->aspect = GS_get_aspect();
160 
161  glGetIntegerv(GL_MATRIX_MODE, &mm);
162  glMatrixMode(GL_PROJECTION);
163  glLoadIdentity();
164  gluPerspective((double).1 * (gv->fov), (double)gd->aspect,
165  (double)gd->nearclip, (double)gd->farclip);
166 
167  glMatrixMode(mm);
168 
169  glLoadIdentity();
170 
171  /* update twist parm */
172  glRotatef((float)(gv->twist / 10.), 0.0, 0.0, 1.0);
173 
174  /* OGLXXX lookat: replace UPx with vector */
175  gluLookAt((double)gv->from_to[FROM][X], (double)gv->from_to[FROM][Y],
176  (double)gv->from_to[FROM][Z], (double)gv->from_to[TO][X],
177  (double)gv->from_to[TO][Y], (double)gv->from_to[TO][Z],
178  (double)up[X], (double)up[Y], (double)up[Z]);
179 
180  /* rotate to get rotation matrix and then save it */
181  if (gv->rotate.do_rot) {
182 
183  glPushMatrix();
184  glLoadMatrixd(gv->rotate.rotMatrix);
185 
186  glRotated(gv->rotate.rot_angle, gv->rotate.rot_axes[0],
187  gv->rotate.rot_axes[1], gv->rotate.rot_axes[2]);
188  glGetDoublev(GL_MODELVIEW_MATRIX, modelMatrix);
189 
190  for (i = 0; i < 16; i++) {
191  gv->rotate.rotMatrix[i] = modelMatrix[i];
192  }
193 
194  glPopMatrix();
195  }
196 
197  gs_get_datacenter(pos);
198  gsd_surf2model(pos);
199  /* translate rotation center to view center, rotate and translate back */
200  glTranslatef(pos[0], pos[1], pos[2]);
201  glMultMatrixd(gv->rotate.rotMatrix);
202  glTranslatef(-pos[0], -pos[1], -pos[2]);
203 
204  /* have to redefine clipping planes when view changes */
205 
207 
208  return;
209 }
210 
211 /*!
212  \brief Check focus
213 
214  \param gv view (geoview)
215  */
217 {
218  float zmax, zmin;
219 
220  GS_get_zrange(&zmin, &zmax, 0);
221 
222  if (gv->infocus) {
223  GS_v3eq(gv->from_to[TO], gv->real_to);
224  gv->from_to[TO][Z] -= zmin;
225  GS_v3mult(gv->from_to[TO], gv->scale);
226  gv->from_to[TO][Z] *= gv->vert_exag;
227 
228  GS_v3normalize(gv->from_to[FROM], gv->from_to[TO]);
229  }
230 
231  return;
232 }
233 
234 /*!
235  \brief Get z-up vector (z-direction)
236 
237  \param gv view (geoview)
238  \param up up vector
239  */
240 void gsd_get_zup(geoview *gv, double *up)
241 {
242  float alpha;
243  float zup[3], fup[3];
244 
245  /* neg alpha OK since sin(-x) = -sin(x) */
246  alpha = (2.0 * atan(1.0)) - acos(gv->from_to[FROM][Z] - gv->from_to[TO][Z]);
247 
248  zup[X] = gv->from_to[TO][X];
249  zup[Y] = gv->from_to[TO][Y];
250 
251  if (sin(alpha)) {
252  zup[Z] = gv->from_to[TO][Z] + 1 / sin(alpha);
253  }
254  else {
255  zup[Z] = gv->from_to[FROM][Z] + 1.0;
256  }
257 
258  GS_v3dir(gv->from_to[FROM], zup, fup);
259 
260  up[X] = fup[X];
261  up[Y] = fup[Y];
262  up[Z] = fup[Z];
263 
264  return;
265 }
266 
267 /*!
268  \brief ADD
269 
270  \param gv view (geoview)
271 
272  \return ?
273  */
275 {
276  float fr_to[2][4];
277  float look_theta, pi;
278  float alpha, beta;
279  float zup[3], yup[3], zupmag, yupmag;
280 
281  pi = 4.0 * atan(1.0);
282 
283  /* *************************************************************** */
284  /* This block of code is used to keep pos z in the up direction,
285  * correcting for SGI system default which is pos y in the up
286  * direction. Involves finding up vectors for both y up and
287  * z up, then determining angle between them. LatLon mode uses y as
288  * up direction instead of z, so no correction necessary. Next rewrite,
289  * we should use y as up for all drawing.
290  */
291  GS_v3eq(fr_to[FROM], gv->from_to[FROM]);
292  GS_v3eq(fr_to[TO], gv->from_to[TO]);
293 
294  /* neg alpha OK since sin(-x) = -sin(x) */
295  alpha = pi / 2.0 - acos(fr_to[FROM][Z] - fr_to[TO][Z]);
296 
297  zup[X] = fr_to[TO][X];
298  zup[Y] = fr_to[TO][Y];
299 
300  if (sin(alpha)) {
301  zup[Z] = fr_to[TO][Z] + 1 / sin(alpha);
302  }
303  else {
304  zup[Z] = fr_to[FROM][Z] + 1.0;
305  }
306 
307  zupmag = GS_distance(fr_to[FROM], zup);
308 
309  yup[X] = fr_to[TO][X];
310  yup[Z] = fr_to[TO][Z];
311 
312  /* neg beta OK since sin(-x) = -sin(x) */
313  beta = pi / 2.0 - acos(fr_to[TO][Y] - fr_to[FROM][Y]);
314 
315  if (sin(beta)) {
316  yup[Y] = fr_to[TO][Y] - 1 / sin(beta);
317  }
318  else {
319  yup[Y] = fr_to[FROM][Y] + 1.0;
320  }
321 
322  yupmag = GS_distance(fr_to[FROM], yup);
323 
324  look_theta = (1800.0 / pi) *
325  acos(((zup[X] - fr_to[FROM][X]) * (yup[X] - fr_to[FROM][X]) +
326  (zup[Y] - fr_to[FROM][Y]) * (yup[Y] - fr_to[FROM][Y]) +
327  (zup[Z] - fr_to[FROM][Z]) * (yup[Z] - fr_to[FROM][Z])) /
328  (zupmag * yupmag));
329 
330  if (fr_to[TO][X] - fr_to[FROM][X] < 0.0) {
331  look_theta = -look_theta;
332  }
333 
334  if (fr_to[TO][Z] - fr_to[FROM][Z] < 0.0) {
335  /* looking down */
336  if (fr_to[TO][Y] - fr_to[FROM][Y] < 0.0) {
337  look_theta = 1800 - look_theta;
338  }
339  }
340  else {
341  /* looking up */
342  if (fr_to[TO][Y] - fr_to[FROM][Y] > 0.0) {
343  look_theta = 1800 - look_theta;
344  }
345  }
346 
347  return ((int)(gv->twist + 1800 + look_theta));
348 }
349 
350 /*!
351  \brief Set current scale
352 
353  \param doexag use z-exaggeration
354  */
355 void gsd_do_scale(int doexag)
356 {
357  float sx, sy, sz;
358  float min, max;
359 
360  GS_get_scale(&sx, &sy, &sz, doexag);
361  gsd_scale(sx, sy, sz);
362  GS_get_zrange(&min, &max, 0);
363  gsd_translate(0.0, 0.0, -min);
364 
365  return;
366 }
367 
368 /*!
369  \brief Convert real to model coordinates
370 
371  \param point[in,out] 3d point (Point3)
372  */
374 {
375  float sx, sy, sz;
376  float min, max, n, s, w, e;
377 
378  GS_get_region(&n, &s, &w, &e);
379  GS_get_scale(&sx, &sy, &sz, 1);
380  GS_get_zrange(&min, &max, 0);
381  point[X] = (point[X] - w) * sx;
382  point[Y] = (point[Y] - s) * sy;
383  point[Z] = (point[Z] - min) * sz;
384 
385  return;
386 }
387 
388 /*!
389  \brief Convert model to real coordinates
390 
391  \param point[in,out] 3d point (x,y,z)
392  */
394 {
395  float sx, sy, sz;
396  float min, max, n, s, w, e;
397 
398  GS_get_region(&n, &s, &w, &e);
399  GS_get_scale(&sx, &sy, &sz, 1);
400  GS_get_zrange(&min, &max, 0);
401  point[X] = (sx ? point[X] / sx : 0.0) + w;
402  point[Y] = (sy ? point[Y] / sy : 0.0) + s;
403  point[Z] = (sz ? point[Z] / sz : 0.0) + min;
404 
405  return;
406 }
407 
408 /*!
409  \brief Convert model to surface coordinates
410 
411  \param gs surface (geosurf)
412  \param point 3d point (Point3)
413  */
414 void gsd_model2surf(geosurf *gs, Point3 point)
415 {
416  float min, max, sx, sy, sz;
417 
418  /* so far, only one geographic "region" allowed, so origin of
419  surface is same as origin of model space, but will need to provide
420  translations here to make up the difference, so not using gs yet */
421 
422  if (gs) {
423  /* need to undo z scaling & translate */
424  GS_get_scale(&sx, &sy, &sz, 1);
425  GS_get_zrange(&min, &max, 0);
426 
427  point[Z] = (sz ? point[Z] / sz : 0.0) + min;
428 
429  /* need to unscale x & y */
430  point[X] = (sx ? point[X] / sx : 0.0);
431  point[Y] = (sy ? point[Y] / sy : 0.0);
432  }
433 
434  return;
435 }
436 
437 /*!
438  \brief Convert surface to model coordinates
439 
440  \param point 3d point (Point3)
441  */
443 {
444  float min, max, sx, sy, sz;
445 
446  /* need to undo z scaling & translate */
447  GS_get_scale(&sx, &sy, &sz, 1);
448  GS_get_zrange(&min, &max, 0);
449 
450  point[Z] = (sz ? (point[Z] - min) * sz : 0.0);
451 
452  /* need to unscale x & y */
453  point[X] = (sx ? point[X] * sx : 0.0);
454  point[Y] = (sy ? point[Y] * sy : 0.0);
455 
456  return;
457 }
458 
459 /*!
460  \brief Convert surface to real coordinates
461 
462  \param gs surface (geosurf)
463  \param[in,out] point 3d point (Point3)
464  */
465 void gsd_surf2real(geosurf *gs, Point3 point)
466 {
467  if (gs) {
468  point[X] += gs->ox;
469  point[Y] += gs->oy;
470  }
471 
472  return;
473 }
474 
475 /*!
476  \brief Convert real to surface coordinates
477 
478  \param gs surface (geosurf)
479  \param[in,out] point 3d point (Point3)
480  */
481 void gsd_real2surf(geosurf *gs, Point3 point)
482 {
483  if (gs) {
484  point[X] -= gs->ox;
485  point[Y] -= gs->oy;
486  }
487 
488  return;
489 }
void GS_v3mult(float *, float)
Multiple vectors.
Definition: gs_util.c:229
int GS_get_region(float *, float *, float *, float *)
Get 2D region extent.
Definition: gs2.c:156
void GS_set_draw(int)
Sets which buffer to draw to.
Definition: gs2.c:2459
void gsd_scale(float, float, float)
Multiply the current matrix by a general scaling matrix.
Definition: gsd_prim.c:525
double GS_get_aspect(void)
Get aspect value.
Definition: gs2.c:3447
int GS_get_zrange(float *, float *, int)
Get z-extent for all loaded surfaces.
Definition: gs2.c:2685
int GS_v3normalize(float *, float *)
Change v2 so that v1v2 is a unit vector.
Definition: gs_util.c:321
int gs_get_datacenter(float *)
Get data center point.
Definition: gs.c:1230
void GS_v3eq(float *, float *)
Copy vector values.
Definition: gs_util.c:178
void gsd_color_func(unsigned int)
Set current color.
Definition: gsd_prim.c:698
void gsd_translate(float, float, float)
Multiply the current matrix by a translation matrix.
Definition: gsd_prim.c:539
void GS_ready_draw(void)
Definition: gs2.c:2485
void gsd_linewidth(short)
Set width of rasterized lines.
Definition: gsd_prim.c:267
void GS_get_scale(float *, float *, float *, int)
Get axis scale.
Definition: gs2.c:3236
int GS_v3dir(float *, float *, float *)
Get a normalized direction from v1 to v2, store in v3.
Definition: gs_util.c:351
float GS_distance(float *, float *)
Calculate distance.
Definition: gs_util.c:141
void gsd_update_cplanes(void)
Update cplaces.
Definition: gsd_cplane.c:86
#define min(x, y)
Definition: draw2.c:29
#define max(x, y)
Definition: draw2.c:30
void gsd_real2model(Point3 point)
Convert real to model coordinates.
Definition: gsd_views.c:373
void gsd_surf2real(geosurf *gs, Point3 point)
Convert surface to real coordinates.
Definition: gsd_views.c:465
void gsd_check_focus(geoview *gv)
Check focus.
Definition: gsd_views.c:216
int gsd_zup_twist(geoview *gv)
ADD.
Definition: gsd_views.c:274
int gsd_get_los(float(*vect)[3], short sx, short sy)
ADD.
Definition: gsd_views.c:40
void gsd_real2surf(geosurf *gs, Point3 point)
Convert real to surface coordinates.
Definition: gsd_views.c:481
void gsd_surf2model(Point3 point)
Convert surface to model coordinates.
Definition: gsd_views.c:442
void gsd_model2real(Point3 point)
Convert model to real coordinates.
Definition: gsd_views.c:393
void gsd_model2surf(geosurf *gs, Point3 point)
Convert model to surface coordinates.
Definition: gsd_views.c:414
void gsd_do_scale(int doexag)
Set current scale.
Definition: gsd_views.c:355
void gsd_get_zup(geoview *gv, double *up)
Get z-up vector (z-direction)
Definition: gsd_views.c:240
void gsd_set_view(geoview *gv, geodisplay *gd)
Set view.
Definition: gsd_views.c:146
#define X
Definition: ogsf.h:140
float Point3[3]
Definition: ogsf.h:205
#define Z
Definition: ogsf.h:142
#define GSD_FRONT
Definition: ogsf.h:104
#define Y
Definition: ogsf.h:141
#define GSD_BACK
Definition: ogsf.h:105
#define TO
Definition: ogsf.h:145
@ FROM
Definition: sqlp.tab.h:67
Definition: ogsf.h:256
double ox
Definition: ogsf.h:263
double oy
Definition: ogsf.h:263
float aspect
Definition: ogsf.h:480
float nearclip
Definition: ogsf.h:480
float farclip
Definition: ogsf.h:480
double rot_angle
Definition: ogsf.h:462
GLdouble rotMatrix[16]
Definition: ogsf.h:464
int do_rot
Definition: ogsf.h:461
double rot_axes[3]
Definition: ogsf.h:463
Definition: ogsf.h:467
float vert_exag
Definition: ogsf.h:474
int twist
Definition: ogsf.h:473
float real_to[4]
Definition: ogsf.h:474
struct georot rotate
Definition: ogsf.h:472
float from_to[2][4]
Definition: ogsf.h:471
int fov
Definition: ogsf.h:473
float scale
Definition: ogsf.h:475
int infocus
Definition: ogsf.h:470