c
c     Copyright (c) 1999,2000, Compaq Computer Corporation
c
c     See the file "license.terms" for information on usage and redistribution
c     of this file, and for a DISCLAIMER OF ALL WARRANTIES.
c
      subroutine surf_boundingbox(window_id,eye_radius,eye_psi,
     1                            eye_theta, lfirst)
c
      implicit none
c
c     Draw bounding box around surface quasi-3D plot
c
c     A stripped down version of the SURF code. Just enough to setup
c     to plot the axis bounding box in 3D.
c
      integer window_id
      real eye_radius, eye_psi, eye_theta
      logical lfirst
c
      include 'const.h'
      include 'sysdep.cmn'
      include 'surf.cmn'
      include 'plot.cmn'
c
      real    window(4)
      real    x(5), x_2d(5), y(5), y_2d(5), z(5)
      real    x_bb1(5), x_bb2(5), x_bb3(2), x_bb4(2), x_bb5(2), x_bb6(2)
      real    y_bb1(5), y_bb2(5), y_bb3(2), y_bb4(2), y_bb5(2), y_bb6(2)
      real    viewport(4)
      integer status, surf_xform_gks, n, i, j
      integer xform_number
c
      integer temp_ws_id, lookup_ws_id
c
      integer mode
c
      data window/0.,1.0,0.0,0.9/
      data xform_number/PLOT_XFORM/
      save x_bb1, x_bb2, x_bb3, x_bb4, x_bb5, x_bb6
      save y_bb1, y_bb2, y_bb3, y_bb4, y_bb5, y_bb6
c
c     Start up window
      temp_ws_id = lookup_ws_id(window_id)
      if( temp_ws_id .eq. PLT$EMPTY_WINDOW_SLOT ) return  !  Bad window id
c
c     Activate plot device, do not clear it, do not start segment on
c     workstations
      call activate_plot_device(window_id, .FALSE., .FALSE.)
      call pop_vs_window(window_id)
      call init_device_depend
c
c     Assign default values (from surf_set_default/surf_general.for)
c     View specification data 
      surf_view%xupvec   = 0.0     ! UP VECTOR is the Z axis
      surf_view%yupvec   = 0.0
      surf_view%zupvec   = 1.0
      surf_view%distance = 1.0
c
      surf_view%umin = -1.25
      surf_view%umax =  1.25
      surf_view%vmin = -1.25
      surf_view%vmax =  1.25
c
      surf_global_xfactor = 1.0
      surf_global_yfactor = 1.0
      surf_global_xoffset = 0.0
      surf_global_yoffset = 0.0
c
c     Use 100000 and -100000 to define clipping in the w direction....
c     NO actual clipping!
      surf_view%wmin     = -100000.
      surf_view%wmax     = 100000.
c
      do j=1,4
         do i=1,4
            surf_unity_xform_3d(i,j) = 0.0
         end do
         surf_unity_xform_3d(j,j) = 1.0
      end do
c
c     Define the new viewing transformation... (from surf.for)
      call surf_spherical_to_cartesian(surf_eye%xposition,
     +       surf_eye%yposition, surf_eye%zposition,
     +       eye_radius, eye_psi, eye_theta)
      surf_view%xnormal = 0.5 - surf_eye%xposition
      surf_view%ynormal = 0.5 - surf_eye%yposition
      surf_view%znormal = 0.5 - surf_eye%zposition
c     Scale u-v plane by eye radius to effect zoom */
c     Magic number 8.66 is radius of default position (5,5,5)
c     used in perspective transformation. In parallel transform,
c     we scale by 8.66 to normalize radius and keep it compatible
c     with previous versions of SURF.
      surf_view%umin = surf_view%umin*surf_eye%radius/8.66
      surf_view%umax = surf_view%umax*surf_eye%radius/8.66
      surf_view%vmin = surf_view%vmin*surf_eye%radius/8.66
      surf_view%vmax = surf_view%vmax*surf_eye%radius/8.66
c
c     Set 3D transformations (from surf_set_view/surf_plot3d.for)
      call surf_set_window_3d (surf_view%umin,surf_view%umax,
     +                         surf_view%vmin,surf_view%vmax,
     +                         surf_view%wmin,surf_view%wmax)
      call surf_init_xform_3d
      call surf_calculate_ndc_transform_3d
c     Port only supports parallel projection
      call surf_set_parallel_transform_3d
c
c     First we set the window transformation to map to the plot area.
c     (from surf_set_plotwindow_gks/surf_general.for)
c      
      viewport(1) = window(1)*surf_global_xfactor + surf_global_xoffset
      viewport(2) = window(2)*surf_global_xfactor + surf_global_xoffset
      viewport(3) = window(3)*ratio_y*surf_global_yfactor +
     +                        surf_global_yoffset
      viewport(4) = window(4)*ratio_y*surf_global_yfactor +
     +                        surf_global_yoffset
      status = surf_xform_gks(xform_number,window,viewport)
      if(status .ne. SURF_NORMAL) then
         call surf_error('INVALID TRANSFORMATION FOR PLOT AREA',
     1                   status,1)
         return
      endif
c
c     Turn on complement writing mode
      mode = 1
      call pm_escape_set_writing_mode(mode)
c
c     If not first time, then redraw old bounding box to erase it
      if( .NOT. lfirst ) then
         n = 5
         call gks$polyline(n,x_bb1,y_bb1)
         call gks$polyline(n,x_bb2,y_bb2)
         n = 2
         call gks$polyline(n,x_bb3,y_bb3)
         call gks$polyline(n,x_bb4,y_bb4)
         call gks$polyline(n,x_bb5,y_bb5)
         call gks$polyline(n,x_bb6,y_bb6)
      endif
c
c     Draw bounding box as polylines with (0,0,0) to (1,1,1) as box
c     corners. All surf coordinates are normalized so that the axes
c     boundaries all run from 0.0 to 1.0.
c
      x(1) = 0.0
      x(2) = 1.0
      x(3) = 1.0
      x(4) = 0.0
      x(5) = 0.0
      y(1) = 0.0
      y(2) = 0.0
      y(3) = 1.0
      y(4) = 1.0
      y(5) = 0.0
      z(1) = 0.0
      z(2) = 0.0
      z(3) = 0.0
      z(4) = 0.0
      z(5) = 0.0
      n = 5
      call surf_xform_pline_3d(n,x,y,z,x_2d,y_2d)
      call gks$polyline(n,x_2d,y_2d)
      do 100 i=1,n
         x_bb1(i) = x_2d(i)
         y_bb1(i) = y_2d(i)
100   continue
c
      z(1) = 1.0
      z(2) = 1.0
      z(3) = 1.0
      z(4) = 1.0
      z(5) = 1.0
      call surf_xform_pline_3d(n,x,y,z,x_2d,y_2d)
      call gks$polyline(n,x_2d,y_2d)
      do 101 i=1,n
         x_bb2(i) = x_2d(i)
         y_bb2(i) = y_2d(i)
101   continue
c
      n = 2
      x(1) = 0.0
      x(2) = 0.0
      y(1) = 0.0
      y(2) = 0.0
      z(1) = 0.0
      z(2) = 1.0
      call surf_xform_pline_3d(n,x,y,z,x_2d,y_2d)
      call gks$polyline(n,x_2d,y_2d)
      do 102 i=1,n
         x_bb3(i) = x_2d(i)
         y_bb3(i) = y_2d(i)
102   continue
c
      x(1) = 1.0
      x(2) = 1.0
      y(1) = 0.0
      y(2) = 0.0
      call surf_xform_pline_3d(n,x,y,z,x_2d,y_2d)
      call gks$polyline(n,x_2d,y_2d)
      do 103 i=1,n
         x_bb4(i) = x_2d(i)
         y_bb4(i) = y_2d(i)
103   continue
c
      x(1) = 1.0
      x(2) = 1.0
      y(1) = 1.0
      y(2) = 1.0
      call surf_xform_pline_3d(n,x,y,z,x_2d,y_2d)
      call gks$polyline(n,x_2d,y_2d)
      do 104 i=1,n
         x_bb5(i) = x_2d(i)
         y_bb5(i) = y_2d(i)
104   continue
c
      x(1) = 0.0
      x(2) = 0.0
      y(1) = 1.0
      y(2) = 1.0
      call surf_xform_pline_3d(n,x,y,z,x_2d,y_2d)
      call gks$polyline(n,x_2d,y_2d)
      do 105 i=1,n
         x_bb6(i) = x_2d(i)
         y_bb6(i) = y_2d(i)
105   continue
c
c     Turn off complement writing mode
c
      mode = 0
      call pm_escape_set_writing_mode(mode)
c
      call gks$update_ws(ws_id, GKS$K_PERFORM_FLAG)
      call deactivate_plot_device
c
      return
      end
c
      subroutine surf_xform_pline_3d(n,x,y,z,x_2d,y_2d)
      implicit none
c
      integer n
      real x(n), y(n), z(n), x_2d(n), y_2d(n)
c
c     Transforms a point in 3-D space to one in the view plane.
c
c     (1) change to NDC
c     (2) perform z (w -- depth) clipping - leave the other
c         two axes to the 2D package, after perspective xform.
c     (3) apply perspective transformation
c
      include 'const.h'
      include 'surf.cmn'
      include 'plot3d.cmn'
c
      real     start(4)
      real s(4)          ! points, clipped
      real ss(4)         ! points, after perspective
      integer    i
c
      do 100 i=1,n
         start(1)  = x(i)
         start(2)  = y(i)
         start(3)  = z(i)
c      
c        Change to NDC
         call surf_transform_point_3d (start, surf_ndc_xform_3d, s)
c      
c        Depth clipping omitted, since we don't use it anyway
c        (front/backplane_ndc are set very large). jvf
c      
c        Now transform into 2D space
         call surf_transform_point_3d(s, surf_perspective_xform_3d, ss)
c      
         x_2d(i)  = ss(1)
         y_2d(i)  = ss(2)
100   continue
c
      return
      end

