c      main program
       parameter(nx=8,ny=2,nx2=nx*ny,ico=7)
       dimension p(nx2)
c
       pi=3.1415926
c
       do 10 j=1,ny
       do 10 i=1,nx
       i1=j+2*(i-1)
       p(i1)=float(i1)/float(nx2)
   10  continue
c      
c
       call initvrml
c
       call facebe
c
       call colorp(nx2,p,ico)
c
       call coordbe
c
       call coordp(nx2,p)
c
       call coorden(nx)
       call faceen
c
       stop
       end
c      subroutine initvrml
       subroutine initvrml
       write(10,10)
   10  format("#VRML V2.0 utf8")
       return
       end
       subroutine colorp(nx,p,ico)
       dimension p(1)
       r=1.0
       g=1.0
       b=1.0
       do 100 i=1,nx
       x1=p(i)
c
       a0=1.20
       a0=1.00
       if(ico.eq.2) go to 32
       if(ico.eq.3) go to 33
       if(ico.eq.4) go to 34
       if(ico.eq.5) go to 35
       if(ico.eq.6) go to 36
       if(ico.eq.7) go to 37
       r=a0-(a0-0.5)*x1
       g=a0-(a0-0.5)*x1
       b=a0-(a0-0.5)*x1
       go to 40
   32  continue
       r=a0-a0*x1
       g=1.0
       b=a0-a0*x1
       go to 40
   33  continue
       r=1.0
       g=a0-a0*x1
       b=a0-a0*x1
       go to 40
   34  continue
       a0=2.0
       y1=x1*(a0+4.0)
       r=a0-y1
       if(y1.ge.a0) r=y1-(a0+1.0)
       g=a0+3.0-y1
       b=a0+1.0-y1
       if(y1.ge.a0+2.0) b=y1-(a0+3.0)
       go to 40
   35  continue
       a0=2.0
       y1=x1*(a0+4.0)
       b=a0-y1
       if(y1.ge.a0) b=y1-(a0+1.0)
       r=a0+3.0-y1
       g=a0+1.0-y1
       if(y1.ge.a0+2.0) g=y1-(a0+3.0)
       go to 40
   36  continue
       a0=2.0
       y1=x1*(a0+4.0)
       r=a0-y1
       if(y1.ge.a0) r=y1-(a0+1.0)
       b=a0+3.0-y1
       g=a0+1.0-y1
       if(y1.ge.a0+2.0) g=y1-(a0+3.0)
       go to 40
   37  continue
       x1=10.0*x1
       b=x1
       if(x1.ge.3.0) b=4.0-x1
       if(x1.ge.6.0) b=0.25*(x1-6.0)
       g=x1-2.0
       if(x1.ge.5.0) g=6.0-x1
       if(x1.ge.6.0) g=0.25*(x1-6.0)
       r=x1
       if(x1.ge.1.0) r=2.0-x1
       if(x1.ge.4.0) r=x1-4.0
       go to 40
c
   40  continue
       r=amax1(r,0.0)
       g=amax1(g,0.0)
       b=amax1(b,0.0)
       r=amin1(r,1.0)
       g=amin1(g,1.0)
       b=amin1(b,1.0)
c
       if(i.ne.nx) write(10,20) r,g,b
       if(i.eq.nx) write(10,30) r,g,b
   20  format(f3.1,1x,f3.1,1x,f3.1,",")
   30  format(f3.1,1x,f3.1,1x,f3.1)
  100  continue
       return
       end
       subroutine coordp(nx,p)
       dimension p(1)
       x=1.0
       y=1.0
       z=1.0
       xb=0.0
       yb=0.0
       dx=1.0
       dy=1.0
       do 100 i=1,nx
       i1=(i-1)/2
       i2=i-i1*2
       x=xb+dx*float(i1)
       y=yb+dy*float(i2-1)
       z=0.0
       if(i.ne.nx) write(10,20) x,y,z
       if(i.eq.nx) write(10,30) x,y,z
   20  format(f3.1,1x,f3.1,1x,f3.1,",")
   30  format(f3.1,1x,f3.1,1x,f3.1)
  100  continue
       return
       end
       subroutine facebe
       write(10,21)
   21  format("Shape {")
       write(10,22)
   22  format("geometry IndexedFaceSet {")
       write(10,23)
   23  format("color Color {")
       write(10,24)
   24  format("color [")
       return
       end
       subroutine coordbe
       write(10,21)
   21  format("]")
       write(10,22)
   22  format("}")
       write(10,23)
   23  format("coord Coordinate {")
       write(10,24)
   24  format("point [")
       return
       end
       subroutine coorden(nx)
       dimension ia(2000)
       write(10,21)
   21  format("]")
       write(10,22)
   22  format("}")
       write(10,23)
   23  format("colorPerVertex TRUE")
       write(10,24)
   24  format("coordIndex [")
c
       do 10 i=1,nx
       i1=1+5*(i-1)
       i2=2+5*(i-1)
       i3=3+5*(i-1)
       i4=4+5*(i-1)
       i5=5+5*(i-1)
       ia(i1)=0+2*(i-1)
       ia(i2)=1+2*(i-1)
       ia(i3)=3+2*(i-1)
       ia(i4)=2+2*(i-1)
       ia(i5)=-1
   10  continue
c
       na2=5*(nx-1)-1
       write(10,102) (ia(i),i=1,na2)
  102  format(4(4i2,i3),2i2,3i3,220i3,2i3,2i4,i3,
     1        449(4i4,i3),2i4,2i5,i3,3000(4i5,i3))
c
       return
       end
       subroutine faceen
       write(10,21)
   21  format("]")
       write(10,22)
   22  format("solid FALSE")
       write(10,23)
       write(10,23)
   23  format("}")
       return
       end