program Koch_snowflake
      implicit none
      integer pgopen,n,big_small,n_a,n_b
      real mhalf_sqrt3
      parameter ( mhalf_sqrt3 = (-1)*sqrt(3.0)/2 )
      character*10 string
      write(*,*) 'Big window or small (2 = big ; 1 = small) :'
      read (*,*) big_small
c     if (big_small.eq.2) big_small = big_small*0.75 
      write(*,*) 'How many levels to draw ?'
      read (*,*) n
c     write(string,'(i2)') n   ! This will also works.
      write(string,100) n
 100  format (i2)
      string = 'n ='//string

      if (pgopen('/xwin').le.0) stop
      call pgpap(6.0*big_small,1.1)
      call pgenv(0.0,1.0,-0.9,0.3,1,-2)
c     call pgsch(2.0)
      call pglab(string,'','Koch Snowflake')

      n_a = n
      n_b = n
      call draw_koch( (0.0,0.0), (1.0,0.0), n)
      call draw_koch( (1.0,0.0), (0.5,mhalf_sqrt3), n_a)
      call draw_koch( (0.5,mhalf_sqrt3), (0.0,0.0), n_b)
      

      call pgclos

      end


      subroutine draw_koch (z1,z2,n)
      implicit none
      integer n,n1,n2,n3,n4
      complex z1,z2,za,zb,zc,dz,phase
      real length,pi_3,x(2),y(2)
      parameter (pi_3=3.14159268/3.0) 
      if (n.gt.0) then
        n = n - 1
        n1 = n
        n2 = n
        n3 = n
        n4 = n
c       write(*,*) 'The new n is : ',n
        dz = (z2-z1)/3
        za = z1 + dz
        zc = z1 + 2*dz
        phase = (1.0,0.0)*cos(pi_3) + (0.0,1.0)*sin(pi_3)
        zb = za + dz*phase
        call draw_koch (z1,za,n1)
        call draw_koch (za,zb,n2)
        call draw_koch (zb,zc,n3)
        call draw_koch (zc,z2,n4)
      else
        x(1) = real(z1)
        x(2) = real(z2)
        y(1) = imag(z1)
        y(2) = imag(z2)
        call pgline (2,x,y)
c       write(*,*) 'A line has been drawn between', z1,' and ',z2
      endif
      end