~hwkrus/f03gl/trunk

! This program plots a function of two variables.  The function, called

! func_to_plot, is an external procedure at the end of the file.  

! This begins with the same module used in the modview example, followed by

! another module for plotting the function, called function_plotter.

! You might want to change default initial settings in modules modview and

! function_plotter.

! William F. Mitchell

! william.mitchell@nist.gov

! Mathematical and Computational Sciences Division

! National Institute of Standards and Technology

! August, 1999

!---------------------------------------------------------------------------

module view_modifier

! This module provides facilities to modify the view in an OpenGL window.

! The mouse buttons and keyboard arrow keys can be used to zoom, pan,

! rotate and change the scale.  A menu or submenu can be used to select which

! buttons perform which function and to reset the view to the initial settings.

! This is limited to one window.

! William F. Mitchell

! william.mitchell@nist.gov

! Mathematical and Computational Sciences Division

! National Institute of Standards and Technology

! April, 1998

! To use this module:

!

! 1) put a USE view_modifier statement in any program unit that calls a

!    procedure in this module

!

! 2) set the initial operation assignments, view and scale below the

!    "Initial configuration" comment below

!

! 3) call view_modifier_init after glutCreateWindow

!    This is a function that returns integer(kind=glcint) menuid.  The menuid

!    is the ID returned by glutCreateMenu.  You can either use the view_modifier

!    menu as your menu by calling glutAttachMenu immediately after

!    view_modifier_init, as in

!       menuid = view_modifier_init()

!       call glutAttachMenu(GLUT_RIGHT_BUTTON)

!    or by using the menuid to attach a submenu to your own menu, as in

!       call glutAddSubMenu("View Modifier",menuid)

!

! 4) in any callback functions that update the display, put

!       call reset_view

!    as the first executable statement

!

! Note that view_modifier_init sets the callback functions for glutMouseFunc,

! glutMotionFunc and glutSpecialFunc, so don't call these yourself

!

! The menu allows you to select what operation is attached to the left and

! middle mouse buttons and arrow keys, reset to the initial view, and quit.

! The right mouse button should be used for the menu.

use opengl_gl

use opengl_glu

use opengl_glut

implicit none

private

public :: view_modifier_init, reset_view

private :: ZOOM, PAN, ROTATE, SCALEX, SCALEY, SCALEZ, RESET, ABOVE, QUIT, &

           PI, &

           left_button_func, middle_button_func, arrow_key_func, &

           init_lookat, init_lookfrom, &

           init_xscale_factor, init_yscale_factor, init_zscale_factor, &

           angle, shift, xscale_factor, yscale_factor, zscale_factor, &

           moving_left, moving_middle, begin_left, begin_middle, &

           cart2sphere, sphere2cart, cart3D_plus_cart3D, cart3D_minus_cart3D, &

           reset_to_init, mouse, motion, arrows, &

           menu_handler, set_left_button, set_middle_button, set_arrow_keys

integer(kind=glcint), parameter :: ZOOM = 1, PAN = 2, ROTATE = 3, SCALEX = 4, &

                      SCALEY = 5, SCALEZ = 6

integer(kind=glcint), parameter :: RESET = 10, ABOVE = 11, QUIT = 12

real(kind=gldouble), parameter :: PI = 3.141592653589793_gldouble

type, private :: cart2D ! 2D cartesian coordinates

   real(kind=gldouble) :: x, y

end type cart2D

type, private :: cart3D ! 3D cartesian coordinates

   real(kind=gldouble) :: x, y, z

end type cart3D

type, private :: sphere3D ! 3D spherical coordinates

   real(kind=gldouble) :: theta, phi, rho

end type sphere3D

type(cart2D), save :: angle

type(cart3D), save :: shift

real(kind=gldouble), save :: xscale_factor, yscale_factor, zscale_factor

logical, save :: moving_left, moving_middle

type(cart2D), save :: begin_left, begin_middle

interface operator(+)

   module procedure cart3D_plus_cart3D

end interface

interface operator(-)

   module procedure cart3D_minus_cart3D

end interface

! ------- Initial configuration -------

! Set the initial operation performed by each button and the arrow keys.

! The operations are ZOOM, PAN, ROTATE, SCALEX, SCALEY, and SCALEZ

integer, save ::   left_button_func = ROTATE, &

                 middle_button_func = ZOOM, &

                     arrow_key_func = PAN

! Set the initial view as the point you are looking at, the point you are

! looking from, and the scale factors

type(cart3D), parameter :: &

     init_lookat = cart3D(0.5_gldouble, 0.5_gldouble, 0.0_gldouble), &

   init_lookfrom = cart3D(5.0_gldouble, 10.0_gldouble, 2.5_gldouble)

real(kind=gldouble), parameter :: &

   init_xscale_factor = 1.0_gldouble, &

   init_yscale_factor = 1.0_gldouble, &

   init_zscale_factor = 1.0_gldouble

! -------- end of Initial configuration ------

contains

!          -------------

subroutine reset_to_init

!          -------------

! This resets the view to the initial configuration

type(sphere3D) :: slookfrom

slookfrom = cart2sphere(init_lookfrom-init_lookat)

angle%x = -180.0_gldouble*slookfrom%theta/PI - 90.0_gldouble

angle%y = -180.0_gldouble*slookfrom%phi/PI

shift%x = 0.0_gldouble

shift%y = 0.0_gldouble

shift%z = -slookfrom%rho

xscale_factor = init_xscale_factor

yscale_factor = init_yscale_factor

zscale_factor = init_zscale_factor

call glutPostRedisplay

return

end subroutine reset_to_init

!          ---------------

subroutine view_from_above

!          ---------------

! This sets the view to be from straight above

type(sphere3D) :: slookfrom

slookfrom = cart2sphere(cart3D(0.0,0.0,1.0))

angle%x = -180.0_gldouble*slookfrom%theta/PI

angle%y = -180.0_gldouble*slookfrom%phi/PI

call glutPostRedisplay

return

end subroutine view_from_above

!          ----------

subroutine reset_view

!          ----------

! This routine resets the view to the current orientation and scale

call glMatrixMode(GL_MODELVIEW)

call glPopMatrix

call glPushMatrix

call glTranslated(shift%x, shift%y, shift%z)

call glRotated(angle%x, 0.0_gldouble, 0.0_gldouble, 1.0_gldouble)

call glRotated(angle%y, cos(PI*angle%x/180.0_gldouble), &

               -sin(PI*angle%x/180.0_gldouble), 0.0_gldouble)

call glTranslated(-init_lookat%x, -init_lookat%y, -init_lookat%z)

call glScaled(xscale_factor,yscale_factor,zscale_factor)

return

end subroutine reset_view

!          -----

subroutine mouse(button, state, x, y) bind(c,name="")

!          -----

integer(kind=glcint), value :: button, state, x, y

! This gets called when a mouse button changes

  if (button == GLUT_LEFT_BUTTON .and. state == GLUT_DOWN) then

    moving_left = .true.

    begin_left = cart2D(x,y)

  endif

  if (button == GLUT_LEFT_BUTTON .and. state == GLUT_UP) then

    moving_left = .false.

  endif

  if (button == GLUT_MIDDLE_BUTTON .and. state == GLUT_DOWN) then

    moving_middle = .true.

    begin_middle = cart2D(x,y)

  endif

  if (button == GLUT_MIDDLE_BUTTON .and. state == GLUT_UP) then

    moving_middle = .false.

  endif

end subroutine mouse

!          ------

subroutine motion(x, y) bind(c,name="")

!          ------

integer(kind=glcint), value :: x, y

! This gets called when the mouse moves

integer :: button_function

type(cart2D) :: begin

real(kind=gldouble) :: factor

! Determine and apply the button function

if (moving_left) then

   button_function = left_button_func

   begin = begin_left

else if(moving_middle) then

   button_function = middle_button_func

   begin = begin_middle

end if

select case(button_function)

case (ZOOM)

   if (y < begin%y) then

      factor = 1.0_gldouble/(1.0_gldouble + .002_gldouble*(begin%y-y))

   else if (y > begin%y) then

      factor = 1.0_gldouble + .002_gldouble*(y-begin%y)

   else

      factor = 1.0_gldouble

   end if

   shift%z = factor*shift%z

case (PAN)

   shift%x = shift%x + .01*(x - begin%x)

   shift%y = shift%y - .01*(y - begin%y)

case (ROTATE)

   angle%x = angle%x + (x - begin%x)

   angle%y = angle%y + (y - begin%y)

case (SCALEX)

   if (y < begin%y) then

      factor = 1.0_gldouble + .002_gldouble*(begin%y-y)

   else if (y > begin%y) then

      factor = 1.0_gldouble/(1.0_gldouble + .002_gldouble*(y-begin%y))

   else

      factor = 1.0_gldouble

   end if

   xscale_factor = xscale_factor * factor

case (SCALEY)

   if (y < begin%y) then

      factor = 1.0_gldouble + .002_gldouble*(begin%y-y)

   else if (y > begin%y) then

      factor = 1.0_gldouble/(1.0_gldouble + .002_gldouble*(y-begin%y))

   else

      factor = 1.0_gldouble

   end if

   yscale_factor = yscale_factor * factor

case (SCALEZ)

   if (y < begin%y) then

      factor = 1.0_gldouble + .002_gldouble*(begin%y-y)

   else if (y > begin%y) then

      factor = 1.0_gldouble/(1.0_gldouble + .002_gldouble*(y-begin%y))

   else

      factor = 1.0_gldouble

   end if

   zscale_factor = zscale_factor * factor

end select

! update private variables and redisplay

if (moving_left) then

   begin_left = cart2D(x,y)

else if(moving_middle) then

   begin_middle = cart2D(x,y)

endif

if (moving_left .or. moving_middle) then

   call glutPostRedisplay

endif

return

end subroutine motion

!          ------

subroutine arrows(key, x, y) bind(c,name="")

!          ------

integer(glcint), value :: key, x, y

! This routine handles the arrow key operations

real(kind=gldouble) :: factor

select case(arrow_key_func)

case(ZOOM)

   select case(key)

   case(GLUT_KEY_DOWN)

      factor = 1.0_gldouble + .02_gldouble

   case(GLUT_KEY_UP)

      factor = 1.0_gldouble/(1.0_gldouble + .02_gldouble)

   case default

      factor = 1.0_gldouble

   end select

   shift%z = factor*shift%z

case(PAN)

   select case(key)

   case(GLUT_KEY_LEFT)

      shift%x = shift%x - .02

   case(GLUT_KEY_RIGHT)

      shift%x = shift%x + .02

   case(GLUT_KEY_DOWN)

      shift%y = shift%y - .02

   case(GLUT_KEY_UP)

      shift%y = shift%y + .02

   end select

case(ROTATE)

   select case(key)

   case(GLUT_KEY_LEFT)

      angle%x = angle%x - 1.0_gldouble

   case(GLUT_KEY_RIGHT)

      angle%x = angle%x + 1.0_gldouble

   case(GLUT_KEY_DOWN)

      angle%y = angle%y + 1.0_gldouble

   case(GLUT_KEY_UP)

      angle%y = angle%y - 1.0_gldouble

   end select

case(SCALEX)

   select case(key)

   case(GLUT_KEY_DOWN)

      factor = 1.0_gldouble/(1.0_gldouble + .02_gldouble)

   case(GLUT_KEY_UP)

      factor = 1.0_gldouble + .02_gldouble

   case default

      factor = 1.0_gldouble

   end select

   xscale_factor = xscale_factor * factor

case(SCALEY)

   select case(key)

   case(GLUT_KEY_DOWN)

      factor = 1.0_gldouble/(1.0_gldouble + .02_gldouble)

   case(GLUT_KEY_UP)

      factor = 1.0_gldouble + .02_gldouble

   case default

      factor = 1.0_gldouble

   end select

   yscale_factor = yscale_factor * factor

case(SCALEZ)

   select case(key)

   case(GLUT_KEY_DOWN)

      factor = 1.0_gldouble/(1.0_gldouble + .02_gldouble)

   case(GLUT_KEY_UP)

      factor = 1.0_gldouble + .02_gldouble

   case default

      factor = 1.0_gldouble

   end select

   zscale_factor = zscale_factor * factor

end select

call glutPostRedisplay

return

end subroutine arrows

!          ------------

subroutine menu_handler(value) bind(c,name="")

!          ------------

integer(kind=glcint), value :: value

! This routine handles the first level entries in the menu

select case(value)

case(RESET)

   call reset_to_init

case(ABOVE)

   call view_from_above

case(QUIT)

   stop

end select

return

end subroutine menu_handler

!          ---------------

subroutine set_left_button(value) bind(c,name="")

!          ---------------

integer(kind=glcint), value :: value

! This routine sets the function of the left button as given by menu selection

left_button_func = value

return

end subroutine set_left_button

!          -----------------

subroutine set_middle_button(value)  bind(c,name="")

!          -----------------

integer(kind=glcint), value :: value

! This routine sets the function of the middle button as given by menu selection

middle_button_func = value

return

end subroutine set_middle_button

!          --------------

subroutine set_arrow_keys(value)  bind(c,name="")

!          --------------

integer(kind=glcint), value :: value

! This routine sets the function of the arrow keys as given by menu selection

arrow_key_func = value

return

end subroutine set_arrow_keys

!        ------------------

function view_modifier_init() result(menuid)

!        ------------------

integer(kind=glcint) :: menuid

! This initializes the view modifier variables and sets initial view.

! It should be called immediately after glutCreateWindow

integer(kind=glcint) :: button_left, button_middle, arrow_keys

! set the callback functions

call glutMouseFunc(mouse)

call glutMotionFunc(motion)

call glutSpecialFunc(arrows)

! create the menu

button_left = glutCreateMenu(set_left_button)

call glutAddMenuEntry("rotate"//char(0),ROTATE)

call glutAddMenuEntry("zoom"//char(0),ZOOM)

call glutAddMenuEntry("pan"//char(0),PAN)

call glutAddMenuEntry("scale x"//char(0),SCALEX)

call glutAddMenuEntry("scale y"//char(0),SCALEY)

call glutAddMenuEntry("scale z"//char(0), SCALEZ)

button_middle = glutCreateMenu(set_middle_button)

call glutAddMenuEntry("rotate"//char(0),ROTATE)

call glutAddMenuEntry("zoom"//char(0),ZOOM)

call glutAddMenuEntry("pan"//char(0),PAN)

call glutAddMenuEntry("scale x"//char(0),SCALEX)

call glutAddMenuEntry("scale y"//char(0),SCALEY)

call glutAddMenuEntry("scale z"//char(0), SCALEZ)

arrow_keys = glutCreateMenu(set_arrow_keys)

call glutAddMenuEntry("rotate"//char(0),ROTATE)

call glutAddMenuEntry("zoom"//char(0),ZOOM)

call glutAddMenuEntry("pan"//char(0),PAN)

call glutAddMenuEntry("scale x"//char(0),SCALEX)

call glutAddMenuEntry("scale y"//char(0),SCALEY)

call glutAddMenuEntry("scale z"//char(0), SCALEZ)

menuid = glutCreateMenu(menu_handler)

call glutAddSubMenu("left mouse button"//char(0),button_left)

call glutAddSubMenu("middle mouse button"//char(0),button_middle)

call glutAddSubMenu("arrow keys"//char(0),arrow_keys)

call glutAddMenuEntry("reset to initial view"//char(0),RESET)

call glutAddMenuEntry("view from above"//char(0),ABOVE)

call glutAddMenuEntry("quit"//char(0),QUIT)

! set the perspective

call glMatrixMode(GL_PROJECTION)

call gluPerspective(10.0_gldouble, 1.0_gldouble, 0.1_gldouble, 200.0_gldouble)

! set the initial view

call glMatrixMode(GL_MODELVIEW)

call glPushMatrix

call reset_to_init

return

end function view_modifier_init

!        -----------

function sphere2cart(spoint) result(cpoint)

!        -----------

type(sphere3D), intent(in) :: spoint

type(cart3D) :: cpoint

! This converts a 3D point from spherical to cartesean coordinates

real(kind=gldouble) :: t,p,r

t=spoint%theta

p=spoint%phi

r=spoint%rho

cpoint%x = r*cos(t)*sin(p)

cpoint%y = r*sin(t)*sin(p)

cpoint%z = r*cos(p)

return

end function sphere2cart

!        -----------

function cart2sphere(cpoint) result(spoint)

!        -----------

type(cart3D), intent(in) :: cpoint

type(sphere3D) :: spoint

! This converts a 3D point from cartesean to spherical coordinates

real(kind=gldouble) :: x,y,z

x=cpoint%x

y=cpoint%y

z=cpoint%z

spoint%rho = sqrt(x*x+y*y+z*z)

if (x==0.0_gldouble .and. y==0.0_gldouble) then

   spoint%theta = 0.0_gldouble

else

   spoint%theta = atan2(y,x)

end if

if (spoint%rho == 0.0_gldouble) then

   spoint%phi = 0.0_gldouble

else

   spoint%phi = acos(z/spoint%rho)

endif

return

end function cart2sphere

!        ------------------

function cart3D_plus_cart3D(cart1,cart2) result(cart3)

!        ------------------

type(cart3D), intent(in) :: cart1, cart2

type(cart3D) :: cart3

! Compute the sum of two 3D cartesean points

cart3%x = cart1%x + cart2%x

cart3%y = cart1%y + cart2%y

cart3%z = cart1%z + cart2%z

return

end function cart3D_plus_cart3D

!        -------------------

function cart3D_minus_cart3D(cart1,cart2) result(cart3)

!        -------------------

type(cart3D), intent(in) :: cart1, cart2

type(cart3D) :: cart3

! Compute the difference of two 3D cartesean points

cart3%x = cart1%x - cart2%x

cart3%y = cart1%y - cart2%y

cart3%z = cart1%z - cart2%z

return

end function cart3D_minus_cart3D

end module view_modifier

!---------------------------------------------------------------------------

module function_plotter

use opengl_gl

use opengl_glut

use view_modifier

implicit none

private

public :: display, draw_func, menu_handler2, make_menu

! symbolic constants

integer, parameter :: surfgrid_toggle = 1, &

                      surfsolid_toggle = 2, &

                      contour_toggle = 3, &

                      quit_selected = 4

integer, parameter :: set_nx = 1, &

                      set_ny = 2, &

                      set_ncontour = 3, &

                      set_contour_val = 4, &

                      set_xrange = 5, &

                      set_yrange = 6, &

                      reset_params = 7

integer, parameter :: black_contour = 1, &

                      rainbow_contour = 2

integer, parameter :: white_surface = 1, &

                      red_surface = 2, &

                      rainbow_surface = 3

! Default initial settings

integer, parameter :: init_ngridx = 40, &

                      init_ngridy = 40, &

                      init_num_contour = 20, &

                      init_contour_color = black_contour, &

                      init_surface_color = rainbow_surface

real(GLDOUBLE), parameter :: init_minx = 0.0_GLDOUBLE, &

                             init_maxx = 1.0_GLDOUBLE, &

                             init_miny = 0.0_GLDOUBLE, &

                             init_maxy = 1.0_GLDOUBLE

logical, parameter :: init_draw_surface_grid = .false., &

                      init_draw_surface_solid = .true., &

                      init_draw_contour = .true.

! Current settings

integer :: ngridx = init_ngridx, &

           ngridy = init_ngridy, &

           num_contour = init_num_contour, &

           contour_color = init_contour_color, &

           surface_color = init_surface_color

real(GLDOUBLE) :: minx = init_minx, &

                  maxx = init_maxx, &

                  miny = init_miny, &

                  maxy = init_maxy, &

                  minz = 0.0_GLDOUBLE, &

                  maxz = 0.0_GLDOUBLE

logical :: draw_surface_grid = init_draw_surface_grid, &

           draw_surface_solid = init_draw_surface_solid, &

           draw_contour = init_draw_contour, &

           contour_values_given = .false.

real(GLDOUBLE), allocatable :: actual_contours(:)

contains

subroutine display() bind(c)

! This gets called when the display needs to be redrawn

call reset_view

call glClear(ior(GL_COLOR_BUFFER_BIT,GL_DEPTH_BUFFER_BIT))

call glCallList(1)

call glutSwapBuffers

return

end subroutine display

subroutine draw_func

real(GLDOUBLE) :: gridx(0:ngridx),gridy(0:ngridy),zval(0:ngridy,2)

integer :: i,j,k,cont

real(GLDOUBLE) :: x1,x2,x3,xt,y1,y2,y3,yt,z1,z2,z3,zt

real(GLDOUBLE) :: frac,xcross1,xcross2,ycross1,ycross2

real(GLDOUBLE) :: contour_value(num_contour)

real(GLFLOAT) :: color(4), normal(3), &

                 red(4) = (/1.0,0.0,0.0,1.0/), &

                 black(4) = (/0.0,0.0,0.0,1.0/), &

                 white(4) = (/1.0,1.0,1.0,1.0/)

real(GLDOUBLE), external :: func_to_plot

! prepare to make a new display list

call reset_view

call glDeleteLists(1_gluint, 1_glsizei)

call glNewList(1_gluint, gl_compile_and_execute)

! set the grid points

gridx = (/ ((minx + i*(maxx-minx)/ngridx),i=0,ngridx) /)

gridy = (/ ((miny + i*(maxy-miny)/ngridy),i=0,ngridy) /)

! if this is the first call and either rainbow coloring of a solid surface

! or contours are being drawn, minz and maxz need to be set

if (minz == 0.0_GLDOUBLE .and. maxz == 0.0_GLDOUBLE .and. &

    ( (draw_surface_solid .and. surface_color == rainbow_surface) .or. &

      draw_contour ) ) then

   do i=0,ngridx

      do j=0,ngridy

         z1 = func_to_plot(gridx(i),gridy(j))

         minz = min(z1,minz)

         maxz = max(z1,maxz)

      end do

   end do

endif

! draw the solid surface

if (draw_surface_solid) then

   call glPolygonMode(gl_front_and_back, gl_fill)

   call glBegin(gl_triangles)

! set the color for a red or white surface.  For white, set the lighting

! such that there is uniform brightness.

   if (surface_color == red_surface) then

      call glMaterialfv(gl_front_and_back,gl_ambient_and_diffuse,red)

   elseif (surface_color == white_surface) then

      call glDisable(gl_light0)

      call glLightModelfv(gl_light_model_ambient, (/1.,1.,1.,1./))

      call glMaterialfv(gl_front_and_back,gl_ambient_and_diffuse,white)

   endif

! compute the function values for the first line of points

   do j=0,ngridy

      zval(j,2) = func_to_plot(gridx(0),gridy(j))

   end do

! for each x grid interval...

   do i=1,ngridx

! copy left side function values from the right side of the previous interval

      zval(:,1) = zval(:,2)

! compute the function values for the right side of the interval

      do j=0,ngridy

         zval(j,2) = func_to_plot(gridx(i),gridy(j))

      end do

      minz = min(minz,minval(zval))

      maxz = max(maxz,maxval(zval))

! for each y grid interval ...

      do j=1,ngridy

! add two triangles to the display list

! Before each triangle, set the normal.  Before each vertex, set the

! color if we're coloring by height

         normal = normcrossprod((/gridx(i-1),gridx(i),gridx(i)/), &

                                (/gridy(j-1),gridy(j-1),gridy(j)/), &

                                (/zval(j-1,1),zval(j-1,2),zval(j,2)/))

         call glNormal3fv(normal)

         if (surface_color == rainbow_surface) then

            call get_rainbow(zval(j-1,1),minz,maxz,color)

            call glMaterialfv(gl_front_and_back,gl_ambient_and_diffuse,color)

         endif

         call glvertex3d(gridx(i-1),gridy(j-1),zval(j-1,1))

         if (surface_color == rainbow_surface) then

            call get_rainbow(zval(j-1,2),minz,maxz,color)