Geografic-Information-System/lib/ogsf/GS2.c

/*!
  \file GS2.c
 
  \brief OGSF library - loading and manipulating surfaces (higher level functions)
 
  GRASS OpenGL gsurf OGSF Library 

  Plans for handling color maps:
  NOW:
  if able to load as unsigned char, make lookup table containing palette
  otherwise, load directly as packed color, set lookup = NULL
  MAYBE LATER:
  if able to load as POSITIVE short, make lookup table containing palette
  - may want to calculate savings first (ie,  numcells > 32768)
  (not exactly, it's Friday & time to go home - figure it later)
  otherwise, load directly as packed color, set lookup = NULL
  MESSY! - need to fix up!
 
  (C) 1999-2008 by the GRASS Development Team
 
  This program is free software under the 
  GNU General Public License (>=v2). 
  Read the file COPYING that comes with GRASS
  for details.
  
  \author Bill Brown USACERL (1993)
  \author Pierre de Mouveaux <p_de_mouveaux hotmail.com> (updated October 1999)
  \author Doxygenized by Martin Landa <landa.martin gmail.com> (May 2008)
*/

#include <stdlib.h>
#include <string.h>
#include <math.h>

#include <grass/config.h>

#if defined(OPENGL_X11) || defined(OPENGL_WINDOWS)
#include <GL/gl.h>
#include <GL/glu.h>
#elif defined(OPENGL_AQUA)
#include <OpenGL/gl.h>
#include <OpenGL/glu.h>
#endif

#include <grass/gis.h>
#include <grass/gstypes.h>
#include <grass/glocale.h>

#include "gsget.h"
#include "rowcol.h"
#include "rgbpack.h"

/* Hack to make NVIZ2.2 query functions.("What's Here" and "Look at")
 * to work.
 * Uses gs_los_intersect1() instead of gs_los_intersect().
 * Pierre de Mouveaux - 31 oct. 1999. p_de_mouveaux@hotmail.com.
 */
#define NVIZ_HACK 1

int gsd_getViewport(GLint *, GLint *);

/* array of surface ids */
static int Surf_ID[MAX_SURFS];
static int Next_surf = 0;
static int SDref_surf = 0;
/* attributes array */
static float Default_const[MAX_ATTS];
static float Default_nulls[MAX_ATTS];
/* largest dimension */
static float Longdim;
/* N, S, W, E */
static float Region[4];
static geoview Gv;
static geodisplay Gd;
static struct Cell_head wind;
static int Buffermode;
static int Numlights = 0;
static int Modelshowing = 0;

void void_func(void)
{
    return;
}

/*!
  \brief Initialize OGSF library

  Get region settings - wind

  Set Region (NSWE array) and compute scale
*/
void GS_libinit(void)
{
    static int first = 1;

    G_get_set_window(&wind);

    Region[0] = wind.north;
    Region[1] = wind.south;
    Region[2] = wind.west;
    Region[3] = wind.east;

    /* scale largest dimension to GS_UNIT_SIZE */
    if ((wind.east - wind.west) > (wind.north - wind.south)) {
	Longdim = (wind.east - wind.west);
    }
    else {
	Longdim = (wind.north - wind.south);
    }

    Gv.scale = GS_UNIT_SIZE / Longdim;

    Cxl_func = void_func;
    Swap_func = void_func;

    if (first) {
	gs_init();
    }

    first = 0;

    return;
}

/*!
  \brief Get largest dimension

  \param[out] dim dimension

  \return 1
*/
int GS_get_longdim(float *dim)
{
    *dim = Longdim;

    G_debug(3, "GS_get_longdim(): dim=%g", *dim);

    return (1);
}

/*!
  \brief Get 2D region extent

  \param[out] n,s,w,e extent values

  \return 1
*/
int GS_get_region(float *n, float *s, float *w, float *e)
{
    *n = Region[0];
    *s = Region[1];
    *w = Region[2];
    *e = Region[3];

    return (1);
}

/*!
  \brief Set default attributes for map objects

  \param defs attributes array (dim MAX_ATTS)
  \param null_defs null attributes array (dim MAX_ATTS)
*/
void GS_set_att_defaults(float *defs, float *null_defs)
{
    int i;

    G_debug (3, "GS_set_att_defaults");

    for (i = 0; i < MAX_ATTS; i++) {
	Default_const[i] = defs[i];
	Default_nulls[i] = null_defs[i];
    }

    return;
}

/*!
  Check if surface exists

  \param id surface id

  \return 0 not found
  \return 1 found
*/
int GS_surf_exists(int id)
{
    int i, found = 0;

    G_debug(3, "GS_surf_exists");


    if (NULL == gs_get_surf(id)) {
	return (0);
    }

    for (i = 0; i < Next_surf && !found; i++) {
	if (Surf_ID[i] == id) {
	    found = 1;
	}
    }

    return (found);
}

/*!
  \brief Add new surface

  Note that origin has 1/2 cell added to represent center of cells
  because library assumes that east - west = (cols - 1) * ew_res,
  since left and right columns are on the edges.

  \return surface id
  \return -1 on error (MAX_SURFS exceded)
*/
int GS_new_surface(void)
{
    geosurf *ns;

    if (Next_surf < MAX_SURFS) {
	ns = gs_get_new_surface();
	gs_init_surf(ns, wind.west + wind.ew_res / 2.,
		     wind.south + wind.ns_res / 2., wind.rows, wind.cols,
		     wind.ew_res, wind.ns_res);
	gs_set_defaults(ns, Default_const, Default_nulls);

	/* make default shine current */
	gs_set_att_src(ns, ATT_SHINE, CONST_ATT);

	Surf_ID[Next_surf] = ns->gsurf_id;
	++Next_surf;

	G_debug(3, "GS_new_surface(): id=%d", ns->gsurf_id);
	
	return (ns->gsurf_id);
    }



    return (-1);
}

/*!
  \brief Add new model light

  \return light model id
  \return -1 on error (MAX_LIGHTS exceded)
*/
int GS_new_light(void)
{
    static int first = 1;
    int i;

    if (first) {
	first = 0;

	for (i = 0; i < MAX_LIGHTS; i++) {
	    Gv.lights[i].position[X] = Gv.lights[i].position[Y] = 0.0;
	    Gv.lights[i].position[Z] = 1.0;
	    Gv.lights[i].position[W] = 0.0;	/* infinite */
	    Gv.lights[i].color[0] = Gv.lights[i].color[1] =
		Gv.lights[i].color[2] = 1.0;
	    Gv.lights[i].ambient[0] = Gv.lights[i].ambient[1] =
		Gv.lights[i].ambient[2] = 0.2;
	    Gv.lights[i].shine = 32.0;
	}

	gsd_init_lightmodel();
    }

    if (Numlights < MAX_LIGHTS) {
	gsd_deflight(Numlights + 1, &(Gv.lights[Numlights]));
	gsd_switchlight(Numlights + 1, 1);

	return (++Numlights);
    }

    return (-1);
}

/*!
  \brief Set light position

  \bug I think lights array doesnt match sgi_light array
  
  \param num light id (starts with 1)
  \param xpos,ypos,zpos coordinates (model)
  \param local local coordinate (for viewport)
*/
void GS_setlight_position(int num, float xpos, float ypos, float zpos,
			  int local)
{
    if (num) {
	num -= 1;
	if (num < Numlights) {
	    Gv.lights[num].position[X] = xpos;
	    Gv.lights[num].position[Y] = ypos;
	    Gv.lights[num].position[Z] = zpos;
	    Gv.lights[num].position[W] = (float) local;

	    gsd_deflight(num + 1, &(Gv.lights[num]));
	}
    }

    return;
}


/*!
  \brief Get light position

  \param num light id (starts with 1?)
  \param[out] xpos,ypos,zpos coordinates
  \param[out] local ?
*/
void GS_getlight_position(int num, float *xpos, float *ypos, float *zpos,
			  int *local)
{
    if (num) {
	num -= 1;
	if (num < Numlights) {
	    *xpos = Gv.lights[num].position[X];
	    *ypos = Gv.lights[num].position[Y];
	    *zpos = Gv.lights[num].position[Z];
	    *local = (int) Gv.lights[num].position[W];

	}
    }

    return;
}

/*!
  \brief Set light color
  
  \param num light id (starts with 1?)
  \param red,green,blue color values (from 0.0 to 1.0)
*/
void GS_setlight_color(int num, float red, float green, float blue)
{
    if (num) {
	num -= 1;
	if (num < Numlights) {
	    Gv.lights[num].color[0] = red;
	    Gv.lights[num].color[1] = green;
	    Gv.lights[num].color[2] = blue;

	    gsd_deflight(num + 1, &(Gv.lights[num]));
	}
    }

    return;
}

/*!
  \brief Get light color

  \param num light id (starts with 1?)
  \param[out] red,green,blue color values
*/
void GS_getlight_color(int num, float *red, float *green, float *blue)
{
    if (num) {
	num -= 1;
	if (num < Numlights) {
	    *red = Gv.lights[num].color[0];
	    *green = Gv.lights[num].color[1];
	    *blue = Gv.lights[num].color[2];
	}
    }

    return;
}

/*!
  \brief Set light ambient

  Red, green, blue from 0.0 to 1.0

  \param num light id (starts with 1?)
  \param red,green,blue color values
*/
void GS_setlight_ambient(int num, float red, float green, float blue)
{
    if (num) {
	num -= 1;
	if (num < Numlights) {
	    Gv.lights[num].ambient[0] = red;
	    Gv.lights[num].ambient[1] = green;
	    Gv.lights[num].ambient[2] = blue;

	    gsd_deflight(num + 1, &(Gv.lights[num]));
	}
    }

    return;
}

/*!
  \brief Get light ambient

  \param num light id (starts with 1?)
  \param[out] red,green,blue color values
*/
void GS_getlight_ambient(int num, float *red, float *green, float *blue)
{
    if (num) {
	num -= 1;
	if (num < Numlights) {
	    *red = Gv.lights[num].ambient[0];
	    *green = Gv.lights[num].ambient[1];
	    *blue = Gv.lights[num].ambient[2];
	}
    }

    return;
}


/*!
  \brief Switch off all lights
*/
void GS_lights_off(void)
{
    int i;

    for (i = 0; i < Numlights; i++) {
	gsd_switchlight(i + 1, 0);
    }

    return;
}

/*!
  \brief Switch on all lights
*/
void GS_lights_on(void)
{
    int i;

    for (i = 0; i < Numlights; i++) {
	gsd_switchlight(i + 1, 1);
    }

    return;
}

/*!
  \brief Switch on/off light

  \param num light id (starts with 1?)
  \param on non-zero for 'on' otherwise 'off'
*/
void GS_switchlight(int num, int on)
{
    if (num) {
	num -= 1;

	if (num < Numlights) {
	    gsd_switchlight(num + 1, on);
	}
    }

    return;
}

/*!
  \brief Check if transparency is set

  \return 0 transparency not set
  \return 1 transparency is set
*/
int GS_transp_is_set(void)
{
    return (gs_att_is_set(NULL, ATT_TRANSP) || (FC_GREY == gsd_getfc()));
}

/*!
  \brief Retrieves coordinates for lighting model position, at center of view

  \param pos[out] coordinates
*/
void GS_get_modelposition1(float pos[])
{
    /* TODO: Still needs work to handle other cases */
    /* this is a quick hack to get lighting adjustments debugged */
    /*
       GS_v3dir(Gv.from_to[FROM], Gv.from_to[TO], center);
       GS_v3mult(center, 1000);
       GS_v3add(center, Gv.from_to[FROM]);
     */

    gs_get_datacenter(pos);
    gs_get_data_avg_zmax(&(pos[Z]));

    G_debug (1, "GS_get_modelposition1(): model position: %f %f %f",
	     pos[X], pos[Y], pos[Z]);

    return;
}

/*!
  \brief Retrieves coordinates for lighting model position, at center of view

  Position at nearclip * 2: tried nearclip + siz, but since need to
  know position to calculate size, have two dependent variables
  (nearclip * 2) from eye.

  \param siz[out] size
  \param pos[out] coordinates (X, Y, Z)
*/
void GS_get_modelposition(float *siz, float *pos)
{
    float dist, near_h, dir[3];

    dist = 2. * Gd.nearclip;

    near_h = 2.0 * tan(4.0 * atan(1.) * Gv.fov / 3600.) * dist;
    *siz = near_h / 8.0;

    /* prevent clipping - would only happen if fov > ~127 degrees, at
       fov = 2.0 * atan(2.0) */

    if (*siz > Gd.nearclip) {
	*siz = Gd.nearclip;
    }

    GS_v3dir(Gv.from_to[FROM], Gv.from_to[TO], dir);

    pos[X] = Gv.from_to[FROM][X] + dir[X] * dist;
    pos[Y] = Gv.from_to[FROM][Y] + dir[Y] * dist;
    pos[Z] = Gv.from_to[FROM][Z] + dir[Z] * dist;

    return;
}


/*!
  \brief Set decoration, north arrow ??

  \todo scale used to calculate len of arrow still needs work
  needs go function that returns center / eye distance
  gsd_get_los function is not working correctly ??

  \param pt point value in true world coordinates (?)
  \param id surface id
  \param[out] pos2 output coordinates
*/
void GS_set_Narrow(int *pt, int id, float *pos2)
{
    geosurf *gs;
    float x, y, z;
    GLdouble modelMatrix[16], projMatrix[16];
    GLint viewport[4];

    if (GS_get_selected_point_on_surface(pt[X], pt[Y], &id, &x, &y, &z)) {
        gs = gs_get_surf(id);
        if (gs) {
            z = gs->zmax;
            pos2[X] = (float)x - gs->ox + gs->x_trans;
            pos2[Y] = (float)y - gs->oy + gs->y_trans;
            pos2[Z] = (float)z + gs->z_trans;

            return;
        }
    }
    else {
        gs = gs_get_surf(id);
        
        /* Need to get model matrix, etc 
         * to run gluUnProject
         */
        gsd_pushmatrix();
        gsd_do_scale(1);
        glGetDoublev(GL_MODELVIEW_MATRIX, modelMatrix);
        glGetDoublev(GL_PROJECTION_MATRIX, projMatrix);
        glGetIntegerv(GL_VIEWPORT, viewport);

        if (gs) {
        GLdouble out_near[3], out_far[3];
        GLdouble factor;
        GLdouble out[3];

            z = (float)gs->zmax + gs->z_trans;

            gluUnProject((GLdouble)pt[X], (GLdouble)pt[Y], (GLdouble)0.,
                            modelMatrix, projMatrix, viewport,
                            &out_near[X], &out_near[Y], &out_near[Z]);
            gluUnProject((GLdouble)pt[X], (GLdouble)pt[Y], (GLdouble)1.,
                            modelMatrix, projMatrix, viewport,
                            &out_far[X], &out_far[Y], &out_far[Z]);
                            
            glPopMatrix();

            factor = (out_near[Z]-z) / (out_near[Z]-out_far[Z]);

            out[X] = out_near[X]-((out_near[X]-out_far[X])*factor);
            out[Y] = out_near[Y]-((out_near[Y]-out_far[Y])*factor);
            out[Z] = z;

            pos2[X] = (float)out[X];
            pos2[Y] = (float)out[Y];
            pos2[Z] = (float)out[Z];

            return;

        }
    }
    return;
}

/*!
  \brief ADD

  \param id surface id
  \param pt point, X, Y value in true world coordinates
*/
void GS_draw_X(int id, float *pt)
{
    geosurf *gs;
    Point3 pos;
    float siz;

    if ((gs = gs_get_surf(id))) {
	GS_get_longdim(&siz);
	siz /= 200.;
	pos[X] = pt[X] - gs->ox;
	pos[Y] = pt[Y] - gs->oy;
	_viewcell_tri_interp(gs, pos);

	gsd_pushmatrix();

	gsd_do_scale(1);
	gsd_translate(gs->x_trans, gs->y_trans, gs->z_trans);
	gsd_linewidth(1);

	if (CONST_ATT == gs_get_att_src(gs, ATT_TOPO)) {
	    pos[Z] = gs->att[ATT_TOPO].constant;
	    gs = NULL;		/* tells gpd_obj to use given Z val */
	}

	gpd_obj(gs, Gd.bgcol, siz, ST_GYRO, pos);
	gsd_flush();

	gsd_popmatrix();
    }

    return;
}

/*!
  \brief Draw line on surface

  \param id surface id
  \param x1,y1,x2,y2 line nodes
*/
void GS_draw_line_onsurf(int id, float x1, float y1, float x2, float y2)
{
    float p1[2], p2[2];
    geosurf *gs;

    if ((gs = gs_get_surf(id))) {
	p1[X] = x1 - gs->ox;
	p1[Y] = y1 - gs->oy;
	p2[X] = x2 - gs->ox;
	p2[Y] = y2 - gs->oy;

	gsd_pushmatrix();

	gsd_do_scale(1);
	gsd_translate(gs->x_trans, gs->y_trans, gs->z_trans);
	gsd_linewidth(1);

	gsd_color_func(GS_default_draw_color());
	gsd_line_onsurf(gs, p1, p2);

	gsd_popmatrix();
	gsd_flush();
    }

    return;
}

/*!
  \brief Draw multiline on surface

  Like above but limits points in line to n or points found in segment,
  whichever is smaller.

  \param id surface id
  \param x1,y1,x2,y2 line nodes
  \param lasp ?
  \paran n ?
  \return number of points used
*/
int GS_draw_nline_onsurf(int id, float x1, float y1, float x2, float y2,
			 float *lasp, int n)
{
    float p1[2], p2[2];
    geosurf *gs;
    int ret = 0;

    if ((gs = gs_get_surf(id))) {
	p1[X] = x1 - gs->ox;
	p1[Y] = y1 - gs->oy;
	p2[X] = x2 - gs->ox;
	p2[Y] = y2 - gs->oy;

	gsd_pushmatrix();

	gsd_do_scale(1);
	gsd_translate(gs->x_trans, gs->y_trans, gs->z_trans);
	gsd_linewidth(1);
	gsd_color_func(GS_default_draw_color());
	ret = gsd_nline_onsurf(gs, p1, p2, lasp, n);
	gsd_surf2real(gs, lasp);

	gsd_popmatrix();
	gsd_flush();
    }

    return (ret);
}

/*!
  \brief Draw flow-line on surace

  This is slow - should be moved to gs_ but GS_ good for testing
  and useful for app programmer

  \param id surface id
  \param x,y coordinates of flow-line
*/
void GS_draw_flowline_at_xy(int id, float x, float y)
{
    geosurf *gs;
    float nv[3], pdir[2], mult;
    float p1[2], p2[2], next[2];
    int i = 0;

    if ((gs = gs_get_surf(id))) {
	p1[X] = x;
	p1[Y] = y;
	/* multiply by 1.5 resolutions to ensure a crossing ? */
	mult = .1 * (VXRES(gs) > VYRES(gs) ? VXRES(gs) : VYRES(gs));

	GS_coordpair_repeats(p1, p1, 50);

	while (1 == GS_get_norm_at_xy(id, p1[X], p1[Y], nv)) {
	    if (nv[Z] == 1.0) {
		if (pdir[X] == 0.0 && pdir[Y] == 0.0) {
		    break;
		}

		p2[X] = p1[X] + (pdir[X] * mult);
		p2[Y] = p1[Y] + (pdir[Y] * mult);
	    }
	    else {
		/* use previous direction */
		GS_v2norm(nv);
		p2[X] = p1[X] + (nv[X] * mult);
		p2[Y] = p1[Y] + (nv[Y] * mult);
		pdir[X] = nv[X];
		pdir[Y] = nv[Y];
	    }

	    if (i > 2000) {
		break;
	    }

	    if (GS_coordpair_repeats(p1, p2, 0)) {
		break;
	    }

	    /* Think about this: */
	    /* degenerate line means edge or level edge ? */
	    /* next is filled with last point drawn */
	    if (2 > GS_draw_nline_onsurf(id, p1[X], p1[Y],
					 p2[X], p2[Y], next, 3)) {
		break;
	    }

	    p1[X] = next[X];
	    p1[Y] = next[Y];
	}

	G_debug(3, "GS_draw_flowline_at_xy(): dir: %f %f",
		nv[X], nv[Y]);
    }

    return;
}

/*!
  \brief Draw fringe around data at selected corners

  \param id surface id
  \param clr ?
  \param elev elevation value
  \param where ?
*/
void GS_draw_fringe(int id, unsigned long clr, float elev, int *where)
{
  geosurf *gs;

  if ((gs = gs_get_surf(id)))
             gsd_display_fringe(gs, clr, elev, where);

}


/*!
  \brief Draw legend

  \todo add legend from list option
  make font loading more flexible

  \param name legend name
  \param fontbase font-base
  \param size ? 
  \param flags legend flags
  \param range values range
  \param pt ?
*/
int GS_draw_legend(char *name, GLuint fontbase, int size, int *flags,
	       float *range, int *pt)
{
    int list_no;

    list_no = gsd_put_legend(name, fontbase, size, flags, range, pt);

    return (list_no);
}

/*!
  \brief Draw pre-defined list

  Uses glFlush() to ensure all drawing is complete
  before returning

  \param list_id list id
*/
void GS_draw_list(GLuint list_id)
{
    gsd_calllist(list_id);
    glFlush();
    return;
}

/*!
  \brief Draw all glLists

  Uses glFlush() to ensure all drawing is complete
  before returning
*/
void GS_draw_all_list(void)
{
    gsd_calllists(0); /* not sure if 0 is right - MN */
    glFlush();
    return;
}

/*!
  \brief Delete pre-defined list

  \param list_id list id
*/
void GS_delete_list(GLuint list_id)
{
    gsd_deletelist(list_id, 1);

    return;
}

/*!
  \brief Draw lighting model
*/
void GS_draw_lighting_model1(void)
{
    static float center[3];
    float tcenter[3];

    if (!Modelshowing) {
	GS_get_modelposition1(center);
    }

    GS_v3eq(tcenter, center);

    gsd_zwritemask(0x0);
    gsd_backface(1);

    gsd_colormode(CM_AD);
    gsd_shademodel(DM_GOURAUD);
    gsd_pushmatrix();
    gsd_do_scale(1);

    if (Gv.vert_exag) {
	tcenter[Z] *= Gv.vert_exag;
	gsd_scale(1.0, 1.0, 1. / Gv.vert_exag);
    }

    gsd_drawsphere(tcenter, 0xDDDDDD, (float) (Longdim / 10.));
    gsd_popmatrix();
    Modelshowing = 1;

    gsd_backface(0);
    gsd_zwritemask(0xffffffff);

    return;
}

/*!
  \brief Draw lighting model

  Just turn off any cutting planes and draw it just outside near
  clipping plane, since lighting is infinite now
*/
void GS_draw_lighting_model(void)
{
    static float center[3], size;
    float tcenter[3], tsize;
    int i, wason[MAX_CPLANES];

    gsd_get_cplanes_state(wason);

    for (i = 0; i < MAX_CPLANES; i++) {
	if (wason[i]) {
	    gsd_cplane_off(i);
	}
    }


    if (!Modelshowing) {
	GS_get_modelposition(&size, center);
    }

    GS_v3eq(tcenter, center);
    tsize = size;

    gsd_zwritemask(0x0);
    gsd_backface(1);

    gsd_colormode(CM_DIFFUSE);
    gsd_shademodel(DM_GOURAUD);
    gsd_pushmatrix();
    gsd_drawsphere(tcenter, 0xDDDDDD, tsize);
    gsd_popmatrix();
    Modelshowing = 1;

    gsd_backface(0);
    gsd_zwritemask(0xffffffff);

    for (i = 0; i < MAX_CPLANES; i++) {
	if (wason[i]) {
	    gsd_cplane_on(i);
	}
    }

    gsd_flush();

    return;
}

/*!
  \brief Update current mask

  May be called to update total mask for a surface at convenient times
  instead of waiting until ready to redraw surface

  \param id surface id

  \return ?
*/
int GS_update_curmask(int id)
{
    geosurf *gs;

    gs = gs_get_surf(id);
    return (gs_update_curmask(gs));
}

/*!
  \brief Check if point is masked ?

  \param id surface id
  \param pt point

  \return 1 masked
  \return 0 not masked
  \return -1 on error, invalid surface id
*/
int GS_is_masked(int id, float *pt)
{
    geosurf *gs;
    Point3 tmp;

    if ((gs = gs_get_surf(id))) {
	tmp[X] = pt[X] - gs->ox;
	tmp[Y] = pt[Y] - gs->oy;

	return (gs_point_is_masked(gs, tmp));
    }

    return (-1);
}

/*!
  \brief Unset ?
*/
void GS_unset_SDsurf(void)
{
    gsdiff_set_SDref(NULL);
    SDref_surf = 0;

    return;
}

/*!
  \brief Set ?

  \param id surface id

  \return 1 on success
  \return 0 on error, invalid surface id
*/
int GS_set_SDsurf(int id)
{
    geosurf *gs;

    if ((gs = gs_get_surf(id))) {
	gsdiff_set_SDref(gs);
	SDref_surf = id;

	return (1);
    }

    return (0);
}

/*!
  \brief Set ?
  
  \param scale scale value

  \return 1
*/
int GS_set_SDscale(float scale)
{
    gsdiff_set_SDscale(scale);

    return (1);
}

/*!
  \brief Get ?

  \param[out] id ?

  \return 1 on success
  \return 0 on error
*/
int GS_get_SDsurf(int *id)
{
    geosurf *gs;

    if ((gs = gsdiff_get_SDref())) {
	*id = SDref_surf;

	return (1);
    }

    return (0);
}

/*!
  \brief Get ?

  \param[out] scale value

  \return 1
*/
int GS_get_SDscale(float *scale)
{
    *scale = gsdiff_get_SDscale();

    return (1);
}

/*!
  \brief Update normals

  \param id surface id

  \return ?
*/
int GS_update_normals(int id)
{
    geosurf *gs;

    gs = gs_get_surf(id);

    return (gs_calc_normals(gs));
}

/*!
  \brief Get attributes

  \param id surface id
  \param att
  \param[out] set
  \param[out] constant
  \param[out] mapname

  \return 1 on success
  \return -1 on error (invalid surface id)
*/
int GS_get_att(int id, int att, int *set, float *constant, char *mapname)
{
    int src;
    geosurf *gs;

    gs = gs_get_surf(id);
    if (gs) {
	if (-1 != (src = gs_get_att_src(gs, att))) {
	    *set = src;

	    if (src == CONST_ATT) {
		*constant = gs->att[att].constant;
	    }
	    else if (src == MAP_ATT) {
		strcpy(mapname, gsds_get_name(gs->att[att].hdata));
	    }

	    return (1);
	}

	return (-1);
    }

    return (-1);
}

#define CATSREADY
#ifdef CATSREADY
/*!
  \brief Print categories on given position

  prints "no data" or a description (i.e., "coniferous forest") to catstr
  Usually call after GS_get_selected_point_on_surface
  att_src must be MAP_ATT

  \param id surface id
  \param att
  \param catstr cat string (must be allocated, dim?)
  \param x,y real coordinates

  \return -1 if no category info or point outside of window
  \return 1 on success
*/
int GS_get_cat_at_xy(int id, int att, char *catstr, float x, float y)
{
    int offset, drow, dcol, vrow, vcol;
    float ftmp, pt[3];
    typbuff *buff;
    geosurf *gs;

    catstr[0] = '\0';
    gs = gs_get_surf(id);

    if (NULL == gs) {
	return (-1);
    }

    pt[X] = x;
    pt[Y] = y;

    gsd_real2surf(gs, pt);
    if (gs_point_is_masked(gs, pt)) {
	return (-1);
    }

    if (!in_vregion(gs, pt)) {
	return (-1);
    }

    if (MAP_ATT != gs_get_att_src(gs, att)) {
	sprintf(catstr, "no category info");
	return (-1);
    }

    buff = gs_get_att_typbuff(gs, att, 0);

    vrow = Y2VROW(gs, pt[Y]);
    vcol = X2VCOL(gs, pt[X]);
    drow = VROW2DROW(gs, vrow);
    dcol = VCOL2DCOL(gs, vcol);

    offset = DRC2OFF(gs, drow, dcol);


    if (GET_MAPATT(buff, offset, ftmp)) {
	return
	    (Gs_get_cat_label
	     (gsds_get_name(gs->att[att].hdata), drow, dcol, catstr));
    }

    sprintf(catstr, "no data");

    return (1);
}

/*!
  \brief Get surface normal at x,y (real coordinates)

  Usually call after GS_get_selected_point_on_surface()

  \param id surface id
  \param x,y real coordinates
  \param[out] nv surface normal
  
  \return -1 if point outside of window or masked
  \return 1 on success
*/
int GS_get_norm_at_xy(int id, float x, float y, float *nv)
{
    int offset, drow, dcol, vrow, vcol;
    float pt[3];
    geosurf *gs;

    gs = gs_get_surf(id);

    if (NULL == gs) {
	return (-1);
    }

    if (gs->norm_needupdate) {
	gs_calc_normals(gs);
    }

    pt[X] = x;
    pt[Y] = y;

    gsd_real2surf(gs, pt);
    if (gs_point_is_masked(gs, pt)) {
	return (-1);
    }

    if (!in_vregion(gs, pt)) {
	return (-1);
    }

    vrow = Y2VROW(gs, pt[Y]);
    vcol = X2VCOL(gs, pt[X]);
    drow = VROW2DROW(gs, vrow);
    dcol = VCOL2DCOL(gs, vcol);

    offset = DRC2OFF(gs, drow, dcol);

    if (gs->norms) {
	FNORM(gs->norms[offset], nv);
    }
    else {
	/* otherwise must be a constant */
	nv[0] = 0.0;
	nv[1] = 0.0;
	nv[2] = 1.0;
    }

    return (1);
}

/*!
  \brief Get RGB color  at given point

  Colors are translated to rgb and returned as Rxxx Gxxx Bxxx
  Usually call after GS_get_selected_point_on_surface()
  
  Prints "NULL" or the value (i.e., "921.5") to valstr

  Usually call after GS_get_selected_point_on_surface()

  \param id surface id
  \param att
  \param[out] valstr value string (allocated, dim?)
  \param x,y real coordinates

  \return -1 if point outside of window or masked
  \return 1 on success
*/
int GS_get_val_at_xy(int id, int att, char *valstr, float x, float y)
{
    int offset, drow, dcol, vrow, vcol;
    float ftmp, pt[3];
    typbuff *buff;
    geosurf *gs;

    valstr = '\0';
    gs = gs_get_surf(id);

    if (NULL == gs) {
	return (-1);
    }

    pt[X] = x;
    pt[Y] = y;

    gsd_real2surf(gs, pt);

    if (gs_point_is_masked(gs, pt)) {
	return (-1);
    }

    if (!in_vregion(gs, pt)) {
	return (-1);
    }

    if (CONST_ATT == gs_get_att_src(gs, att)) {
	if (att == ATT_COLOR) {
	    int r, g, b, i;
	    i = gs->att[att].constant;
	    sprintf(valstr, "R%d G%d B%d",
		    INT_TO_RED(i, r), INT_TO_GRN(i, g), INT_TO_BLU(i, b));
	}
	else {
	    sprintf(valstr, "%f", gs->att[att].constant);
	}

	return (1);
    }
    else if (MAP_ATT != gs_get_att_src(gs, att)) {
	return (-1);
    }

    buff = gs_get_att_typbuff(gs, att, 0);

    vrow = Y2VROW(gs, pt[Y]);
    vcol = X2VCOL(gs, pt[X]);
    drow = VROW2DROW(gs, vrow);
    dcol = VCOL2DCOL(gs, vcol);

    offset = DRC2OFF(gs, drow, dcol);

    if (GET_MAPATT(buff, offset, ftmp)) {
	if (att == ATT_COLOR) {
	    int r, g, b, i;
	    i = gs_mapcolor(gs_get_att_typbuff(gs, ATT_COLOR, 0),
			    &(gs->att[ATT_COLOR]), offset);
	    sprintf(valstr, "R%d G%d B%d",
		    INT_TO_RED(i, r), INT_TO_GRN(i, g), INT_TO_BLU(i, b));
	}
	else {
	    sprintf(valstr, "%f", ftmp);
	}

	return (1);
    }

    sprintf(valstr, "NULL");

    return (1);
}

#endif

/*!
  \brief Unset attribute

  \param id surface id
  \param att attribute id

  \return ?
*/
int GS_unset_att(int id, int att)
{
    geosurf *gs;

    gs = gs_get_surf(id);
    gs->mask_needupdate = 1;

    return (gs_set_att_src(gs, att, NOTSET_ATT));
}

/*!
  \brief Set attribute constant

  \param id surface id
  \param att attribute id
  \param constant value

  \return ?
*/
int GS_set_att_const(int id, int att, float constant)
{
    geosurf *gs;
    int ret;

    gs = gs_get_surf(id);
    ret = (gs_set_att_const(gs, att, constant));

    Gs_update_attrange(gs, att);

    return (ret);
}

/*!
  \brief Set mask mode

  Mask attribute special: constant is set to indicate invert or no

  \param id surface id
  \param mode id

  \return mode id
  \return -1 on error (invalid surface id)
*/
int GS_set_maskmode(int id, int mode)
{
    geosurf *gs;

    gs = gs_get_surf(id);

    if (gs) {
	gs->att[ATT_MASK].constant = mode;
	gs->mask_needupdate = 1;

	return (mode);
    }

    return (-1);
}

/*!
  \brief Get mask mode

  \param id surface id
  \param[out] mode id

  \return 1 on success
  \return -1 on error (invalid surface id)
*/
int GS_get_maskmode(int id, int *mode)
{
    geosurf *gs;

    gs = gs_get_surf(id);

    if (gs) {
	*mode = gs->att[ATT_MASK].constant;

	return (1);
    }

    return (-1);
}

/*!
  \brief Set client data

  \param id surface id
  \param clientd pointer to client data struct

  \return 1 on success
  \return -1 on error (invalid surface id)
*/
int GS_Set_ClientData(int id, void *clientd)
{
    geosurf *gs;

    gs = gs_get_surf(id);
    if (gs) {
	gs->clientdata = clientd;

	return (1);
    }

    return (-1);
}

/*!
  \brief Get client data

  \param id surface id

  \return pointer to client data
  \return NULL on error
*/
void *GS_Get_ClientData(int id)
{
    geosurf *gs;

    gs = gs_get_surf(id);
    if (gs) {
	return (gs->clientdata);
    }

    return (NULL);
}

/*!
  \brief Get number of surfaces

  \return number of surfaces
*/
int GS_num_surfs(void)
{
    return (gs_num_surfaces());
}

/*!
  \brief Get surface list
  
  Must be freed when not neeed!

  \param[out] numsurf number of available surfaces

  \return pointer to surface array
  \return NULL on error
*/
int *GS_get_surf_list(int *numsurfs)
{
    int i, *ret;

    *numsurfs = Next_surf;

    if (Next_surf) {
	ret = (int *) G_malloc(Next_surf * sizeof(int));

	for (i = 0; i < Next_surf; i++) {
	    ret[i] = Surf_ID[i];
	}

	return (ret);
    }

    return (NULL);
}

/*!
  \brief Delete surface

  \param id surface id

  \return 1 on success
  \return -1 on error
*/
int GS_delete_surface(int id)
{
    int i, j, found = 0;

    G_debug(3, "GS_delete_surface");

    if (GS_surf_exists(id)) {
	gs_delete_surf(id);

	for (i = 0; i < Next_surf && !found; i++) {
	    if (Surf_ID[i] == id) {
		found = 1;

		for (j = i; j < Next_surf; j++) {
		    Surf_ID[j] = Surf_ID[j + 1];
		}
	    }
	}

	gv_update_drapesurfs();

	if (found) {
	    --Next_surf;
	    return (1);
	}
    }

    return (-1);
}


/*!
  \brief Load raster map as attribute

  \param id surface id
  \param filename filename
  \param att attribute descriptor

  \return -1 on error (invalid surface id)
  \return ?
*/
int GS_load_att_map(int id, char *filename, int att)
{
    geosurf *gs;
    unsigned int changed;
    unsigned int atty;
    char *mapset;
    struct Cell_head rast_head;
    int reuse = 0, begin, hdata, ret, neg = 0, has_null = 0;
    typbuff *tbuff;

    G_debug (3, "GS_load_att_map(): map=%s", filename);

    gs = gs_get_surf(id);

    if (NULL == gs) {
	return (-1);
    }

    gs->mask_needupdate = (ATT_MASK == att || ATT_TOPO == att ||
			   (gs->nz_topo && ATT_TOPO == att) ||
			   (gs->nz_color && ATT_COLOR == att));

    gs_set_att_src(gs, att, MAP_ATT);

    /* Check against maps already loaded in memory   */
    /* if to be color attribute:
       - if packed color for another surface, OK to reuse
       - if unchanged, ok to reuse IF it's of type char (will have lookup)
     */
    begin = hdata = 1;

    /* Get MAPSET to ensure names are fully qualified */
    mapset = G_find_cell2(filename, "");
    if (mapset == NULL) {
	/* Check for valid filename */
	G_warning ("Raster map <%s> not found",
		   filename);
	return -1;
    }
    filename = G_fully_qualified_name(filename, mapset);

    /* Check to see if map is in Region */
    G_get_cellhd(filename, mapset, &rast_head);
    if (rast_head.north <= wind.south ||
		    rast_head.south >= wind.north ||
		    rast_head.east <= wind.west ||
		    rast_head.west >= wind.east)
    {
	G_fatal_error(_("Raster map <%s> is outside of current region. Load failed."),
		      filename);
    }

    while (!reuse && (0 < hdata)) {
	changed = CF_COLOR_PACKED;
	atty = ATTY_FLOAT | ATTY_CHAR | ATTY_INT | ATTY_SHORT | ATTY_MASK;

	if (0 < (hdata = gsds_findh(filename, &changed, &atty, begin))) {
	    
	    G_debug(3, "GS_load_att_map(): %s already has data handle %d.CF=%x",
		    filename, hdata, changed);

	    /* handle found */
	    if (ATT_COLOR == att) {
		if ((changed == CF_COLOR_PACKED) ||
		    (!changed && atty == ATTY_CHAR)) {
		    reuse = 1;
		}
	    }
	    else if (atty == ATTY_MASK && att != ATT_MASK) {
		reuse = 0;
		/* should also free mask data & share new - but need backward
		   reference? */
	    }
	    else if (!changed) {
		reuse = 1;
	    }
	}

	begin = 0;
    }

    if (reuse) {
	gs->att[att].hdata = hdata;
	gs_set_att_type(gs, att, atty);	/* ?? */

	/* free lookup  & set to NULL! */
	if (atty == ATTY_INT) {
	    if (gs->att[att].lookup) {
		free(gs->att[att].lookup);
		gs->att[att].lookup = NULL;
	    }
	}
	/* TODO: FIX THIS stuff with lookup sharing! */

	G_debug(3, "GS_load_att_map(): %s is being reused. hdata=%d",
		filename, hdata);
    }
    else {
	G_debug(3, "GS_load_att_map(): %s not loaded in correct form - loading now",
		filename);
	
	/* not loaded - need to get new dataset handle */
	gs->att[att].hdata = gsds_newh(filename);

	tbuff = gs_get_att_typbuff(gs, att, 1);

	/* TODO: Provide mechanism for loading certain attributes at
	   specified sizes, allow to scale or cap, or scale non-zero */
	if (ATT_MASK == att) {
	    atty = ATTY_MASK;
	}
	else {
	    atty = Gs_numtype(filename, &neg);
	    filename = G_fully_qualified_name(filename, mapset);
	}

#ifdef MAYBE_LATER
	if (att == ATT_COLOR && atty == ATTY_SHORT) {
	    atty = (neg ? ATTY_INT : ATTY_SHORT);
	}
#endif

	if (att == ATT_COLOR && atty == ATTY_SHORT) {
	    atty = ATTY_INT;
	}

	if (0 > gs_malloc_att_buff(gs, att, ATTY_NULL)) {
	    G_fatal_error (_("GS_load_att_map(): Out of memory. Unable to load map"));
	}

	switch (atty) {
	case ATTY_MASK:
	    if (0 > gs_malloc_att_buff(gs, att, ATTY_MASK)) {
		G_fatal_error (_("GS_load_att_map(): Out of memory. Unable to load map"));
	    }

	    ret = Gs_loadmap_as_bitmap(&wind, filename, tbuff->bm);
	    filename = G_fully_qualified_name(filename, mapset);

	    break;
	case ATTY_CHAR:
	    if (0 > gs_malloc_att_buff(gs, att, ATTY_CHAR)) {
		G_fatal_error (_("GS_load_att_map(): Out of memory. Unable to load map"));
	    }

	    ret = Gs_loadmap_as_char(&wind, filename, tbuff->cb,
				     tbuff->nm, &has_null);
	    filename = G_fully_qualified_name(filename, mapset);

	    break;
	case ATTY_SHORT:
	    if (0 > gs_malloc_att_buff(gs, att, ATTY_SHORT)) {
		G_fatal_error (_("GS_load_att_map(): Out of memory. Unable to load map"));
	    }

	    ret = Gs_loadmap_as_short(&wind, filename, tbuff->sb,
				      tbuff->nm, &has_null);
	    filename = G_fully_qualified_name(filename, mapset);

	    break;
	case ATTY_FLOAT:
	    if (0 > gs_malloc_att_buff(gs, att, ATTY_FLOAT)) {
		G_fatal_error (_("GS_load_att_map(): Out of memory. Unable to load map"));
	    }

	    ret = Gs_loadmap_as_float(&wind, filename, tbuff->fb,
				      tbuff->nm, &has_null);
	    filename = G_fully_qualified_name(filename, mapset);

	    break;
	case ATTY_INT:
	default:
	    if (0 > gs_malloc_att_buff(gs, att, ATTY_INT)) {
		G_fatal_error (_("GS_load_att_map(): Out of memory. Unable to load map"));
	    }

	    ret = Gs_loadmap_as_int(&wind, filename, tbuff->ib,
				    tbuff->nm, &has_null);
	    filename = G_fully_qualified_name(filename, mapset);

	    break;

	}			/* Done with switch */

	if (ret == -1) {
	    gsds_free_data_buff(gs->att[att].hdata, ATTY_NULL);
	    return (-1);
	}

	G_debug(4, "  has_null=%d", has_null);

	if (!has_null) {
	    gsds_free_data_buff(gs->att[att].hdata, ATTY_NULL);
	}
	else {
	    gs_update_curmask(gs);
	}

    }				/* end if not reuse */

    if (ATT_COLOR == att) {
#ifdef MAYBE_LATER
	if (ATTY_INT == atty) {
	    Gs_pack_colors(filename, tbuff->ib, gs->rows, gs->cols);
	    gsds_set_changed(gs->att[att].hdata, CF_COLOR_PACKED);
	    gs->att[att].lookup = NULL;
	    filename = G_fully_qualified_name(filename, mapset);
	}
	else {
	    gs_malloc_lookup(gs, att);
	    Gs_build_lookup(filename, gs->att[att].lookup);
	    filename = G_fully_qualified_name(filename, mapset);
	}
#else

	if (ATTY_CHAR == atty) {
	    if (!gs->att[att].lookup) {
		/* might already exist if reusing */
		gs_malloc_lookup(gs, att);
		Gs_build_256lookup(filename, gs->att[att].lookup);
		filename = G_fully_qualified_name(filename, mapset);
	    }
	}
	else if (ATTY_FLOAT == atty) {
	    if (!reuse) {
		if (0 > gs_malloc_att_buff(gs, att, ATTY_INT)) {
		    G_fatal_error (_("GS_load_att_map(): Out of memory. Unable to load map"));
		}

		Gs_pack_colors_float(filename, tbuff->fb, tbuff->ib,
				     gs->rows, gs->cols);
		gsds_set_changed(gs->att[att].hdata, CF_COLOR_PACKED);
		gsds_free_data_buff(gs->att[att].hdata, ATTY_FLOAT);
		gs->att[att].lookup = NULL;
		filename = G_fully_qualified_name(filename, mapset);
	    }
	}
	else {
	    if (!reuse) {
		Gs_pack_colors(filename, tbuff->ib, gs->rows, gs->cols);
		gsds_set_changed(gs->att[att].hdata, CF_COLOR_PACKED);
		gs->att[att].lookup = NULL;
		filename = G_fully_qualified_name(filename, mapset);
	    }
	}
#endif
    }

    if (ATT_TOPO == att) {
	gs_init_normbuff(gs);
	/* S_DIFF: should also check here to see if this surface is a
	   reference surface for scaled differences, if so update references
	   to it */
    }

    if (ret < 0) {
	G_warning (_("Loading failed"));
    }

    if (-1 == Gs_update_attrange(gs, att)) {
	G_warning(_("Error finding range"));
    }

    return (ret);
}

/*!
  \brief Draw surface

  \param id surface id
*/
void GS_draw_surf(int id)
{
    geosurf *gs;

    G_debug(3, "GS_draw_surf(): id=%d", id);

    gs = gs_get_surf(id);
    if (gs) {
	gsd_shademodel(gs->draw_mode & DM_GOURAUD);

	if (gs->draw_mode & DM_POLY) {
	    gsd_surf(gs);
	}

	if (gs->draw_mode & DM_WIRE) {
	    gsd_wire_surf(gs);
	}

	/* TODO: write wire/poly draw routines */
	if (gs->draw_mode & DM_WIRE_POLY) {
	    gsd_surf(gs);
	    gsd_wire_surf(gs);
	}
    }

    return;
}

/*!
  \brief Draw surface wire

  Overrides draw_mode for fast display

  \param id surface id
*/
void GS_draw_wire(int id)
{
    geosurf *gs;

    G_debug(3, "GS_draw_wire(): id=%d", id);

    gs = gs_get_surf(id);

    if (gs) {
	gsd_wire_surf(gs);
    }

    return;
}

/*!
  \brief Draw all wires

  Overrides draw_mode for fast display
*/
void GS_alldraw_wire(void)
{
    geosurf *gs;
    int i;

    for (i = 0; i < Next_surf; i++) {
	if ((gs = gs_get_surf(Surf_ID[i]))) {
	    gsd_wire_surf(gs);
	}
    }

    return;
}

/*!
  \brief Draw all surfaces
*/
void GS_alldraw_surf(void)
{
    int i;

    for (i = 0; i < Next_surf; i++) {
	GS_draw_surf(Surf_ID[i]);
    }

    return;
}

/*!
  \brief Set Z exag for surface

  \param id surface id
  \param exag z-exag value
*/
void GS_set_exag(int id, float exag)
{
    geosurf *gs;

    G_debug(3, "GS_set_exag");

    gs = gs_get_surf(id);

    if (gs) {
	if (gs->z_exag != exag) {
	    gs->norm_needupdate = 1;
	}

	gs->z_exag = exag;
    }

    return;
}

/*!
  \brief Set global z-exag value
  
  \param exag exag value to be set up
*/
void GS_set_global_exag(float exag)
{

    G_debug(3, "GS_set_global_exag");

    Gv.vert_exag = exag;
    /* GL_NORMALIZE */
    /* Only need to update norms gs_norms.c
     * if exag is used in norm equation which
     * it is not! If GL_NORMALIZE is disabled
     * will need to include.
    gs_setall_norm_needupdate();
    */

    return;
}

/*!
  \brief Get global z-exag value

  \return value
*/
float GS_global_exag(void)
{
    G_debug(3, "GS_global_exag(): %g", Gv.vert_exag);

    return (Gv.vert_exag);
}

/*!
  \brief Set wire color

  \todo error-handling

  \param id surface id
  \param colr color value
*/
void GS_set_wire_color(int id, int colr)
{
    geosurf *gs;

    G_debug(3, "GS_set_wire_color");

    gs = gs_get_surf(id);

    if (gs) {
	gs->wire_color = colr;
    }

    return;
}

/*!
  \brief Get wire color

  \param id surface id
  \param[out] colr color value

  \return 1 on success
  \return -1 on error
*/
int GS_get_wire_color(int id, int *colr)
{
    geosurf *gs;

    gs = gs_get_surf(id);

    if (gs) {
	*colr = gs->wire_color;

	return (1);
    }

    return (-1);
}

/*!
  \brief Set all draw-modes

  \param mode mode id

  \return 0 on success
  \return -1 on error
*/
int GS_setall_drawmode(int mode)
{
    int i;

    for (i = 0; i < Next_surf; i++) {
	if (0 != GS_set_drawmode(Surf_ID[i], mode)) {
	    return (-1);
	}
    }

    return (0);
}

/*!
  \brief Set draw mode

  \param id surface id
  \param mode mode type(s)

  \return 0 on success
  \return -1 on error (invalid surface id)
*/
int GS_set_drawmode(int id, int mode)
{
    geosurf *gs;

    G_debug(3, "GS_set_drawmode(): id=%d mode=%d",
	    id, mode);

    gs = gs_get_surf(id);

    if (gs) {
	gs->draw_mode = mode;

	return (0);
    }

    return (-1);
}

/*!
  \brief Get draw mode

  \param id surface id
  \param[out] mode mode id

  \return 1 on success
  \return -1 on error (invalid surface id)
*/
int GS_get_drawmode(int id, int *mode)
{
    geosurf *gs;

    gs = gs_get_surf(id);

    if (gs) {
	*mode = gs->draw_mode;

	return (1);
    }

    return (-1);
}

/*!
  \brief Set no-zero ?
  
  \param id surface id
  \param att attribute id
  \param mode mode id
*/
void GS_set_nozero(int id, int att, int mode)
{
    geosurf *gs;

    G_debug(3, "GS_set_nozero");

    gs = gs_get_surf(id);

    if (gs) {
	if (att == ATT_TOPO) {
	    gs->nz_topo = mode;
	    gs->mask_needupdate = 1;
	}

	if (att == ATT_COLOR) {
	    gs->nz_color = mode;
	    gs->mask_needupdate = 1;
	}
    }

    return;
}

/*!
  \brief Get no-zero ?
  
  \param id surface id
  \param att attribute id
  \param[out] mode mode id

  \return -1 on error (invalid surface id)
  \return 1 on success
*/
int GS_get_nozero(int id, int att, int *mode)
{
    geosurf *gs;

    G_debug(3, "GS_set_nozero");

    gs = gs_get_surf(id);

    if (gs) {
	if (att == ATT_TOPO) {
	    *mode = gs->nz_topo;
	}
	else if (att == ATT_COLOR) {
	    *mode = gs->nz_color;
	}
	else {
	    return (-1);
	}

	return (1);
    }

    return (-1);
}

/*!
  \brief Set all draw resolutions

  \param xres,yres x/y resolution value
  \param xwire,ywire x/y wire value

  \return 0 on success
  \return -1 on error
*/
int GS_setall_drawres(int xres, int yres, int xwire, int ywire)
{
    int i;

    for (i = 0; i < Next_surf; i++) {
	if (0 != GS_set_drawres(Surf_ID[i], xres, yres, xwire, ywire)) {
	    return (-1);
	}
    }

    return (0);
}

/*!
  \brief Set draw resolution for surface

  \param id surface id
  \param xres,yres x/y resolution value
  \param xwire,ywire x/y wire value

  \return -1 on error
  \return 0 on success
*/
int GS_set_drawres(int id, int xres, int yres, int xwire, int ywire)
{
    geosurf *gs;

    G_debug(3, "GS_set_drawres() id=%d xyres=%d/%d xywire=%d/%d",
	    id, xres, yres, xwire, ywire);

    if (xres < 1 || yres < 1 || xwire < 1 || ywire < 1) {
	return (-1);
    }

    gs = gs_get_surf(id);

    if (gs) {
	if (gs->x_mod != xres || gs->y_mod != yres) {
	    gs->norm_needupdate = 1;
	}

	gs->x_mod = xres;
	gs->y_mod = yres;
	gs->x_modw = xwire;
	gs->y_modw = ywire;
    }

    return (0);
}

/*!
  \brief Get draw resolution of surface

  \param id surface id
  \param[out] xres,yres x/y resolution value
  \param[out] xwire,ywire x/y wire value
*/
void GS_get_drawres(int id, int *xres, int *yres, int *xwire, int *ywire)
{
    geosurf *gs;

    G_debug(3, "GS_get_drawres");

    gs = gs_get_surf(id);

    if (gs) {
	*xres = gs->x_mod;
	*yres = gs->y_mod;
	*xwire = gs->x_modw;
	*ywire = gs->y_modw;
    }

    return;
}

/*!
  \brief Get dimension of surface

  \param id surface id
  \param[out] rows,cols number of rows/cols
*/
void GS_get_dims(int id, int *rows, int *cols)
{
    geosurf *gs;

    gs = gs_get_surf(id);

    if (gs) {
	*rows = gs->rows;
	*cols = gs->cols;
    }

    return;
}

/*!
  \brief Get exag-value guess

  Use no_zero range because if zero IS data, then range won't be that
  much off (it's just a GUESS, after all), but if zero is NO data, could
  drastically affect guess

  \param id surface id
  \param[out] exag exag value

  \return 1 on success
  \return -1 on error
*/
int GS_get_exag_guess(int id, float *exag)
{
    geosurf *gs;
    float guess;

    gs = gs_get_surf(id);
    guess = 1.0;

    /* if gs is type const return guess = 1.0 */
    if (CONST_ATT == gs_get_att_src(gs, ATT_TOPO))
    {
	return(1);
    }

    if (gs) {
	if (gs->zrange_nz == 0.0) {
	    *exag = 0.0;

	    return (1);
	}

	G_debug (3, "GS_get_exag_guess(): %f %f",
		 gs->zrange_nz, Longdim);

	while (gs->zrange_nz * guess / Longdim >= .25) {
	    guess *= .1;

	    G_debug (3, "GS_get_exag_guess(): %f",
		     guess);
	}

	while (gs->zrange_nz * guess / Longdim < .025) {
	    guess *= 10.;

	    G_debug (3, "GS_get_exag_guess(): %f",
		     guess);
	}

	*exag = guess;

	return (1);
    }

    return (-1);
}

/*!
  \brief Get Z extents for all loaded surfaces

  Treating zeros as "no data"

  \param[out] min min value
  \param[out] max max value
*/
void GS_get_zrange_nz(float *min, float *max)
{
    int i, first = 1;
    geosurf *gs;

    for (i = 0; i < Next_surf; i++) {
	if ((gs = gs_get_surf(Surf_ID[i]))) {
	    if (first) {
		first = 0;
		*min = gs->zmin_nz;
		*max = gs->zmax_nz;
	    }

	    if (gs->zmin_nz < *min) {
		*min = gs->zmin_nz;
	    }

	    if (gs->zmax_nz > *max) {
		*max = gs->zmax_nz;
	    }
	}
    }

    G_debug (3, "GS_get_zrange_nz(): min=%g max=%g",
             *min, *max);

    return;
}

/*!
  \brief Set trans ?

  \param id surface id
  \param xtrans,ytrans,ztrans real trans coordinates
*/
void GS_set_trans(int id, float xtrans, float ytrans, float ztrans)
{
    geosurf *gs;

    G_debug(3, "GS_set_trans");

    gs = gs_get_surf(id);

    if (gs) {
	gs->x_trans = xtrans;
	gs->y_trans = ytrans;
	gs->z_trans = ztrans;
    }

    return;
}

/*!
  \brief Get trans ?

  \param id surface id
  \param[out] xtrans,ytrans,ztrans real trans coordinates
*/
void GS_get_trans(int id, float *xtrans, float *ytrans, float *ztrans)
{
    geosurf *gs;

    G_debug(3, "GS_get_trans");

    gs = gs_get_surf(id);

    if (gs) {
	*xtrans = gs->x_trans;
	*ytrans = gs->y_trans;
	*ztrans = gs->z_trans;
    }

    return;
}


/*!
  \brief Get default draw color

  \return color value
*/
unsigned int GS_default_draw_color(void)
{

    G_debug(3, "GS_default_draw_color");

    return ((unsigned int) Gd.bgcol);
}

/*!
  \brief Get background color

  \return color value
*/
unsigned int GS_background_color(void)
{
    return ((unsigned int) Gd.bgcol);
}

/*!
  \brief Sets which buffer to draw to

  \param where GSD_BOTH, GSD_FRONT, GSD_BACK
*/
void GS_set_draw(int where)
{
    Buffermode = where;

    switch (where) {
    case GSD_BOTH:
	gsd_frontbuffer(1);
	gsd_backbuffer(1);

	break;
    case GSD_FRONT:
	gsd_frontbuffer(1);
	gsd_backbuffer(0);

	break;
    case GSD_BACK:
    default:
	gsd_frontbuffer(0);
	gsd_backbuffer(1);

	break;
    }

    return;
}

/*
  \brief Ready to draw
*/
void GS_ready_draw(void)
{

    G_debug(3, "GS_ready_draw");

    gsd_set_view(&Gv, &Gd);

    return;
}

/*!
  \brief Draw done, swap buffers
*/
void GS_done_draw(void)
{

    G_debug(3, "GS_done_draw");

    if (GSD_BACK == Buffermode) {
	gsd_swapbuffers();
    }

    gsd_flush();

    return;
}

/*!
  \brief Set focus

  \param realto real coordinates to
*/
void GS_set_focus(float *realto)
{

    G_debug(3, "GS_set_focus(): %f,%f,%f",
	    realto[0], realto[1], realto[2]);

    Gv.infocus = 1;
    GS_v3eq(Gv.real_to, realto);

    gsd_set_view(&Gv, &Gd);

    return;
}

/*!
  \brief Set real focus

  \param realto real coordinates to
*/
void GS_set_focus_real(float *realto)
{

    G_get_set_window(&wind);
    realto[X] = realto[X] - wind.west - (wind.ew_res/2.);
    realto[Y] = realto[Y] - wind.south - (wind.ns_res/2.);

    Gv.infocus = 1;
    GS_v3eq(Gv.real_to, realto);

    gsd_set_view(&Gv, &Gd);

    return;
}


/*!
  \brief Get focus

  OK to call with NULL argument if just want to check state

  \param realto real coordinates to

  \return ?
*/
int GS_get_focus(float *realto)
{

    G_debug(3, "GS_get_focus");

    if (Gv.infocus) {
	if (realto) {
	    GS_v3eq(realto, Gv.real_to);
	}
    }

    return (Gv.infocus);
}

/*!
  \brief Set focus to map center

  \param id surface id
*/
void GS_set_focus_center_map(int id)
{
    float center[3];
    geosurf *gs;

    G_debug(3, "GS_set_focus_center_map");

    gs = gs_get_surf(id);

    if (gs) {
	center[X] = (gs->xmax - gs->xmin) / 2.;
	center[Y] = (gs->ymax - gs->ymin) / 2.;
	center[Z] = (gs->zmax_nz + gs->zmin_nz) / 2.;

	/* not yet working
	   buff = gs_get_att_typbuff(gs, ATT_TOPO, 0);
	   offset = gs->rows*gs->cols/2 + gs->cols/2;
	   if (buff)
	   {
	   if (GET_MAPATT(buff, offset, tmp))
	   {
	   center[Z] = tmp;
	   }
	   }
	 */

	GS_set_focus(center);
    }
}

/*!
  \brief Move viewpoint 

  \param pt 'from' model coordinates
*/
void GS_moveto(float *pt)
{
    float ft[3];

    G_debug(3, "GS_moveto(): %f,%f,%f", pt[0], pt[1], pt[2]);

    if (Gv.infocus) {
	GS_v3eq(Gv.from_to[FROM], pt);
	/*
	GS_v3eq(Gv.from_to[TO], Gv.real_to);
	*/
	GS_v3normalize(Gv.from_to[FROM], Gv.from_to[TO]);
	/* update inclination, look_dir if we're keeping these */
    }
    else {
	GS_v3eq(ft, Gv.from_to[TO]);
	GS_v3sub(ft, Gv.from_to[FROM]);
	GS_v3eq(Gv.from_to[FROM], pt);
	GS_v3eq(Gv.from_to[TO], pt);
	GS_v3add(Gv.from_to[TO], ft);
    }

    return;
}

/*!
  \brief Move position to (real)

  \param pt point real coordinates
*/
void GS_moveto_real(float *pt)
{
    gsd_real2model(pt);
    GS_moveto(pt);

    return;
}

/*!
  \brief Get z-extent for a single surface

  \param id surface id
  \param[out] min min z-value
  \param[out] max max z-value
  \param[out] mid middle z-value

  \return -1 on error (invalid surface id)
  \return ?
*/
int GS_get_zextents(int id, float *min, float *max, float *mid)
{
    geosurf *gs;

    if (NULL == (gs = gs_get_surf(id))) {
	return (-1);
    }

    G_debug(3, "GS_get_zextents(): id=%d", id);

    return (gs_get_zextents(gs, min, max, mid));
}

/*!
  \brief Get z-extent for all loaded surfaces

  \param[out] min min z-value
  \param[out] max max z-value
  \param doexag use z-exaggeration

  \return 1 on success
  \return -1 on error
*/
int GS_get_zrange(float *min, float *max, int doexag)
{
    int ret_surf, ret_vol;
    float surf_min, surf_max;
    float vol_min, vol_max;
    
    ret_surf = gs_get_zrange(&surf_min, &surf_max);
    ret_vol = gvl_get_zrange(&vol_min, &vol_max);
    
    if (ret_surf > 0 && ret_vol > 0) {
	*min = (surf_min < vol_min) ? surf_min : vol_min;
	*max = (surf_max < vol_max) ? surf_max : vol_max;
    } else if (ret_surf > 0) {
	*min = surf_min;
	*max = surf_max;
    } else if (ret_vol > 0) {
	*min = vol_min;
	*max = vol_max;
    }
    
    if (doexag) {
	*min *= Gv.vert_exag;
	*max *= Gv.vert_exag;
    }
    
    G_debug(3, "GS_get_zrange(): min=%g max=%g",
	    *min, *max);
    return ((ret_surf > 0 || ret_vol > 0) ? (1) :(-1));
}

/*!
  \brief Get viewpoint 'from' position

  \param[out] fr from model coordinates
*/
void GS_get_from(float *fr)
{
    GS_v3eq(fr, Gv.from_to[FROM]);

    G_debug(3, "GS_get_from(): %f,%f,%f", fr[0], fr[1], fr[2]);

    return;
}

/*!
  \brief Get viewpoint 'from' real coordinates

  \param[out] fr 'from' real coordinates
*/
void GS_get_from_real(float *fr)
{
    GS_v3eq(fr, Gv.from_to[FROM]);
    gsd_model2real(fr);

    return;
}

/*!
  \brief Get 'to' real coordinates

  \param[out] to 'to' real coordinates
*/
void GS_get_to_real(float *to)
{
float realto[3];

    G_get_set_window(&wind);
    GS_get_focus(realto);
    to[X] = realto[X] + wind.west + (wind.ew_res/2.);
    to[Y] = realto[Y] + wind.south + (wind.ns_res/2.);
    to[Z] = realto[Z];

    return;
}


/*!
  \brief Get zoom setup
  
  \param[out] a,b,c,d current viewport settings
  \param[out] maxx,maxy max viewport size
*/
void GS_zoom_setup(int *a, int *b, int *c, int *d, int *maxx, int *maxy)
{
    GLint tmp[4];
    GLint num[2];
    gsd_getViewport(tmp, num);
    *a = tmp[0];
    *b = tmp[1];
    *c = tmp[2];
    *d = tmp[3];
    *maxx = num[0];
    *maxy = num[1];

    return;
}

/*!
  \brief Get 'to' model coordinates

  \todo need set_to? - just use viewdir?

  \param[out] to 'to' model coordinates
*/
void GS_get_to(float *to)
{
    G_debug(3, "GS_get_to");

    GS_v3eq(to, Gv.from_to[TO]);

    return;
}

/*!
  \brief Get viewdir

  \param[out] dir viewdir value
*/
void GS_get_viewdir(float *dir)
{
    GS_v3dir(Gv.from_to[FROM], Gv.from_to[TO], dir);

    return;
}

/*!
  \brief Set viewdir
  
  Automatically turns off focus

  \param dir viewdir value
*/
void GS_set_viewdir(float *dir)
{
    float tmp[3];

    GS_v3eq(tmp, dir);
    GS_v3norm(tmp);
    GS_v3eq(Gv.from_to[TO], Gv.from_to[FROM]);
    GS_v3add(Gv.from_to[TO], tmp);

    GS_set_nofocus();
    gsd_set_view(&Gv, &Gd);

    return;
}

/*!
  \brief Set field of view

  \param fov fov value
*/
void GS_set_fov(int fov)
{
    Gv.fov = fov;

    return;
}

/*!
  \brief Get fied of view

  \return field of view, in 10ths of degrees
*/
int GS_get_fov(void)
{
    return (Gv.fov);
}

/*!
  \brief Get twist value

  10ths of degrees off twelve o'clock
*/
int GS_get_twist(void)
{
    return (Gv.twist);
}

/*!
  \brief Set viewpoint twist value

  10ths of degrees off twelve o'clock
  
  \params t tenths of degrees clockwise from 12:00.
*/
void GS_set_twist(int t)
{
    Gv.twist = t;

    return;
}

/*!
  \brief Unset focus
*/
void GS_set_nofocus(void)
{
    G_debug(3, "GS_set_nofocus");

    Gv.infocus = 0;

    return;
}

/*!
  \brief Set focus

  Make sure that the center of view is set
*/
void GS_set_infocus(void)
{
    G_debug(3, "GS_set_infocus");

    Gv.infocus = 1;

    return;
}

/*!
  \brief Set viewport

  \param left,right,bottom,top viewport extent values
*/
void GS_set_viewport(int left, int right, int bottom, int top)
{
    G_debug(3, "GS_set_viewport(): left=%d, right=%d, "
	    "bottom=%d, top=%d", left, right, bottom, top);

    gsd_viewport(left, right, bottom, top);

    return;
}

/*!
  \brief Send screen coords sx and sy, lib traces through surfaces; sets
  new center to point of nearest intersection.

  If no intersection, uses line of sight with length of current view
  ray (eye to center) to set new center.

  Reset center of view to screen coordinates sx, sy.

  \param sx,sy screen coordinates

  \return 1 on success
  \return 0 on error (invalid surface id)
*/
int GS_look_here(int sx, int sy)
{
    float x, y, z, len, los[2][3];
    Point3 realto, dir;
    int id;
    geosurf *gs;

    if (GS_get_selected_point_on_surface(sx, sy, &id, &x, &y, &z)) {
	gs = gs_get_surf(id);
	if (gs) {
	    realto[X] = x - gs->ox + gs->x_trans;
	    realto[Y] = y - gs->oy + gs->y_trans;
	    realto[Z] = z + gs->z_trans;
	    GS_set_focus(realto);

	    return (1);
	}
    }
    else {
	if (gsd_get_los(los, (short) sx, (short) sy)) {
	    len = GS_distance(Gv.from_to[FROM], Gv.real_to);
	    GS_v3dir(los[FROM], los[TO], dir);
	    GS_v3mult(dir, len);
	    realto[X] = Gv.from_to[FROM][X] + dir[X];
	    realto[Y] = Gv.from_to[FROM][Y] + dir[Y];
	    realto[Z] = Gv.from_to[FROM][Z] + dir[Z];
	    GS_set_focus(realto);

	    return (1);
	}
    }

    return (0);
}

/*!
  \brief Get selected point of surface

  Given screen coordinates sx and sy, find closest intersection of
  view ray with surfaces and return coordinates of intersection in x, y,
  z, and identifier of surface in id.

  \param sx,sy screen coordinates
  \param[out] id surface id
  \param[out] x,y,z point on surface (model coordinates?)

  \returns 0 if no intersections found
  \return number of intersections
*/
int GS_get_selected_point_on_surface(int sx, int sy, int *id, float *x,
				     float *y, float *z)
{
    float los[2][3], find_dist[MAX_SURFS], closest;
    Point3 point, tmp, finds[MAX_SURFS];
    int surfs[MAX_SURFS], i, iclose, numhits = 0;
    geosurf *gs;

    /* returns surface-world coords */
    gsd_get_los(los, (short) sx, (short) sy);

    if (!gs_setlos_enterdata(los)) {
	G_debug (3, "gs_setlos_enterdata(los): returns false");
	return (0);
    }

    for (i = 0; i < Next_surf; i++) {
	G_debug (3, "id=%d", i);

	gs = gs_get_surf(Surf_ID[i]);

	/* los_intersect expects surf-world coords (xy transl, no scaling) */

#if NVIZ_HACK
	if (gs_los_intersect1(Surf_ID[i], los, point)) {
#else
	if (gs_los_intersect(Surf_ID[i], los, point)) {
#endif
	    if (!gs_point_is_masked(gs, point)) {
		GS_v3eq(tmp, point);
		tmp[X] += gs->x_trans;
		tmp[Y] += gs->y_trans;
		tmp[Z] += gs->z_trans;
		find_dist[numhits] = GS_distance(los[FROM], tmp);
		gsd_surf2real(gs, point);
		GS_v3eq(finds[numhits], point);
		surfs[numhits] = Surf_ID[i];
		numhits++;
	    }
	}
    }

    for (i = iclose = 0; i < numhits; i++) {
	closest = find_dist[iclose];

	if (find_dist[i] < closest) {
	    iclose = i;
	}
    }

    if (numhits) {
	*x = finds[iclose][X];
	*y = finds[iclose][Y];
	*z = finds[iclose][Z];
	*id = surfs[iclose];
    }

    G_debug(3, "NumHits %d, next %d", numhits, Next_surf);

    return (numhits);
}

/*!
  \brief Set cplace rotation

  \param num cplace id
  \param dx,dy,dz rotation values
*/
void GS_set_cplane_rot(int num, float dx, float dy, float dz)
{
    gsd_cplane_setrot(num, dx, dy, dz);

    return;
}

/*!
  \brief Set cplace trans
  
  \param num cplace id
  \param dx,dy,dz rotation values
*/
void GS_set_cplane_trans(int num, float dx, float dy, float dz)
{
    gsd_cplane_settrans(num, dx, dy, dz);

    return;
}


/*!
  \brief Draw cplace

  \param num cplace id
*/
void GS_draw_cplane(int num)
{
    geosurf *gsurfs[MAX_SURFS];
    int nsurfs;

    nsurfs = gs_num_surfaces();
    if (2 == nsurfs) {
	/* testing */
	gs_getall_surfaces(gsurfs);
	gsd_draw_cplane_fence(gsurfs[0], gsurfs[1], num);
    }
    else {
	gsd_draw_cplane(num);
    }

    return;
}

/*!
  \brief Draw cplace fence ?

  \param hs1,hs2
  \param num cplane id

  \return 0 on error
  \return 1 on success
*/
int GS_draw_cplane_fence(int hs1, int hs2, int num)
{
    geosurf *gs1, *gs2;

    if (NULL == (gs1 = gs_get_surf(hs1))) {
	return (0);
    }

    if (NULL == (gs2 = gs_get_surf(hs2))) {
	return (0);
    }

    gsd_draw_cplane_fence(gs1, gs2, num);

    return (1);
}

/*!
  \brief Draw all cplace fences ?
*/
void GS_alldraw_cplane_fences(void)
{
    int onstate[MAX_CPLANES], i;

    gsd_get_cplanes_state(onstate);

    for (i = 0; i < MAX_CPLANES; i++) {
	if (onstate[i]) {
	    GS_draw_cplane_fence(Surf_ID[0], Surf_ID[1], i);
	}
    }

    return;
}

/*!
  \brief Set cplace

  \param num cplane id
*/
void GS_set_cplane(int num)
{
    gsd_cplane_on(num);

    return;
}

/*!
  \brief Unset clip place (turn off)

  \param num cplane id
*/
void GS_unset_cplane(int num)
{
    gsd_cplane_off(num);

    return;
}

/*!
  \brief Get axis scale

  \param sx,sy,sz x/y/z scale values
  \param doexag use vertical exaggeration
*/
void GS_get_scale(float *sx, float *sy, float *sz, int doexag)
{
    float zexag;

    zexag = doexag ? Gv.vert_exag : 1.;
    *sx = *sy = Gv.scale;
    *sz = Gv.scale * zexag;

    return;
}

/*!
  \brief Set fence color

  \param mode mode id
*/
void GS_set_fencecolor(int mode)
{
    gsd_setfc(mode);

    return;
}

/*!
  \brief Get fence color

  \return color value
*/
int GS_get_fencecolor(void)
{
    return gsd_getfc();
}

/*!
  \brief Measure distance "as the ball rolls" between two points on surface

  \param hs
  \param x1,y1,x2,y2 distance line nodes
  \param dist ?
  \param use_exag use exag ?

  \return 0 on error or if one or more points is not in region,
  \return distance following terrain
*/
int GS_get_distance_alongsurf(int hs, float x1, float y1, float x2, float y2,
			      float *dist, int use_exag)
{
    geosurf *gs;
    float p1[2], p2[2];

    if (NULL == (gs = gs_get_surf(hs))) {
	return (0);
    }

    p1[X] = x1;
    p1[Y] = y1;
    p2[X] = x2;
    p2[Y] = y2;
    gsd_real2surf(gs, p1);
    gsd_real2surf(gs, p2);

    return (gs_distance_onsurf(gs, p1, p2, dist, use_exag));
}

/*!
  \brief Save 3d view

  \param vname view file name
  \param surfid surface id

  \return ?
*/
int GS_save_3dview(char *vname, int surfid)
{
    return (Gs_save_3dview(vname, &Gv, &Gd, &wind, gs_get_surf(surfid)));
}

/*!
  \brief Load 3d view

  \param vname view file name
  \param surfid surface id

  \return ?
*/
int GS_load_3dview(char *vname, int surfid)
{

    return (Gs_load_3dview(vname, &Gv, &Gd, &wind, gs_get_surf(surfid)));

    /* what to do about lights - I guess, delete all &
       create any that exist in 3dview file */
}

/************************************************************************
* Following routines use Graphics Library
************************************************************************/

/*!
  \brief Init viewpoint

  \todo allow to set center?
*/
void GS_init_view(void)
{
    static int first = 1;

    G_debug(3, "GS_init_view");

    if (first) {
	first = 0;
	glMatrixMode(GL_MODELVIEW);

	/* OGLXXX doublebuffer: use GLX_DOUBLEBUFFER in attriblist */
	/* glxChooseVisual(*dpy, screen, *attriblist); */
	/* OGLXXX
	 * ZMIN not needed -- always 0.
	 * ZMAX not needed -- always 1.
	 * getgdesc other posiblilties:
	 *      glxGetConfig();
	 *      glxGetCurrentContext();
	 *      glxGetCurrentDrawable();
	 * GLint gdtmp;
	 * getgdesc other posiblilties:
	 *      glxGetConfig();
	 *      glxGetCurrentContext();
	 *      glxGetCurrentDrawable();
	 * GLint gdtmp;
	 * glDepthRange params must be scaled to [0, 1]
	 */
	glDepthRange(0.0, 1.0);
	glEnable(GL_DEPTH_TEST);
	glDepthFunc(GL_LEQUAL);
    /* } */

    /* replace these with something meaningful */
    Gv.fov = 450;
    Gv.twist = 0;
    Gv.from_to[FROM][X] = Gv.from_to[FROM][Y] =
	Gv.from_to[FROM][Z] = GS_UNIT_SIZE / 2.;

    Gv.from_to[TO][X] = GS_UNIT_SIZE / 2.;
    Gv.from_to[TO][Y] = GS_UNIT_SIZE / 2.;
    Gv.from_to[TO][Z] = 0.;
    Gv.from_to[TO][W] = Gv.from_to[FROM][W] = 1.;

    Gv.real_to[W] = 1.;
    Gv.vert_exag = 1.;

    GS_v3eq(Gv.real_to, Gv.from_to[TO]);
    GS_v3normalize(Gv.from_to[FROM], Gv.from_to[TO]);

    /*
    Gd.nearclip = 50;
    Gd.farclip = 10000.;
    */
    Gd.nearclip = 10.;
    Gd.farclip = 10000.;
    Gd.aspect = (float) GS_get_aspect();

    GS_set_focus(Gv.real_to);
 }

    return;
}

/*!
  \brief Clear view

  \param col color value
*/
void GS_clear(int col)
{
    G_debug(3, "GS_clear");

    col = col | 0xFF000000;

    /* OGLXXX
     * change glClearDepth parameter to be in [0, 1]
     * ZMAX not needed -- always 1.
     * getgdesc other posiblilties:
     *      glxGetConfig();
     *      glxGetCurrentContext();
     *      glxGetCurrentDrawable();
     * GLint gdtmp;
     */
    glClearDepth(1.0);
    glClearColor(((float) ((col) & 0xff)) / 255.,
		 (float) ((col) >> 8 & 0xff) / 255.,
		 (float) ((col) >> 16 & 0xff) / 255.,
		 (float) ((col) >> 24 & 0xff) / 255.);
    glClear(GL_DEPTH_BUFFER_BIT | GL_COLOR_BUFFER_BIT);

    Gd.bgcol = col;
    Modelshowing = 0;
    gsd_flush();

    return;
}

/*!
  \brief Get aspect value

  \return aspect value
*/
double GS_get_aspect(void)
{
    int left, right, bottom, top;
    GLint tmp[4];

    /* OGLXXX
     * get GL_VIEWPORT:
     * You can probably do better than this.
     */
    glGetIntegerv(GL_VIEWPORT, tmp);
    left = tmp[0];
    right = tmp[0] + tmp[2] - 1;
    bottom = tmp[1];
    top = tmp[1] + tmp[3] - 1;

    G_debug(3, "GS_get_aspect(): left=%d, right=%d, top=%d, bottom=%d",
	    left, right, top, bottom);

    return ((double) (right - left) / (top - bottom));
}

/*!
  \brief Check for transparency

  Disabled.

  \return 1
*/
int GS_has_transparency(void)
{
    /* OGLXXX
     * getgdesc other posiblilties:
     *      glxGetConfig();
     *      glxGetCurrentContext();
     *      glxGetCurrentDrawable();
     * GLint gdtmp;
     * blending is ALWAYS supported.
     * This function returns whether it is enabled.
     * return((glGetIntegerv(GL_BLEND, &gdtmp), gdtmp));
     */

    return (1);
}