program fall_polt_drag_dagger_2
      real k1,m
      real dgx(4),dgy(4),dgx0(4),dgy0(4),dgxold(4),dgyold(4)

      dgx0(1)=0.0
      dgy0(1)=0.1
      dgx0(2)=0.6
      dgy0(2)=0.0
      dgx0(3)=0.0
      dgy0(3)=-0.1
      dgx0(4)=dgx0(1)
      dgy0(4)=dgy0(1)


      write(*,*)'initial position ?'
      read(*,*)x
      xini=x
      write(*,*)'initial hight ?'
      read(*,*)y
      write(*,*)'Vx ?'
      read(*,*)vx
      write(*,*)'Vy ?'
      read(*,*)vy
      write(*,*)'air drag coefficient?'
      read(*,*)k1
      write(*,*)'Mass of the object?'
      read(*,*)m

      if(pgopen('/xwin')<=0)stop
      call pgenv(0.0,10.0,0.0,10.0,1,0)

      v=sqrt(vx**2+vy**2)
      do i=1,4
        dgx(i)=dgx0(i)*vx/v + dgy0(i)*(-1)*vy/v
        dgy(i)=dgx0(i)*vy/v + dgy0(i)*vx/v
      end do

      do i=1,4
        dgx(i)=dgx(i)+x
        dgy(i)=dgy(i)+y
      end do

      do i=1,4
        dgxold(i)=dgx0(i)
        dgyold(i)=dgy0(i)
      end do

      dt=0.00001
      do while(y>0)
        
        v=sqrt(vx**2+vy**2)
        do i=1,4
          dgx(i)=dgx0(i)*vx/v + dgy0(i)*(-1)*vy/v
          dgy(i)=dgx0(i)*vy/v + dgy0(i)*vx/v
        end do        

        f_d=(-1)*k1*v
        f_dx=f_d*vx/v
        f_dy=f_d*vy/v
        ax=f_dx/m
        ay=-9.8+f_dy/m
        vx=vx+ax*dt
        vy=vy+ay*dt
        x=x+vx*dt
        y=y+vy*dt
        do i=1,4
          dgx(i)=dgx(i)+x
          dgy(i)=dgy(i)+y
        end do
        call pgbbuf
        call pgsci(0)
        call pgline(4,dgxold,dgyold)
        call pgsci(1)
        call pgline(4,dgx,dgy)
        call pgebuf
        do i=1,4
          dgxold(i)=dgx(i)
          dgyold(i)=dgy(i)
        end do
      end do

      write(*,*)'Final hight =',y,'meters.'
      write(*,*)'Final distace =',x-xini,'meters.'

      call pgclos

      end