GRASS Programmer's Manual  6.4.2(2012)
gsd_views.c
Go to the documentation of this file.
00001 
00019 #include <grass/config.h>
00020 
00021 #if defined(OPENGL_X11) || defined(OPENGL_WINDOWS)
00022 #include <GL/gl.h>
00023 #include <GL/glu.h>
00024 #elif defined(OPENGL_AQUA)
00025 #include <OpenGL/gl.h>
00026 #include <OpenGL/glu.h>
00027 #endif
00028 
00029 #include <grass/gstypes.h>
00030 #include "math.h"
00031 
00040 int gsd_get_los(float (*vect)[3], short sx, short sy)
00041 {
00042     double fx, fy, fz, tx, ty, tz;
00043     GLdouble modelMatrix[16], projMatrix[16];
00044     GLint viewport[4];
00045 
00046     GS_ready_draw();
00047     glPushMatrix();
00048 
00049     gsd_do_scale(1);
00050     glGetDoublev(GL_MODELVIEW_MATRIX, modelMatrix);
00051     glGetDoublev(GL_PROJECTION_MATRIX, projMatrix);
00052     glGetIntegerv(GL_VIEWPORT, viewport);
00053     glPopMatrix();
00054 
00055     /* OGLXXX XXX I think this is backwards gluProject(XXX); */
00056     /* WAS: mapw(Vobj, sx, sy, &fx, &fy, &fz, &tx, &ty, &tz); */
00057     gluUnProject((GLdouble) sx, (GLdouble) sy, 0.0, modelMatrix,
00058                  projMatrix, viewport, &fx, &fy, &fz);
00059     gluUnProject((GLdouble) sx, (GLdouble) sy, 1.0, modelMatrix,
00060                  projMatrix, viewport, &tx, &ty, &tz);
00061     vect[FROM][X] = fx;
00062     vect[FROM][Y] = fy;
00063     vect[FROM][Z] = fz;
00064     vect[TO][X] = tx;
00065     vect[TO][Y] = ty;
00066     vect[TO][Z] = tz;
00067 
00068     /* DEBUG - should just be a dot */
00069     /* OGLXXX frontbuffer: other possibilities include GL_FRONT_AND_BACK */
00070     glDrawBuffer((1) ? GL_FRONT : GL_BACK);
00071     glPushMatrix();
00072     gsd_do_scale(1);
00073     gsd_linewidth(3);
00074     gsd_color_func(0x8888FF);
00075 
00076     /* OGLXXX for multiple, independent line segments: use GL_LINES */
00077     glBegin(GL_LINE_STRIP);
00078     glVertex3fv(vect[FROM]);
00079     glVertex3fv(vect[TO]);
00080     glEnd();
00081 
00082     gsd_linewidth(1);
00083     glPopMatrix();
00084 
00085     /* OGLXXX frontbuffer: other possibilities include GL_FRONT_AND_BACK */
00086     glDrawBuffer((0) ? GL_FRONT : GL_BACK);
00087 
00088     return (1);
00089 }
00090 
00091 
00100 void gsd_set_view(geoview * gv, geodisplay * gd)
00101 {
00102     double up[3];
00103     GLint mm;
00104 
00105     /* will expand when need to check for in focus, ortho, etc. */
00106 
00107     gsd_check_focus(gv);
00108     gsd_get_zup(gv, up);
00109 
00110     gd->aspect = GS_get_aspect();
00111 
00112     glGetIntegerv(GL_MATRIX_MODE, &mm);
00113     glMatrixMode(GL_PROJECTION);
00114     glLoadIdentity();
00115     gluPerspective((double).1 * (gv->fov), (double)gd->aspect,
00116                    (double)gd->nearclip, (double)gd->farclip);
00117 
00118     glMatrixMode(mm);
00119 
00120     glLoadIdentity();
00121 
00122     /* update twist parm */
00123     glRotatef((float)(gv->twist / 10.), 0.0, 0.0, 1.0);
00124 
00125     /* OGLXXX lookat: replace UPx with vector */
00126     gluLookAt((double)gv->from_to[FROM][X], (double)gv->from_to[FROM][Y],
00127               (double)gv->from_to[FROM][Z], (double)gv->from_to[TO][X],
00128               (double)gv->from_to[TO][Y], (double)gv->from_to[TO][Z],
00129               (double)up[X], (double)up[Y], (double)up[Z]);
00130 
00131     /* have to redefine clipping planes when view changes */
00132 
00133     gsd_update_cplanes();
00134 
00135     return;
00136 }
00137 
00143 void gsd_check_focus(geoview * gv)
00144 {
00145     float zmax, zmin;
00146 
00147     GS_get_zrange(&zmin, &zmax, 0);
00148 
00149     if (gv->infocus) {
00150         GS_v3eq(gv->from_to[TO], gv->real_to);
00151         gv->from_to[TO][Z] -= zmin;
00152         GS_v3mult(gv->from_to[TO], gv->scale);
00153         gv->from_to[TO][Z] *= gv->vert_exag;
00154 
00155         GS_v3normalize(gv->from_to[FROM], gv->from_to[TO]);
00156     }
00157 
00158     return;
00159 }
00160 
00167 void gsd_get_zup(geoview * gv, double *up)
00168 {
00169     float alpha;
00170     float zup[3], fup[3];
00171 
00172     /* neg alpha OK since sin(-x) = -sin(x) */
00173     alpha =
00174         (2.0 * atan(1.0)) - acos(gv->from_to[FROM][Z] - gv->from_to[TO][Z]);
00175 
00176     zup[X] = gv->from_to[TO][X];
00177     zup[Y] = gv->from_to[TO][Y];
00178 
00179     if (sin(alpha)) {
00180         zup[Z] = gv->from_to[TO][Z] + 1 / sin(alpha);
00181     }
00182     else {
00183         zup[Z] = gv->from_to[FROM][Z] + 1.0;
00184     }
00185 
00186     GS_v3dir(gv->from_to[FROM], zup, fup);
00187 
00188     up[X] = fup[X];
00189     up[Y] = fup[Y];
00190     up[Z] = fup[Z];
00191 
00192     return;
00193 }
00194 
00202 int gsd_zup_twist(geoview * gv)
00203 {
00204     float fr_to[2][4];
00205     float look_theta, pi;
00206     float alpha, beta;
00207     float zup[3], yup[3], zupmag, yupmag;
00208 
00209     pi = 4.0 * atan(1.0);
00210 
00211     /* *************************************************************** */
00212     /* This block of code is used to keep pos z in the up direction,
00213      * correcting for SGI system default which is pos y in the up
00214      * direction.  Involves finding up vectors for both y up and
00215      * z up, then determining angle between them.  LatLon mode uses y as
00216      * up direction instead of z, so no correction necessary.  Next rewrite,
00217      * we should use y as up for all drawing.
00218      */
00219     GS_v3eq(fr_to[FROM], gv->from_to[FROM]);
00220     GS_v3eq(fr_to[TO], gv->from_to[TO]);
00221 
00222     /* neg alpha OK since sin(-x) = -sin(x) */
00223     alpha = pi / 2.0 - acos(fr_to[FROM][Z] - fr_to[TO][Z]);
00224 
00225     zup[X] = fr_to[TO][X];
00226     zup[Y] = fr_to[TO][Y];
00227 
00228     if (sin(alpha)) {
00229         zup[Z] = fr_to[TO][Z] + 1 / sin(alpha);
00230     }
00231     else {
00232         zup[Z] = fr_to[FROM][Z] + 1.0;
00233     }
00234 
00235     zupmag = GS_distance(fr_to[FROM], zup);
00236 
00237     yup[X] = fr_to[TO][X];
00238     yup[Z] = fr_to[TO][Z];
00239 
00240     /* neg beta OK since sin(-x) = -sin(x) */
00241     beta = pi / 2.0 - acos(fr_to[TO][Y] - fr_to[FROM][Y]);
00242 
00243     if (sin(beta)) {
00244         yup[Y] = fr_to[TO][Y] - 1 / sin(beta);
00245     }
00246     else {
00247         yup[Y] = fr_to[FROM][Y] + 1.0;
00248     }
00249 
00250     yupmag = GS_distance(fr_to[FROM], yup);
00251 
00252     look_theta = (1800.0 / pi) *
00253         acos(((zup[X] - fr_to[FROM][X]) * (yup[X] - fr_to[FROM][X])
00254               + (zup[Y] - fr_to[FROM][Y]) * (yup[Y] - fr_to[FROM][Y])
00255               + (zup[Z] - fr_to[FROM][Z]) * (yup[Z] - fr_to[FROM][Z])) /
00256              (zupmag * yupmag));
00257 
00258     if (fr_to[TO][X] - fr_to[FROM][X] < 0.0) {
00259         look_theta = -look_theta;
00260     }
00261 
00262     if (fr_to[TO][Z] - fr_to[FROM][Z] < 0.0) {
00263         /* looking down */
00264         if (fr_to[TO][Y] - fr_to[FROM][Y] < 0.0) {
00265             look_theta = 1800 - look_theta;
00266         }
00267     }
00268     else {
00269         /* looking up */
00270         if (fr_to[TO][Y] - fr_to[FROM][Y] > 0.0) {
00271             look_theta = 1800 - look_theta;
00272         }
00273     }
00274 
00275     return ((int)(gv->twist + 1800 + look_theta));
00276 }
00277 
00283 void gsd_do_scale(int doexag)
00284 {
00285     float sx, sy, sz;
00286     float min, max;
00287 
00288     GS_get_scale(&sx, &sy, &sz, doexag);
00289     gsd_scale(sx, sy, sz);
00290     GS_get_zrange(&min, &max, 0);
00291     gsd_translate(0.0, 0.0, -min);
00292 
00293     return;
00294 }
00295 
00301 void gsd_real2model(Point3 point)
00302 {
00303     float sx, sy, sz;
00304     float min, max, n, s, w, e;
00305 
00306     GS_get_region(&n, &s, &w, &e);
00307     GS_get_scale(&sx, &sy, &sz, 1);
00308     GS_get_zrange(&min, &max, 0);
00309     point[X] = (point[X] - w) * sx;
00310     point[Y] = (point[Y] - s) * sy;
00311     point[Z] = (point[Z] - min) * sz;
00312 
00313     return;
00314 }
00315 
00321 void gsd_model2real(Point3 point)
00322 {
00323     float sx, sy, sz;
00324     float min, max, n, s, w, e;
00325 
00326     GS_get_region(&n, &s, &w, &e);
00327     GS_get_scale(&sx, &sy, &sz, 1);
00328     GS_get_zrange(&min, &max, 0);
00329     point[X] = (sx ? point[X] / sx : 0.0) + w;
00330     point[Y] = (sy ? point[Y] / sy : 0.0) + s;
00331     point[Z] = (sz ? point[Z] / sz : 0.0) + min;
00332 
00333     return;
00334 }
00335 
00342 void gsd_model2surf(geosurf * gs, Point3 point)
00343 {
00344     float min, max, sx, sy, sz;
00345 
00346     /* so far, only one geographic "region" allowed, so origin of
00347        surface is same as origin of model space, but will need to provide 
00348        translations here to make up the difference, so not using gs yet */
00349 
00350     if (gs) {
00351         /* need to undo z scaling & translate */
00352         GS_get_scale(&sx, &sy, &sz, 1);
00353         GS_get_zrange(&min, &max, 0);
00354 
00355         point[Z] = (sz ? point[Z] / sz : 0.0) + min;
00356 
00357         /* need to unscale x & y */
00358         point[X] = (sx ? point[X] / sx : 0.0);
00359         point[Y] = (sy ? point[Y] / sy : 0.0);
00360     }
00361 
00362     return;
00363 }
00364 
00371 void gsd_surf2real(geosurf * gs, Point3 point)
00372 {
00373     if (gs) {
00374         point[X] += gs->ox;
00375         point[Y] += gs->oy;
00376     }
00377 
00378     return;
00379 }
00380 
00387 void gsd_real2surf(geosurf * gs, Point3 point)
00388 {
00389     if (gs) {
00390         point[X] -= gs->ox;
00391         point[Y] -= gs->oy;
00392     }
00393 
00394     return;
00395 }
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines