!+
PROGRAM GetMeanAerodynamicChord
! ---------------------------------------------------------------------------
! PURPOSE - Compute the mean aerodynamic chord of a wing.
! AUTHOR - Ralph L. Carmichael, Public Domain Aeronautical Software
! REVISION HISTORY
!   DATE  VERS PERSON  STATEMENT OF CHANGES
! 28Feb92  0.5   RLC   Original coding (Fortran 77) (from old fragments)
! 22Oct00  0.6   RLC   Recoded

IMPLICIT NONE

  TYPE::WINGCHORD
    REAL:: y,xle,c
  END TYPE WINGCHORD

  CHARACTER(LEN=*),PARAMETER:: HEADER = "Mean Aerodynamic Chord Calculator"
  CHARACTER(LEN=*),PARAMETER:: AUTHOR = &
    "Ralph Carmichael, Public Domain Aeronautical Software"
  CHARACTER(LEN=*),PARAMETER:: VERSION = "0.5 (28 Oct 2000)"

  TYPE(WINGCHORD),DIMENSION(100):: ch
  INTEGER:: errCode
  CHARACTER(LEN=80):: fileName
  INTEGER,PARAMETER:: IDAT=1
  INTEGER:: k,n
  REAL:: area,mac,ymac,xlemac
  REAL:: y,xle,c

  NAMELIST /CHORD/ y,xle,c
!----------------------------------------------------------------------------

  WRITE(*,*) HEADER
  WRITE(*,*) AUTHOR
  WRITE(*,*) VERSION

  DO
    WRITE(*,*) "Enter the name of the input file: "
    READ(*,'(A)') fileName
    IF (fileName == ' ') STOP
    OPEN(UNIT=IDAT, FILE=FILENAME, IOSTAT=errCode, &
      STATUS='OLD', ACTION='READ', POSITION='REWIND')
    IF (errCode == 0) THEN
      EXIT
    ELSE
      WRITE(*,*) "Unable to open this file. Try again."
    END IF
  END DO

  n=0
  DO
    READ(IDAT,CHORD,IOSTAT=errCode)
    IF (errCode /= 0) EXIT
    n=n+1
    ch(n)%y=y
    ch(n)%xle=xle
    ch(n)%c=c
  END DO

  CALL MACofWing(ch(1:n), area,mac,ymac,xlemac)
  WRITE(*,*) " #        y        xle        c"
  WRITE(*,'(I3,3F12.6)')  (k,ch(k),k=1,n)
  WRITE(*,*) 'area=', area
  WRITE(*,*) 'length of MAC=', mac
  WRITE(*,*) 'y of MAC=',ymac
  WRITE(*,*) 'xle of MAC=', xlemac
  WRITE(*,*) 'xte of MAC=', xlemac+mac
  WRITE(*,*) 'x of c/4 of MAC=', xlemac+0.25*mac

  STOP


CONTAINS

!+
SUBROUTINE MACofWing(c,area,mac,ymac,xlemac)
! ---------------------------------------------------------------------------
! PURPOSE -
  TYPE(WINGCHORD),INTENT(IN),DIMENSION(:):: c
  REAL,INTENT(OUT):: area,mac,ymac,xlemac

  REAL:: aseg,cseg,xseg,yseg
  INTEGER:: k,n
!----------------------------------------------------------------------------
  area=0.0
  mac=0.0
  ymac=0.0
  xlemac=0.0

  n=SIZE(c)
  IF (n < 2) RETURN
!  write(*,*) n

  DO k=2,n
    CALL MACofOneSegment(c(k-1),c(k),aseg,cseg,yseg,xseg)
    area=area+aseg
    mac=mac+cseg*aseg
    ymac=ymac+yseg*aseg
    xlemac=xlemac+xseg*aseg
  END DO

  mac=mac/area
  ymac=ymac/area
  xlemac=xlemac/area

  RETURN
END Subroutine MACofWing   ! ------------------------------------------------

!+
SUBROUTINE MACofOneSegment(c1,c2,area,mac,ymac,xlemac)
! ---------------------------------------------------------------------------
! PURPOSE -
  TYPE(WINGCHORD),INTENT(IN):: c1,c2
  REAL,INTENT(OUT):: area,mac,ymac,xlemac

!  REAL:: croot,ctip
  REAL:: fraction
  REAL:: span
  REAL:: taper,taper1
!----------------------------------------------------------------------------
!  croot=c1%xte - c1%xle
!  ctip=c2%xte - c2%xle
  span=c2%y - c1%y

  area=0.5*(c1%c+c2%c)*span
  taper=c2%c/c1%c
  taper1=taper+1.0
  fraction=(taper+taper1)/(3.0*taper1)
  mac=c1%c*(taper*taper+taper+1.0)/(1.5*taper1)
  ymac=c1%y + fraction*span
  xlemac=c1%xle + fraction*(c2%xle-c1%xle)
!  write(*,*) taper,fraction
  RETURN
END Subroutine MACofOneSegment   ! ------------------------------------------

END Program GetMeanAerodynamicChord   ! =====================================