Subversion Repositories lagranto.ecmwf

Rev

Rev 3 | Show entire file | Ignore whitespace | Details | Blame | Last modification | View Log | RSS feed

Rev 3 Rev 44
Line 734... Line 734...
734 
      real      pi180
 734 
      real      pi180
735 
      parameter (pi180=3.14159/180.)
 735 
      parameter (pi180=3.14159/180.)
736 
 
 736 
 
737 
c     Auixiliary variables
 737 
c     Auixiliary variables
738 
      real     dx,dy
 738 
      real     dx,dy
- 
 
 739 
      real     dlon
739 
 
 740 
 
- 
 
 741 
      dlon = lon1-lon0
- 
 
 742 
      if ( dlon.gt.180.  ) dlon = dlon - 360.
- 
 
 743 
      if ( dlon.lt.-180. ) dlon = dlon + 360.
- 
 
 744
740 
      dx = (lon1-lon0) * cos(pi180*0.5*(lat0+lat1))
 745 
      dx = dlon * cos(pi180*0.5*(lat0+lat1))
741 
      dy = lat1-lat0
 746 
      dy = lat1-lat0
742 
 
 747 
 
743 
      call getangle(1.,0.,dx,dy,head)
 748 
      call getangle(1.,0.,dx,dy,head)
744 
 
 749 
 
745 
      end
 750 
      end
- 
 
 751
- 
 
 752 
c     -------------------------------------------------------------------------
- 
 
 753 
c     Directional change of the trajectory (HEAD)
- 
 
 754 
c     -------------------------------------------------------------------------
- 
 
 755 
 
- 
 
 756 
      subroutine calc_DANGLE (dangle,lon0,lat0,lon1,lat1,lon2,lat2)
- 
 
 757
- 
 
 758 
      implicit none
- 
 
 759 
 
- 
 
 760 
c     Argument declaration
- 
 
 761 
      real  dangle     ! Directiopanl change
- 
 
 762 
      real  lon0,lat0  ! t-1
- 
 
 763 
      real  lon1,lat1  ! t
- 
 
 764 
      real  lon2,lat2  ! t+1
- 
 
 765
- 
 
 766 
c     Physical parameters
- 
 
 767 
      real      pi180
- 
 
 768 
      parameter (pi180=3.14159/180.)
- 
 
 769 
      real      eps
- 
 
 770 
      parameter (eps=0.00001)
- 
 
 771 
 
- 
 
 772 
c     Auixiliary variables
- 
 
 773 
      real     dx1,dy1,dx2,dy2,norm,cross,dlon1,dlon2
- 
 
 774
- 
 
 775 
      dlon1 = lon1 - lon0
- 
 
 776 
      if ( dlon1.gt.180.  ) dlon1 = dlon1 - 360.
- 
 
 777 
      if ( dlon1.lt.-180. ) dlon1 = dlon1 + 360.
- 
 
 778 
      dlon2 = lon2 - lon1
- 
 
 779 
      if ( dlon2.gt.180.  ) dlon2 = dlon2 - 360.
- 
 
 780 
      if ( dlon2.lt.-180. ) dlon2 = dlon2 + 360.
- 
 
 781
- 
 
 782 
      dx1 = dlon1 * cos(pi180*0.5*(lat1+lat0))
- 
 
 783 
      dy1 = lat1 - lat0
- 
 
 784 
      dx2 = dlon2 * cos(pi180*0.5*(lat2+lat1))
- 
 
 785 
      dy2 = lat2 - lat1
- 
 
 786
- 
 
 787 
      norm = sqrt( (dx1**2 + dy1**2) * (dx2**2 + dy2**2) )
- 
 
 788 
 
- 
 
 789 
      if ( norm.gt.eps ) then
- 
 
 790 
         cross = ( dx1 * dy2 - dy1 * dx2 ) / norm
- 
 
 791 
         if ( cross.ge.1. ) then
- 
 
 792 
            dangle = 90.
- 
 
 793 
         elseif (cross.le.-1.) then
- 
 
 794 
            dangle = -90.
- 
 
 795 
         else
- 
 
 796 
            dangle = 1./pi180 * asin( cross )
- 
 
 797 
         endif
- 
 
 798 
      else
- 
 
 799 
         dangle = -999.
- 
 
 800 
      endif
- 
 
 801 
 
- 
 
 802 
      end
746 
 
 803 
 
747 
c
804 
c
748 
c     *************************************************************************
 805 
c     *************************************************************************
749 
c     Auxiliary subroutines and functions
 806 
c     Auxiliary subroutines and functions
750 
c     *************************************************************************
 807 
c     *************************************************************************