c-----------------------------------------------------------
c     subrtn package2(pk2subrtn.f)
c     
c      update  2001-10-12
c
c-----------------------------------------------------------
c-----------------------------------------------------------
       subroutine initvrml
c-----------------------------------------------------------
       write(10,10)
   10  format("#VRML V2.0 utf8")
       return
       end
c-----------------------------------------------------------
       subroutine box1(r,g,b,sz)
c-----------------------------------------------------------
       x=sz
       y=sz
       z=sz
       write(10,12)
   12  format("Transform {")
       write(10,14)
   14  format("children [")
       write(10,16)
   16  format("Shape {")
       write(10,18)
   18  format("appearance Appearance {")
       write(10,20)
   20  format("material Material {")
       write(10,22) r,g,b
   22  format("diffuseColor",f5.2,f5.2,f5.2)
       write(10,24)
   24  format("}",/"}")
       write(10,26)
   26  format("geometry Box {")
       write(10,28) x,y,z
   28  format("size",f5.2,f5.2,f5.2)
       write(10,30)
   30  format("}",/"}",/"]",/"}")
c
       return
       end
c------------------------------------------------------
       subroutine box2(r,g,b,sx,sy,sz)
c-----------------------------------------------------
       write(10,12)
   12  format("Transform {")
       write(10,14)
   14  format("children [")
       write(10,16)
   16  format("Shape {")
       write(10,18)
   18  format("appearance Appearance {")
       write(10,20)
   20  format("material Material {")
       write(10,22) r,g,b
   22  format("diffuseColor",f5.2,f5.2,f5.2)
       write(10,24)
   24  format("}",/"}")
       write(10,26)
   26  format("geometry Box {")
       write(10,28) sx,sy,sz
   28  format("size",f5.2,f5.2,f5.2)
       write(10,30)
   30  format("}",/"}",/"]",/"}")
c
       return
       end
c--------------------------------------------------------
       subroutine cone1k1(rd,h,r,g,b,x,y,z)
c--------------------------------------------------------
       write(10,21) x,y,z
   21  format("Transform {",/"translation ",f5.2,f5.2,f5.2)
       write(10,22)
   22  format("children [")
       write(10,23)
   23  format("Shape {")
       write(10,20)
   20  format("appearance Appearance {")
       write(10,25)
   25  format("material Material {")
       write(10,26) r,g,b
   26  format("diffuseColor",f5.2,f5.2,f5.2)
       write(10,27)
   27  format("}",/"}")
       write(10,28)
   28  format("geometry Cone {")
       write(10,29) rd
   29  format("bottomRadius",f5.2)
       write(10,30) h
   30  format("height",f5.2)
       write(10,31)
   31  format("}",/"}",/"]",/"}")
c
       return
       end
c-----------------------------------------------------
       subroutine cylinder1(rd,h,r,g,b)
c-----------------------------------------------------
       write(10,21)
   21  format("Transform {")
       write(10,22)
   22  format("children [")
       write(10,23)
   23  format("Shape {")
       write(10,20)
   20  format("appearance Appearance {")
       write(10,25)
   25  format("material Material {")
       write(10,26) r,g,b
   26  format("diffuseColor",f5.2,f5.2,f5.2)
       write(10,27)
   27  format("}",/"}")
       write(10,28)
   28  format("geometry Cylinder {")
       write(10,29) rd
   29  format("radius",f5.2)
       write(10,30) h
   30  format("height",f5.2)
       write(10,31)
   31  format("}",/"}",/"]",/"}")
c
       return
       end
c------------------------------------------------------
       subroutine sphere2(r,g,b,rd)
c------------------------------------------------------
       write(10,21)
   21  format("Transform {")
       write(10,23)
   23  format("children [")
       write(10,24)
   24  format("Shape {")
       write(10,25)
   25  format("appearance Appearance {")
       write(10,26)
   26  format("material Material {")
       write(10,27) r,g,b
   27  format("diffuseColor ",f5.2,f5.2,f5.2)
       write(10,28)
   28  format("}",/"}")
       write(10,29)
   29  format("geometry Sphere {")
       write(10,30) rd
   30  format("radius ",f4.2)
       write(10,31)
   31  format("}",/"}",/"]",/"}")
c
       return
       end
c------------------------------------------------------
       subroutine text1(x,y,z,n,r,g,b,sz,ch1)
c------------------------------------------------------
       character ch1(10)*40
       write(10,11)
   11  format("Transform {")
       write(10,12) x,y,z
   12  format("translation ",f6.2,f6.2,f6.2)
c      write(10,13) xr,yr,zr,rd
c  13  format("rotation ",f7.3,f7.3,f7.3,f7.3)
c      write(10,14) xs,ys,zs
c  14  format("scale ",f5.2,f5.2,f5.2)
       write(10,15)
   15  format("children [",/"Shape {",/"appearance Appearance {")
       write(10,16)
   16  format("material Material {")
       write(10,17) r,g,b
   17  format("diffuseColor",f5.2,f5.2,f5.2)
       write(10,18)
   18  format("}",/"}",/"geometry Text {",/"string [")
c
       do 21 i=1,n
       if (i.eq.n) then
       write(10,19) ch1(i)
   19  format('"',a40,'"'," ]")
       else
       write(10,20) ch1(i)
   20  format('"',a40,'"',",")
       end if
   21  continue
c
       write(10,22) sz
   22  format("fontStyle FontStyle {",/"size ",f5.1)  
       write(10,23)
   23  format("}",/"}",/"}",/"]",/"}")
       return
       end
c-------------------------------------------------------
       subroutine line1(last,po,r,g,b)
c------------------------------------------------------
       dimension po(1)
c
       call linebek1(last,r,g,b)
c
       do 10 i=1,last
       i1=3*(i-1)+1
       i2=i1+1
       i3=i2+1
       x=po(i1)
       y=po(i2)
       z=po(i3)
       call linep(x,y,z,i,last)
   10  continue
c
       call lineenk1(last)
c
       return
       end
c------------------------------------------------------
       subroutine linebe(r,g,b)
c------------------------------------------------------
c
       write(10,21)
   21  format("Shape {")
       write(10,22)
   22  format("geometry IndexedLineSet {")
       write(10,23)
   23  format("color Color {")
       write(10,20) r,g,b
   20  format("color [ ",f3.1,1x,f3.1,1x,f3.1," ]")
       write(10,25)
   25  format("}")
       write(10,26)
   26  format("coord Coordinate {")
       write(10,27)
   27  format("point [")
       return
       end
c------------------------------------------------------
       subroutine linebek1(last,r,g,b)
c------------------------------------------------------
c
       write(10,21)
   21  format("Shape {")
       write(10,22)
   22  format("geometry IndexedLineSet {")
       write(10,23)
   23  format("color Color {",/"color [")
c
       do 25 i=1,last
       write(10,24) r,g,b
   24  format(f3.1,1x,f3.1,1x,f3.1)
   25  continue
c
       write(10,26)
   26  format("]")
       write(10,27)
   27  format("}")
       write(10,28)
   28  format("coord Coordinate {")
       write(10,29)
   29  format("point [")
       return
       end
c-------------------------------------------------------
       subroutine linep(x,y,z,ii,last)
c-------------------------------------------------------
       if(ii.lt.last) then
       write(10,10) x,y,z
   10  format(f6.2,1x,f6.2,1x,f6.2,",")
       else
       write(10,20) x,y,z
   20  format(f6.2,1x,f6.2,1x,f6.2)
       end if
       return
       end
c-------------------------------------------------------
       subroutine lineen(last)
c-------------------------------------------------------
       integer ic(3000)
c
       do 10 i=1,3000
       ic(i)=i-1
   10  continue
       write(10,21)
   21  format("]")
       write(10,22)
   22  format("}")
       write(10,23)
   23  format("colorIndex [ 0 ]")
       write(10,24)
   24  format("coordIndex [")
c
       write(10,20) (ic(i),i=1,last)
   20  format(10i2,90i3,900i4/(2000i5))
c
       write(10,25)
   25  format("]")
       write(10,26)
   26  format("colorPerVertex FALSE")
       write(10,27)
   27  format("}")
       write(10,28)
   28  format("}")
       return
       end
c-------------------------------------------------------
       subroutine lineenk1(last)
c-------------------------------------------------------
       integer ic(3000)
c
       do 10 i=1,3000
       ic(i)=i-1
   10  continue
       write(10,21)
   21  format("]")
       write(10,22)
   22  format("}")
       write(10,24)
   24  format("coordIndex [")
c
       write(10,20) (ic(i),i=1,last)
   20  format(10i2,90i3,900i4/(2000i5))
c
       write(10,25)
   25  format("]")
       write(10,27)
   27  format("}")
       write(10,28)
   28  format("}")
       return
       end
c-------------------------------------------------------
       subroutine line2(last,po,sco,eco)
c-------------------------------------------------------
       dimension po(1)
       integer   sco eco
c
       call linebe2
       call colorp1k1(last,sco,eco)
       call coordbek1
c
       do 10 i=1,last
       i1=3*(i-1)+1
       i2=i1+1
       i3=i2+1
       x=po(i1)
       y=po(i2)
       z=po(i3)
       call linep(x,y,z,i,last)
   10  continue
c
       call lineenk1(last)
c
       return
       end
c-------------------------------------------------------
      subroutine linebe2
c-------------------------------------------------------
       write(10,21)
   21  format("Shape {")
       write(10,22)
   22  format("geometry IndexedLineSet {")
       write(10,23)
   23  format("color Color {")
       write(10,24)
   24  format("color [ ")
       return
       end
c--------------------------------------------------------
       subroutine colorp1k1(nx,sco,eco)
c-------------------------------------------------------
       integer   sco eco
       pp=1.0/float(nx)
c
       do 100 i=1,nx
       ppi=pp*float(i-1)
       r=1.0
       g=1.0
       b=1.0
c
       if(sco.eq.2) go to 22
       if(sco.eq.3) go to 23
       if(sco.eq.4) go to 24
       if(sco.eq.5) go to 25
       if(sco.eq.6) go to 26
       if(sco.eq.7) go to 27
c
   22  b=0.0
       go to 30
   23  g= 0.0
       go to 30
   24  g=0.0
       b=0.0
       go to 30
   25  r=0.0
       go to 30
   26  r=0.0
       b=0.0
       go to 30
   27  r=0.0
       g=0.0
c
   30  if(eco.eq.2) go to 32
       if(eco.eq.3) go to 33
       if(eco.eq.4) go to 34
       if(eco.eq.5) go to 35
       if(eco.eq.6) go to 36
       if(eco.eq.7) go to 37
c
       r=r+ppi
       g=g+ppi
       b=b+ppi
   32  r=r+ppi
       g=g+ppi
       b=b-ppi
       go to 40
   33  r=r+ppi
       g=g-ppi
       b=b+ppi
       go to 40
   34  r=r+ppi
       g=g-ppi
       b=b-ppi
       go to 40
   35  r=r-ppi
       g=g+ppi
       b=b+ppi
       go to 40
   36  r=r-ppi
       g=g+ppi
       b=b-ppi
       go to 40
   37  r=r-ppi
       g=g-ppi
       b=b+ppi
c
   40  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,70) r,g,b
       if(i.eq.nx) write(10,80) r,g,b
   70  format(f4.2,1x,f4.2,1x,f4.2,",")
   80  format(f4.2,1x,f4.2,1x,f4.2)
  100  continue
       return
       end
c------------------------------------------------------
       subroutine coordbek1
c------------------------------------------------------
       write(10,25)
   25  format("]",/"}")
       write(10,26)
   26  format("coord Coordinate {")
       write(10,27)
   27  format("point [")
       return
       end
c------------------------------------------------------
       subroutine point1(ln,po,r,g,b)
c-------------------------------------------------------
       dimension po(1)
c
       call pointbe(ln,r,g,b)
c
       do 10 i=1,ln
       i1=3*(i-1)+1
       i2=i1+1
       i3=i2+1
       x=po(i1)
       y=po(i2)
       z=po(i3)
       call linep(x,y,z,i,ln)
   10  continue
c
       call pointen
c
       return
       end
c---------------------------------------------------------
      subroutine pointbe(ln,r,g,b)
c--------------------------------------------------------
       write(10,21)
   21  format("Shape {")
       write(10,22)
   22  format("geometry PointSet {")
       write(10,23)
   23  format("color Color {",/"color [")
c
       do 25 i=1,ln
       write(10,24) r,g,b
   24  format(f3.1,1x,f3.1,1x,f3.1)
   25  continue
c
       write(10,26)
   26  format("]")
       write(10,27)
   27  format("}")
       write(10,28)
   28  format("coord Coordinate {")
       write(10,29)
   29  format("point [")
       return
       end
c--------------------------------------------------------
       subroutine pointen
c--------------------------------------------------------
       write(10,21)
   21  format("]")
       write(10,22)
   22  format("}")
       write(10,27)
   27  format("}")
       write(10,28)
   28  format("}")
       return
       end
c--------------------------------------------------------
       subroutine point2(last,po,sco,eco)
c-------------------------------------------------------
       dimension po(1)
       integer   sco eco
c
       call pointbe2
       call colorp1k1(last,sco,eco)
       call coordbek1
c
       do 10 i=1,last
       i1=3*(i-1)+1
       i2=i1+1
       i3=i2+1
       x=po(i1)
       y=po(i2)
       z=po(i3)
       call linep(x,y,z,i,last)
   10  continue
c
       call pointen
c
       return
       end
c---------------------------------------------------------
      subroutine pointbe2
c--------------------------------------------------------
       write(10,21)
   21  format("Shape {")
       write(10,22)
   22  format("geometry PointSet {")
       write(10,23)
   23  format("color Color {")
       write(10,24)
   24  format("color [ ")
       return
       end
c-------------------------------------------------------
       subroutine face1(ln,po,r,g,b)
c-------------------------------------------------------
       call facebek1(r,g,b)
       call coordp2k1(ln,po)
       call coorden1k1(ln)
       call faceen2
c
       return
       end
c-------------------------------------------------------
       subroutine imagetx1(ln,po,ch1,n)
c--------------------------------------
       call imagebe(ch1,n)
       call coordp2k1(ln,po)
       call imageen(ln)
c
       return
       end
c------------------------------------------------------
       subroutine coordp2k1(ln,po)
c-------------------------------------------------------
       dimension po(1)
c
       do 100 i=1,ln
       do 100 j=1,2
       i1=3*(i-1)+1+3*ln*(j-1)
       i2=i1+1
       i3=i2+1
       x=po(i1)
       y=po(i2)
       z=po(i3)
       call linep2k1(x,y,z,i,j,ln)
  100  continue
c
       return
       end
c------------------------------------------------< 1 color >
      subroutine facebek1(r,g,b)
c------------------------------------------------
       write(10,21)
   21  format("Shape {")
       write(10,22)
   22  format("geometry IndexedFaceSet {")
       write(10,23)
   23  format("color Color {",/"color [ ")
       write(10,20) r,g,b
   20  format(f3.1,1x,f3.1,1x,f3.1)
       write(10,25)
   25  format("]",/"}")
       write(10,26)
   26  format("coord Coordinate {")
       write(10,27)
   27  format("point [")
       return
       end
c----------------------------------------------------
       subroutine faceen2
c----------------------------------------------------
       write(10,21)
   21  format("]")
       write(10,22)
   22  format("solid FALSE")
       write(10,23)
   23  format("}",/"}")
       return
       end
c--------------------------------------------------------
       subroutine linep2k1(x,y,z,i,j,last)
c--------------------------------------------------------
       if(i.eq.last.and.j.eq.2) then
       write(10,10) x,y,z
   10  format(f6.2,1x,f6.2,1x,f6.2)
       else
       write(10,20) x,y,z
   20  format(f6.2,1x,f6.2,1x,f6.2,",")
       end if
       return
       end
c---------------------------------------------------------
       subroutine coorden1k1(ln)
c---------------------------------------------------------
       dimension ia(2000)
       write(10,21)
   21  format("]")
       write(10,22)
   22  format("}")
       write(10,23)
   23  format("colorPerVertex FALSE")
       write(10,24)
   24  format("colorIndex [ ")
       ii=(ln-1)*2
       do 26 i=1,ii
       write(10,25)
   25  format(" 0")
   26  continue
       write(10,27)
   27  format(" ]")
       write(10,29)
   29  format("coordIndex [")
c
       jj=ln*2
       do 100 i=1,jj
       i1=1+4*(i-1)
       i2=2+4*(i-1)
       i3=3+4*(i-1)
       i4=4+4*(i-1)
       ia(i1)=0+i-1
       ia(i2)=1+i-1
       ia(i3)=2+i-1
       ia(i4)=-1
  100  continue
c
       jj=4*(ln-1)*2-1
       write(10,102) (ia(i),i=1,jj)
  102  format(8(3i2,i3),2i2,2i3,i2,3i3,88i3,2i3,i4,i3,i3,2i4,i3,
     1        898(3i4,i3),2i4,2i5,i3,i4,2i5,i3,2000(3i5,i3))
       return
       end
c-------------------------------------------------------
       subroutine triangm(nx,ny,u,p)      
c-------------------------------------------------------
       dimension p(1),u(1)
       ny1=ny-1
c
       do 10 j=1,ny1
       call triangb(nx,j,u,p)
   10  continue  
c
       return
       end
c---------------------------------------------------------
       subroutine triangb(nx,ny,u,p)
c---------------------------------------------------------
       dimension p(1),u(1)
c
       call facebe
       call colorp2(nx,ny,u)
       call coordbe
       call coordp2(nx,ny,p)
       call coorden1(nx)
       call ifaceen 
c
       return
       end
c--------------------------------------------------
       subroutine facebe
c--------------------------------------------------
       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
c-------------------------------------------------------
       subroutine colorp2(nx,ny,u)
c--------------------------------------------------------
       dimension u(1)
c
       do 100 i=1,nx
       do 100 j=1,2
       j1=j+ny-1
c
       i1=3*(i-1)+1+3*nx*(j1-1)
       i2=i1+1
       i3=i1+2
       r=u(i1)
       g=u(i2)
       b=u(i3)
c
       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.eq.nx.and.j.eq.2) then
       write(10,10) r,g,b
   10  format(f6.2,1x,f6.2,1x,f6.2)
       else
       write(10,20) r,g,b
   20  format(f6.2,1x,f6.2,1x,f6.2,",")
       end if
c
  100  continue
       return
       end
c------------------------------------------------------
       subroutine coordbe
c------------------------------------------------------
       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
c------------------------------------------------------
       subroutine coorden1(nx)
c------------------------------------------------------
       dimension ia(9000)
       write(10,21)
   21  format("]")
       write(10,22)
   22  format("}")
       write(10,23)
   23  format("colorPerVertex TRUE")
       write(10,24)
   24  format("coordIndex [")
c
       nx2=nx*2
       do 100 i=1,nx2
       i1=1+4*(i-1)
       i2=2+4*(i-1)
       i3=3+4*(i-1)
       i4=4+4*(i-1)
       ia(i1)=0+i-1
       ia(i2)=1+i-1
       ia(i3)=2+i-1
       ia(i4)=-1
  100  continue
c
       jj=4*(nx-1)*2-1
       write(10,102) (ia(i),i=1,jj)
  102  format(8(3i2,i3),2i2,2i3,i2,3i3,88i3,2i3,i4,i3,i3,2i4,i3,
     1        398(3i4,i3),/500(3i4,i3),/2i4,i5,i3,i4,2i5,i3,
     2        /500(3i5,i3),/500(3i5,i3),/500(3i5,i3),/500(3i5,i3))
       return
       end
c-----------------------------------------------------
       subroutine ifaceen
c-----------------------------------------------------
       write(10,20)
   20  format("]")
       write(10,21)
   21  format("solid FALSE",/"}")
       write(10,22)
   22  format("appearance Appearance {")
       write(10,23)
   23  format("material Material {")
       write(10,24)
   24  format("transparency 0.000000",/"}",/"}")
       write(10,25)
   25  format("}")
       return
       end
c-------------------------------------------------------
       subroutine coordp2(nx,ny,p)
c-------------------------------------------------------
       dimension p(1)
c
       do 100 i=1,nx
       do 100 j=1,2
       j1=j+ny-1
c
       i1=3*(i-1)+1+3*nx*(j1-1)
       i2=i1+1
       i3=i1+2
       x=p(i1)
       y=p(i2)
       z=p(i3)
c
       if(i.eq.nx.and.j.eq.2) then
       write(10,10) x,y,z
   10  format(f9.2,1x,f9.2,1x,f9.2)
       else
       write(10,20) x,y,z
   20  format(f9.2,1x,f9.2,1x,f9.2,",")
       end if
c
  100  continue
       return
       end
c---------------------------------------------------------
       subroutine arrow2(nx,ny,u,p,rd1,rd2,h,cx,cy,cz,br,bg,bb,cr,cg,cb)
c--------------------------------------------------------
       dimension p(1),u(1)
       ny1=ny-1
c
       do 10 j=1,ny1
       call cylinderb(nx,j,u,p,rd1,br,bg,bb) 
   10  continue
c
       jj=ny1
       i=nx
       j=2
       j1=j+jj-1
       i1=3*(i-1)+1+3*nx*(j1-1)
       i2=i1+1
       i3=i1+2
       cx=p(i1)+0.8*rd2
       cy=p(i2)
       cz=p(i3)+rd1
c
       x1=0.0
       y1=0.0
       z1=-1.0
       rr=1.57
c
       call cone1(rd2,h,cr,cg,cb,cx,cy,cz,x1,y1,z1,rr)
       return
       end
c------------------------------------------------------
       subroutine cylinderb(nx,ny,u,p,rd,br,bg,bb)
c          * gradation color *
c------------------------------------------------------
       dimension p(1),u(1)
c
       call facebe
       call colorp3(nx,ny,u,br,bg,bb)
       call coordbe
       call coordp3(nx,ny,p,rd)
       call coorden2(nx)
       call ifaceen
c
       return
       end
c------------------------------------------------------
       subroutine cone1(rd,h,r,g,b,x,y,z,x1,y1,z1,rr)
c------------------------------------------------------
       write(10,20) x,y,z
   20  format("Transform {",/"translation ",f7.2,f7.2,f7.2)
       write(10,21) x1,y1,z1,rr
   21  format("rotation ",f7.2,f7.2,f7.2,f7.2)
       write(10,22)
   22  format("children [")
       write(10,23)
   23  format("Shape {")
       write(10,24)
   24  format("appearance Appearance {")
       write(10,25)
   25  format("material Material {")
       write(10,26) r,g,b
   26  format("diffuseColor",f5.2,f5.2,f5.2)
       write(10,27)
   27  format("}",/"}")
       write(10,28)
   28  format("geometry Cone {")
       write(10,29) rd
   29  format("bottomRadius  ",f5.2)
       write(10,30) h
   30  format("height  ",f5.2)
       write(10,31)
   31  format("}",/"}",/"]",/"}")
c
       return
       end
c--------------------------------------------------------
       subroutine colorp3(nx,ny,u,br,bg,bb)
c           * bottom tuki *
c--------------------------------------------------------
       dimension u(1)
c
       do 100 i=1,nx
       do 100 j=1,2
       j1=j+ny-1
c
       i1=3*(i-1)+1+3*nx*(j1-1)
       i2=i1+1
       i3=i1+2
       r=u(i1)
       g=u(i2)
       b=u(i3)
c
       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
       write(10,10) r,g,b
   10  format(f8.2,1x,f8.2,1x,f8.2,",")
c
  100  continue
c----( bottom chusin color )
       write(10,120) br,bg,bb
  120  format(f8.2,1x,f8.2,1x,f8.2)
c
       return
       end
c--------------------------------------------------
       subroutine coordp3(nx,ny,p,rd)
c          * bottom tuki *
c--------------------------------------------------
       dimension p(1)
c
       do 100 i=1,nx
       do 100 j=1,2
       j1=j+ny-1
c
       i1=3*(i-1)+1+3*nx*(j1-1)
       i2=i1+1
       i3=i1+2
       x=p(i1)
       y=p(i2)
       z=p(i3)
c
       write(10,20) x,y,z
   20  format(f9.2,1x,f9.2,1x,f9.2,",")
c
  100  continue
c
c----( bottom no chusin )
c
       i1=3*nx*2+1
       i2=i1+1
       i3=i1+2
       x=p(i1)
       y=p(i2)
       z=p(i3)
c
       write(10,110) x,y,z
  110  format(f9.2,1x,f9.2,1x,f9.2)
c
       return
       end
c------------------------------------------------------
       subroutine coorden2(nx)
c          * bottom tuki *
c------------------------------------------------------
       dimension ia(9000)
       write(10,21)
   21  format("]")
       write(10,22)
   22  format("}")
       write(10,23)
   23  format("colorPerVertex TRUE")
       write(10,24)
   24  format("coordIndex [")
c
       nx2=nx*2
       do 100 i=1,nx2
       i1=1+4*(i-1)
       i2=2+4*(i-1)
       i3=3+4*(i-1)
       i4=4+4*(i-1)
       ia(i1)=0+i-1
       ia(i2)=1+i-1
       ia(i3)=2+i-1
       ia(i4)=-1
  100  continue
c
       jj=4*(nx-1)*2
       write(10,102) (ia(i),i=1,jj)
  102  format(8(3i2,i3),2i2,2i3,i2,3i3,88i3,2i3,i4,i3,i3,2i4,i3,
     1        398(3i4,i3),/500(3i4,i3),/2i4,i5,i3,i4,2i5,i3,
     2        /500(3i5,i3),/500(3i5,i3),/500(3i5,i3),/500(3i5,i3))
c
       write(10,103) ia(jj-2),ia(jj-1),ia(1)
  103  format(3i5," -1")
       write(10,104) ia(jj-1),ia(1),ia(2)
  104  format(3i5," -1")
c
c---- bottom
       nx2=nx*2
       do 200 i=1,nx
       i1=1+4*(i-1)
       i2=2+4*(i-1)
       i3=3+4*(i-1)
       i4=4+4*(i-1)
       ia(i1)=0+2*(i-1)
       ia(i2)=ia(i1)+2
       ia(i3)=nx2
       ia(i4)=-1
  200  continue
c
       jj=4*(nx-1)
       write(10,202) (ia(i),i=1,jj)
  202  format(500(3i4,i3),/500(3i4,i3),/500(3i5,i3),/500(3i5,i3))
c
       write(10,203) ia(jj-2),ia(1),nx2
  203  format(3i4," -1")
c
c---- cap
       nx2=nx*2
       do 300 i=1,nx
       i1=1+4*(i-1)
       i2=2+4*(i-1)
       i3=3+4*(i-1)
       i4=4+4*(i-1)
       ia(i1)=1+2*(i-1)
       ia(i2)=ia(i1)+2
       ia(i3)=nx2
       ia(i4)=-1
  300  continue
c
       jj=4*(nx-1)
       write(10,302) (ia(i),i=1,jj)
  302  format(500(3i4,i3),/500(3i4,i3),/500(3i5,i3),/500(3i5,i3))
c
       write(10,303) ia(jj-2),ia(1),nx2
  303  format(3i4," -1")
c
       return
       end
c-------------------------------------------------------
       subroutine linemesh(p,n1,rd,r,g,b)
c-------------------------------------------------------
       dimension p(1)
       dimension p1(3*2)
       dimension q(3*12*2)
c
       n11=n1-2
       do 20 i=1,n11
       do 10 j=1,2
c
       i1=3*(i-1)+1+3*(j-1)
       i2=i1+1
       i3=i1+2
c---<< chusin ten set >>
       p1(1+(j-1)*3)=p(i1)
       p1(2+(j-1)*3)=p(i2)
       p1(3+(j-1)*3)=p(i3)
c
       xa=abs(p(i1)-p(i1+3))
       ya=abs(p(i2)-p(i2+3))
       za=abs(p(i3)-p(i3+3))
       amx=amax1(xa,ya,za)  
c
       x=p(i1)
       y=p(i2)
       z=p(i3)
c
       if(amx.eq.xa)  then
               k=1
        else
          if(amx.eq.ya)  then
                 k=2
           else
             if(amx.eq.za)  then
                   k=3
             end if
          end if
       end if
c
c---<< ensyu no point wo motomeru >>
c
       call ensyu(x,y,z,rd,q,j,k)
   10  continue
c
       call linemeshb(12,1,q,rd,r,g,b,k,p1)
   20  continue
c---<< cone >>
       rr=1.57
       if(k.eq.1)   then
            x1=0.0
            y1=0.0
            z1=-1.0
       end if
       if(k.eq.3)   then
            x1=-1.0
            y1=0.0
            z1=0.0
       end if
c
       rd2=rd*2
       h=rd*3
       x=p1(4)
       y=p1(5)
       z=p1(6)
       call cone1(rd2,h,r,g,b,x,y,z,x1,y1,z1,rr)
c
       return
       end
c------------------------------------------------------
       subroutine ensyu(x,y,z,rd,q,j,k)
c------------------------------------------------------
       dimension q(1)
c
       x1=rd/2*1.73205
       y1=rd/2
       j2=36*(j-1)
c--[ max=xa ]
       if(k.eq.1)    then
c (x)
       do 10 i=1,12
       q(1+(i-1)*3+j2)=x
   10  continue
c (y)
       q(2+j2)=y+rd
       q(5+j2)=y+x1
       q(8+j2)=y+y1
       q(11+j2)=y
       q(14+j2)=y-y1
       q(17+j2)=y-x1
       q(20+j2)=y-rd
       q(23+j2)=y-x1
       q(26+j2)=y-y1
       q(29+j2)=y
       q(32+j2)=y+y1
       q(35+j2)=y+x1
c (z)
       q(3+j2)=z
       q(6+j2)=z+y1
       q(9+j2)=z+x1
       q(12+j2)=z+rd
       q(15+j2)=z+x1
       q(18+j2)=z+y1
       q(21+j2)=z
       q(24+j2)=z-y1
       q(27+j2)=z-x1
       q(30+j2)=z-rd
       q(33+j2)=z-x1
       q(36+j2)=z-y1
       go to 99
       end if
c--[ max=ya ]
       if(k.eq.2)    then
c (x)
       q(1+j2)=x
       q(4+j2)=x+y1
       q(7+j2)=x+x1
       q(10+j2)=x+rd
       q(13+j2)=x+x1
       q(16+j2)=x+y1
       q(19+j2)=x
       q(22+j2)=x-y1
       q(25+j2)=x-x1
       q(28+j2)=x-rd
       q(31+j2)=x-x1
       q(34+j2)=x-y1
c (y)
       do 20 i=1,12
       q(2+(i-1)*3+j2)=y
   20  continue
c (z)
       q(2+j2)=z+rd
       q(5+j2)=z+x1
       q(8+j2)=z+y1
       q(11+j2)=z
       q(14+j2)=z-y1
       q(17+j2)=z-x1
       q(20+j2)=z-rd
       q(23+j2)=z-x1
       q(26+j2)=z-y1
       q(29+j2)=z
       q(32+j2)=z+y1
       q(35+j2)=z+x1
       go to 99
       end if
c--[ max=za ]
       if(k.eq.3)    then
c (x)
       q(1+j2)=x
       q(4+j2)=x+y1
       q(7+j2)=x+x1
       q(10+j2)=x+rd
       q(13+j2)=x+x1
       q(16+j2)=x+y1
       q(19+j2)=x
       q(22+j2)=x-y1
       q(25+j2)=x-x1
       q(28+j2)=x-rd
       q(31+j2)=x-x1
       q(34+j2)=x-y1
c (y)
       q(2+j2)=y+rd
       q(5+j2)=y+x1
       q(8+j2)=y+y1
       q(11+j2)=y
       q(14+j2)=y-y1
       q(17+j2)=y-x1
       q(20+j2)=y-rd
       q(23+j2)=y-x1
       q(26+j2)=y-y1
       q(29+j2)=y
       q(32+j2)=y+y1
       q(35+j2)=y+x1
c (z)
       do 30 i=1,12
       q(3*i+j2)=z
   30  continue
       end if
c
   99  continue
       return
       end
c------------------------------------------------------
       subroutine linemeshb(nx,ny,q,rd,r,g,b,k,p1)
c          * 1 color *    
c------------------------------------------------------
       dimension q(1)
c
       call facebe
c
c----( color )
       write(10,10) r,g,b
   10  format(f9.2,1x,f9.2,1x,f9.2)
c
       call coordbe1(nx)
       call coordp4(nx,ny,q,rd,k,p1)
       call coorden3(nx)
       call ifaceen
c
       return
       end
c------------------------------------------------------
       subroutine coordbe1(nx)
c------------------------------------------------------
       write(10,21)
   21  format("]")
       write(10,22)
   22  format("}")
c----( colorindexed )
       nx4=nx*4
       write(10,20)
   20  format("colorIndex [ ")
       write(10,30) (0,i=1,nx4)
   30  format(5000i3)
       write(10,40)
   40  format(" ]")
c
       write(10,23)
   23  format("coord Coordinate {")
       write(10,24)
   24  format("point [")
       return
       end
c------------------------------------------------------
       subroutine coordp4(nx,ny,q,rd,k,p1)
c          * bottom tuki *
c------------------------------------------------------
       dimension q(1)
       dimension p1(1)
c
       do 100 i=1,nx
       do 100 j=1,2
c
       i1=3*(i-1)+1+3*nx*(j-1)
       i2=i1+1
       i3=i1+2
       x=q(i1)
       y=q(i2)
       z=q(i3)
c
       write(10,20) x,y,z
   20  format(f9.2,1x,f9.2,1x,f9.2,",")
c
  100  continue
c
c----( bottom no chusin )
       x=p1(1)
       y=p1(2)
       z=p1(3)
c
       write(10,30) x,y,z
   30  format(f9.2,1x,f9.2,1x,f9.2,",")
c
c----( cap no chusin )
       x=p1(4)
       y=p1(5)
       z=p1(6)
c
       write(10,40) x,y,z
   40  format(f9.2,1x,f9.2,1x,f9.2)
c
       return
       end
c--------------------------------------------------
       subroutine coorden3(nx)
c          * bottom tuki * * 1 color     *
c--------------------------------------------------
       dimension ia(9000)
       write(10,21)
   21  format("]")
       write(10,22)
   22  format("}")
       write(10,23)
   23  format("colorPerVertex FALSE")
       write(10,24)
   24  format("coordIndex [")
c
       nx2=nx*2
       do 100 i=1,nx2
       i1=1+4*(i-1)
       i2=2+4*(i-1)
       i3=3+4*(i-1)
       i4=4+4*(i-1)
       ia(i1)=0+i-1
       ia(i2)=1+i-1
       ia(i3)=2+i-1
       ia(i4)=-1
  100  continue
c
       jj=4*(nx-1)*2
       write(10,102) (ia(i),i=1,jj)
  102  format(8(3i2,i3),2i2,2i3,i2,3i3,88i3,2i3,i4,i3,i3,2i4,i3,
     1        398(3i4,i3),/500(3i4,i3),/2i4,i5,i3,i4,2i5,i3,
     2        /500(3i5,i3),/500(3i5,i3),/500(3i5,i3),/500(3i5,i3))
c
       write(10,103) ia(jj-2),ia(jj-1),ia(1)
  103  format(3i5," -1")
       write(10,104) ia(jj-1),ia(1),ia(2)
  104  format(3i5," -1")
c
c---- bottom
       nx2=nx*2
       do 200 i=1,nx
       i1=1+4*(i-1)
       i2=2+4*(i-1)
       i3=3+4*(i-1)
       i4=4+4*(i-1)
       ia(i1)=0+2*(i-1)
       ia(i2)=ia(i1)+2
       ia(i3)=nx2
       ia(i4)=-1
  200  continue
c
       jj=4*(nx-1)
       write(10,202) (ia(i),i=1,jj)
  202  format(500(3i4,i3),/500(3i4,i3),/500(3i5,i3),/500(3i5,i3))
c
       write(10,203) ia(jj-2),ia(1),nx2
  203  format(3i4," -1")
c
c---- cap
       nx21=nx*2+1
       do 300 i=1,nx
       i1=1+4*(i-1)
       i2=2+4*(i-1)
       i3=3+4*(i-1)
       i4=4+4*(i-1)
       ia(i1)=1+2*(i-1)
       ia(i2)=ia(i1)+2
       ia(i3)=nx21
       ia(i4)=-1
  300  continue
c
       jj=4*(nx-1)
       write(10,302) (ia(i),i=1,jj)
  302  format(500(3i4,i3),/500(3i4,i3),/500(3i5,i3),/500(3i5,i3))
c
       write(10,303) ia(jj-2),ia(1),nx21
  303  format(3i4," -1")
c
       return
       end
c-------------------------------------------------------
       subroutine backd(n1,n2,rd)
c-------------------------------------------------------
       write(10,21)
   21  format("Background {",/"skyColor [ ")
       call colors1(n2)
       write(10,22)
   22  format("skyAngle [ ")
       call angle1(rd,n1)
c
       write(10,31)
   31  format("groundColor [ ")
       call colorg1(n1)
       write(10,32)
   32  format("groundAngle [ ")
       call angle1(rd,n1)
c
       write(10,99)
   99  format("}")
       return
       end
c-------------------------------------------------------
       subroutine colors1(n2)
c------------------------------------------------------
       r=0.0
       do 10 i=1,n2
       g=0.1*float(i)
       b=0.1*float(i)
c
       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
       call colorp(r,g,b)
   10  continue
c
       r=0.0
       do 20 i=1,n2
       g=0.1*float(n2+i)
       b=0.1*float(n2+i)
c
       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 
       call colorp(r,g,b)
   20  continue
c
       write(10,30)
   30  format(" ]")
       return
       end
c----------------------------------------------------
       subroutine colorg1(n1)
c----------------------------------------------------
       do 10 i=1,n1
       r=0.0
       b=0.0
       g=0.1*float(i)
       call colorp(r,g,b)
   10  continue
c
       write(10,30)
   30  format(" 0.20 1.00 0.70")
       write(10,40)
   40  format(" ]")
       return
       end
c-------------------------------------------------------
       subroutine angle1(rd,n1)
c-------------------------------------------------------
       do 30 i=1,n1
       ang=rd*i
       if (i.eq.n1) then
       write(10,22) ang
   22  format(f7.4," ]")
       else
       write(10,23) ang
   23  format(f7.4,",")
       end if
   30  continue
       return
       end
c------------------------------------------------------
       subroutine colorp(r,g,b)
c------------------------------------------------------
       write(10,20) r,g,b
   20  format(f5.2,f5.2,f5.2)
       return
       end
c-------------------------------------------------------
       subroutine triangbk1(nx,ny,u,p)
c-------------------------------------------------------
       dimension p(1),u(1)
c
       call facebe
       call colorp1k2(ii,nx,ny,u)
       call coordbe
       call coordp1(ii,nx,ny,p)
       call coorden1(nx)
       call ifaceen
c
       return
       end
c-------------------------------------------------------
       subroutine colorp1k2(ii,nx,ny,u)
c-------------------------------------------------------
       dimension u(1)
c
       do 100 i=1,nx
       do 100 j=1,2
c
       i1=3*(i-1)+1+3*nx*(j-1)
       i2=i1+1
       i3=i1+2
       r=u(i1)
       g=u(i2)
       b=u(i3)
c
       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.eq.nx.and.j.eq.2) then
       write(10,10) r,g,b
   10  format(f6.2,1x,f6.2,1x,f6.2)
       else
       write(10,20) r,g,b
   20  format(f6.2,1x,f6.2,1x,f6.2,",")
       end if
c
  100  continue
       return
       end
c-------------------------------------------------------
       subroutine coordp1(ii,nx,ny,p)
c-------------------------------------------------------
       dimension p(1)
c
       do 100 i=1,nx
       do 100 j=1,2
c
       i1=3*(i-1)+1+3*nx*(j-1)
       i2=i1+1
       i3=i1+2
       x=p(i1)
       y=p(i2)
       z=p(i3)
c
       if(i.eq.nx.and.j.eq.2) then
       write(10,10) x,y,z
   10  format(f6.2,1x,f6.2,1x,f6.2)
       else
       write(10,20) x,y,z
   20  format(f6.2,1x,f6.2,1x,f6.2,",")
       end if
c
  100  continue
       return
       end
c----------------------------------------------------
       subroutine faceen
c----------------------------------------------------
       write(10,21)
   21  format("]")
       write(10,22)
   22  format("solid FALSE")
       write(10,23)
       write(10,23)
   23  format("}")
       return
       end
c--------------------------------------------------
       subroutine defusem(nx,ny,u,p,x1,y1,z1,x2,y2,z2,r,n)
c--------------------------------------------------
       dimension p(1),u(1)
       ny1=ny-1
c
       do 10 j=1,ny1
       call defusem1(nx,j,u,p,x1,y1,z1,x2,y2,z2,r,n)
   10  continue
c
       return
       end
c--------------------------------------------------
       subroutine defusem1(nx,ny,u,p,x1,y1,z1,x2,y2,z2,r,n)
c--------------------------------------------------
       dimension p(1),u(1)
c
       call defface
       call colorp2(nx,ny,u)
       call coordbe
       call coordp2(nx,ny,p)
       call coorden1(nx)
       call ifaceen
c
       do 10 i=1,n
       call useface(x1,y1,z1,x2,y2,z2,r,i)
   10  continue
c
       return
       end
c--------------------------------------------------
       subroutine defface
c--------------------------------------------------
       write(10,21)
   21  format("DEF Face Shape {")
       write(10,22)
   22  format("geometry IndexedFaceSet {")
       write(10,23)
   23  format("color Color {")
       write(10,24)
   24  format("color [")
       return
       end
c-------------------------------------------------------
       subroutine useface(x1,y1,z1,x2,y2,z2,r,i)
c--------------------------------------------------------
       write(10,21)
   21  format("Transform {")
       write(10,22) x1,y1,z1,r
   22  format("rotation ",f6.2,f6.2,f6.2,f6.2)
       write(10,23) x2,y2*i,z2
   23  format("translation ",f9.2,f9.2,f9.2)
       write(10,24)
   24  format("children USE Face",/"}")
       return
       end
c-------------------------------------------------------
       subroutine triangm1(nx,ny,u,p)
c-------------------------------------------------------
       dimension p(1),u(1)
       ny1=ny-1
c
       call facebe
c
       do 10 j=1,ny1
       call colorp2(nx,j,u)
   10  continue
c
       call coordbe
c
       do 20 j=1,ny1
       call coordp2(nx,j,p)
   20  continue
c
       call coorden5(nx,ny1)
c
       call ifaceen 
c
       return
       end
c------------------------------------------------------
       subroutine coorden5(nx,ny1)
c------------------------------------------------------
       dimension ia(9000)
       write(10,21)
   21  format("]")
       write(10,22)
   22  format("}")
       write(10,23)
   23  format("colorPerVertex TRUE")
       write(10,24)
   24  format("coordIndex [")
c
cc      nx2=(nx-1)*2
       nx2=nx*2
       do 200 j=1,ny1
c
       do 100 i=1,nx2
       i1=1+4*(i-1)
       i2=2+4*(i-1)
       i3=3+4*(i-1)
       i4=4+4*(i-1)
       ia(i1)=(j-1)*2*nx+0+i-1
       ia(i2)=(j-1)*2*nx+1+i-1
       ia(i3)=(j-1)*2*nx+2+i-1
       ia(i4)=-1
  100  continue
c
       jj=4*(nx-1)*2-1
       write(10,102) (ia(i),i=1,jj)
c 102  format(8(3i2,i3),2i2,2i3,i2,3i3,88i3,2i3,i4,i3,i3,2i4,i3,
c    1        398(3i4,i3),/500(3i4,i3),/2i4,i5,i3,i4,2i5,i3,
  102  format(250(3i7,i3),/250(3i7,i3),/250(3i7,i3),/250(3i7,i3),
     1        /250(3i7,i3),/250(3i7,i3),/250(3i7,i3),/250(3i7,i3),
     2        /250(3i7,i3),/250(3i7,i3),/250(3i7,i3),/250(3i7,i3),
     3        /250(3i7,i3),/250(3i7,i3),/250(3i7,i3),/250(3i7,i3),
     4        /250(3i7,i3),/250(3i7,i3),/250(3i7,i3),/250(3i7,i3),
     5        /250(3i7,i3),/250(3i7,i3),/250(3i7,i3),/250(3i7,i3),
     6        /250(3i7,i3),/250(3i7,i3),/250(3i7,i3),/250(3i7,i3),
     7        /250(3i7,i3),/250(3i7,i3),/250(3i7,i3),/250(3i7,i3),
     8        /250(3i7,i3),/250(3i7,i3),/250(3i7,i3),/250(3i7,i3))
c
  200  continue
       return
       end
c-------------------------------------------------------
       subroutine imagebe(chr,n)
c-------------------------------------------------------
       character chr*15,ica*15,ich(15)*1
       equivalence (ica,ich)
       ica=chr
       write(10,21)
   21  format("Shape {")
       write(10,22)
   22  format("appearance Appearance {")
       write(10,23)
   23  format("texture ImageTexture {")
   24  write(10,*) "url",'"',(ich(i),i=1,n),'"'
       write(10,25)
   25  format("}",/"}")
       write(10,26)
   26  format("geometry IndexedFaceSet {")
       write(10,27)
   27  format("coord Coordinate {")
       write(10,28)
   28  format("point [")
       return
       end
c-------------------------------
       subroutine imageen(ln)
c-------------------------------
       dimension ia(2000)
       write(10,21)
   21  format("]")
       write(10,22)
   22  format("}")
       write(10,23)
   23  format("colorPerVertex FALSE")
       write(10,24)
   24  format("coordIndex [")
c
       jj=ln*2
       do 30 i=1,jj
       i1=1+4*(i-1)
       i2=2+4*(i-1)
       i3=3+4*(i-1)
       i4=4+4*(i-1)
       ia(i1)=0+i-1
       ia(i2)=1+i-1
       ia(i3)=2+i-1
       ia(i4)=-1
   30  continue
c
       jj=4*(ln-1)*2-1
       write(10,31) (ia(i),i=1,jj)
   31  format(8(3i2,i3),2i2,2i3,i2,3i3,88i3,2i3,i4,i3,i3,2i4,i3,
     1        898(3i4,i3),2i4,2i5,i3,i4,2i5,i3,2000(3i5,i3))
c
       write(10,32)
   32  format("]")
       write(10,33)
   33  format("solid FALSE")
       write(10,34)
   34  format("}",/"}")
       return
       end
