!+ Subroutine to control fitting of cloud properties in one band.
!
! Method:
!	Each of the quantites to be fitted is considered in turn.
!	A suitable scling of the data is applied for numerical
!	purposes and a conjugate gradient algorithm is called.
!
! Current Owner of Code: J. M. Edwards
!
! History:
!	Version		Date			Comment
!	1.0		12-04-95		First Version under RCS
!						(J. M. Edwards)
!
! Description of Code:
!   FORTRAN 77  with extensions listed in documentation.
!
!- ---------------------------------------------------------------------
      SUBROUTINE CLOUD_FIT_BAND(IERR
     &   , N_DATA, I_FIT_TYPE
     &   , DENSITY_PARTICLE
     &   , VOL_FRACTION, RADIUS_EFFECTIVE
     &   , SCATTERING, ABSORPTION, ASYMMETRY
     &   , CLOUD_PARAMETER_LIST)
!
!
!
      IMPLICIT NONE
!
!
!     INCLUDE HEADERS.
      INCLUDE 'incl/dimensions_pp.ucf'
      INCLUDE 'incl/dimensions_spec.ucf'
      INCLUDE 'incl/def_std_io.icf'
      INCLUDE 'incl/opt_variable.pcf'
      INCLUDE 'incl/cloud_parametrization.pcf'
      INCLUDE 'incl/cloud_parametrization.ncf'
      INCLUDE 'incl/error.pcf'
!
!     DUMMY ARGUMENTS.
      INTEGER 	!, INTENT(OUT)
     &     IERR
!             SUBROUTINE ERROR FLAG
      INTEGER 	!, INTENT(IN)
     &     N_DATA
!             NUMBER OF DATA POINTS
     &   , I_FIT_TYPE
!             TYPE OF FITTING REQUIRED
      REAL*8	!, INTENT(IN)
     &     DENSITY_PARTICLE
!             DENSITY OF PARTICLES
     &   , VOL_FRACTION(N_DATA)
!             VOLUME FRACTION OF ABSORBER
     &   , RADIUS_EFFECTIVE(N_DATA)
!             EFFECTIVE RADIUS OF DROPS
     &   , SCATTERING(N_DATA)
!             SCATTERING IN BAND
     &   , ABSORPTION(N_DATA)
!             ABSORPTION IN BAND
     &   , ASYMMETRY(N_DATA)
!             ASYMMETRY AVERAGED ACROSS BAND
      REAL*8	!, INTENT(OUT)
     &     CLOUD_PARAMETER_LIST(NPD_CLOUD_PARAMETER)
!             PARAMETERS OF CLOUD SCHEME
!
!     LOCAL VARIABLES.
      INTEGER
     &     I
!             LOOP VARIABLE
     &   , N_PARAMETER
!             NUMBER OF FITTING PARAMETERS
      REAL*8
     &     EXTINCTION(NPD_MIE_BLOCK)
!             TOTAL EXTINCTION
     &   , SS_ALBEDO(NPD_MIE_BLOCK)
!             SINGLE SCATTERING ALBEDO
     &   , RADIUS_SCALED(NPD_MIE_BLOCK)
!             DIMENSIONLESS EFFECTIVE RADIUS
     &   , EXTINCTION_SCALED(NPD_MIE_BLOCK)
!             DIMENSIONLESS EXTINCTION
     &   , PARAMETER_LOCAL(NPD_CLOUD_PARAMETER
     &        , NPD_CLOUD_OPTICAL_VARIABLE)
!             PARAMETERS EACH OPTICAL PROPERTY
     &   , RADIUS_MEAN
!             MEAN EFFECTIVE RADIUS
     &   , MEAN_EXTINCTION
!             MEAN EXTINCTION/UNIT VOLUME
     &   , SS_COALBEDO
!             COALBEDO
!
!
!     SUBROUTINES CALLED:
      EXTERNAL
     &     CONJUGATE_GRADIENT_CLOUD, BUILD_PARAMETER_LIST
     &   , EXTINCTION_GRADIENT, SS_ALBEDO_GRADIENT
     &   , ASYMMETRY_GRADIENT, SCATTERING_GRADIENT
     &   , ABSORPTION_GRADIENT
!
!     FUNCTIONS CALLED:
      REAL*8
     &     EXTINCTION_RESIDUAL
!             RESIDUAL FUNCTION FOR EXTINCTION
     &   , SS_ALBEDO_RESIDUAL
!             RESIDUAL FUNCTION FOR SCATTERING
     &   , ASYMMETRY_RESIDUAL
!             RESIDUAL FUNCTION FOR ASYMMETRY
     &   , SCATTERING_RESIDUAL
!             RESIDUAL FUNCTION FOR SCATTERING
     &   , ABSORPTION_RESIDUAL
!             RESIDUAL FUNCTION FOR ABSORPTION
      EXTERNAL
     &     EXTINCTION_RESIDUAL, SS_ALBEDO_RESIDUAL
     &   , ASYMMETRY_RESIDUAL, SCATTERING_RESIDUAL
     &   , ABSORPTION_RESIDUAL
!
!
!     INCLUDE HEADER FILES.
      INCLUDE 'incl/initial.pcf'
      INCLUDE 'incl/cloud_parm_data.dcf'
!
!
!
!     SEPARATE THE PARAMETRIZATION SCHEMES ACCORDING TO THE VARIABLES TO
!     BE FITTED.
      IF ( (I_FIT_TYPE.EQ.IP_SLINGO_SCHRECKER)
     &     .OR.(I_FIT_TYPE.EQ.IP_ACKERMAN_STEPHENS) ) THEN
!        FIT THE EXTINCTION, SINGLE SCATTERING ALBEDO AND ASYMMETRY.
         DO I=1, N_DATA
            EXTINCTION(I)=ABSORPTION(I)+SCATTERING(I)
            SS_ALBEDO(I)=SCATTERING(I)/EXTINCTION(I)
         ENDDO
         DO I=1, N_CLOUD_FIT_PARAMETER(I_FIT_TYPE, IP_EXTINCTION)
            PARAMETER_LOCAL(I, IP_EXTINCTION)=0.0D+00
         ENDDO
         DO I=1, N_CLOUD_FIT_PARAMETER(I_FIT_TYPE, IP_SS_ALBEDO)
            PARAMETER_LOCAL(I, IP_SS_ALBEDO)=0.0D+00
         ENDDO
         DO I=1, N_CLOUD_FIT_PARAMETER(I_FIT_TYPE, IP_ASYMMETRY)
            PARAMETER_LOCAL(I, IP_ASYMMETRY)=0.0D+00
         ENDDO
!
!        SCALE THE QUANTITIES FOR FITTING.
         RADIUS_MEAN=0.0D+00
         SS_COALBEDO=0.0D+00
         MEAN_EXTINCTION=0.0D+00
         DO I=1, N_DATA
            RADIUS_MEAN=RADIUS_MEAN+RADIUS_EFFECTIVE(I)
            SS_COALBEDO=MAX(SS_COALBEDO, 1.0D+00-SS_ALBEDO(I))
            MEAN_EXTINCTION=MEAN_EXTINCTION
     &         +EXTINCTION(I)/VOL_FRACTION(I)
         ENDDO
         RADIUS_MEAN=RADIUS_MEAN/DBLE(N_DATA)
         MEAN_EXTINCTION=MEAN_EXTINCTION/DBLE(N_DATA)
         DO I=1, N_DATA
            EXTINCTION_SCALED(I)=EXTINCTION(I)/(MEAN_EXTINCTION
     &         *VOL_FRACTION(I))
            RADIUS_SCALED(I)=RADIUS_EFFECTIVE(I)/RADIUS_MEAN
            SS_ALBEDO(I)=(1.0D+00-SS_ALBEDO(I))/SS_COALBEDO
         ENDDO
!
!
         N_PARAMETER=N_CLOUD_FIT_PARAMETER(I_FIT_TYPE, IP_EXTINCTION)
!
         WRITE(IU_MONITOR, '(/A ,/A)')
     &      '------------------------------------'
     &      , 'FITTING THE EXTINCTION.'
!
         CALL CONJUGATE_GRADIENT_CLOUD(IERR
     &      , N_DATA, I_FIT_TYPE, N_PARAMETER
     &      , VOL_FRACTION, RADIUS_SCALED
     &      , EXTINCTION_SCALED, PARAMETER_LOCAL(1, IP_EXTINCTION)
     &      , EXTINCTION_RESIDUAL, EXTINCTION_GRADIENT
     &      )
         IF (IERR.NE.I_NORMAL) THEN
            IF (IERR.EQ.I_ABORT_CALCULATION) THEN
!              RECOVER FORM MILD ERRORS OR RETURN FROM FATAL ONES.
               IERR=I_NORMAL
            ELSE
               RETURN
            ENDIF
         ENDIF
!
!        RESCALE PARAMETERS AND CONVERT DIVIDE BY THE DENSITY OF THE
!        PARTICLES TO ALLOW FOR THE SUBSEQUENT USE OF 
!        MASS MIXING RATIOS.
         IF (I_FIT_TYPE.EQ.IP_SLINGO_SCHRECKER) THEN
            PARAMETER_LOCAL(1, IP_EXTINCTION)
     &         =MEAN_EXTINCTION*PARAMETER_LOCAL(1, IP_EXTINCTION)
     &         /DENSITY_PARTICLE
            PARAMETER_LOCAL(2, IP_EXTINCTION)=RADIUS_MEAN
     &         *MEAN_EXTINCTION*PARAMETER_LOCAL(2, IP_EXTINCTION)
     &         /DENSITY_PARTICLE
         ELSE IF (I_FIT_TYPE.EQ.IP_ACKERMAN_STEPHENS) THEN
            PARAMETER_LOCAL(1, IP_EXTINCTION)
     &         =MEAN_EXTINCTION*PARAMETER_LOCAL(1, IP_EXTINCTION)
     &         /DENSITY_PARTICLE
            PARAMETER_LOCAL(2, IP_EXTINCTION)
     &         =MEAN_EXTINCTION*PARAMETER_LOCAL(2, IP_EXTINCTION)
     &         *EXP(-PARAMETER_LOCAL(3, IP_EXTINCTION)*LOG(RADIUS_MEAN))
     &         /DENSITY_PARTICLE
         ENDIF
!
!
         N_PARAMETER=N_CLOUD_FIT_PARAMETER(I_FIT_TYPE, IP_SS_ALBEDO)
!
         WRITE(IU_MONITOR, '(/A ,/A)')
     &      '------------------------------------'
     &      , 'FITTING THE COALBEDO.'
!
         CALL CONJUGATE_GRADIENT_CLOUD(IERR
     &      , N_DATA, I_FIT_TYPE, N_PARAMETER
     &      , VOL_FRACTION, RADIUS_SCALED
     &      , SS_ALBEDO, PARAMETER_LOCAL(1, IP_SS_ALBEDO)
     &      , SS_ALBEDO_RESIDUAL, SS_ALBEDO_GRADIENT
     &      )
         IF (IERR.NE.I_NORMAL) THEN
            IF (IERR.EQ.I_ABORT_CALCULATION) THEN
!              RECOVER FORM MILD ERRORS OR RETURN FROM FATAL ONES.
               IERR=I_NORMAL
            ELSE
               RETURN
            ENDIF
         ENDIF
!
!        RESCALE PARAMETERS.
         IF (I_FIT_TYPE.EQ.IP_SLINGO_SCHRECKER) THEN
            PARAMETER_LOCAL(1, IP_SS_ALBEDO)=SS_COALBEDO
     &         *PARAMETER_LOCAL(1, IP_SS_ALBEDO)
            PARAMETER_LOCAL(2, IP_SS_ALBEDO)=SS_COALBEDO
     &         *PARAMETER_LOCAL(2, IP_SS_ALBEDO)/RADIUS_MEAN
         ELSE IF (I_FIT_TYPE.EQ.IP_ACKERMAN_STEPHENS) THEN
            PARAMETER_LOCAL(1, IP_SS_ALBEDO)=SS_COALBEDO
     &         *PARAMETER_LOCAL(1, IP_SS_ALBEDO)
            PARAMETER_LOCAL(2, IP_SS_ALBEDO)=SS_COALBEDO
     &         *PARAMETER_LOCAL(2, IP_SS_ALBEDO)
     &         *EXP(-PARAMETER_LOCAL(3, IP_SS_ALBEDO)*LOG(RADIUS_MEAN))
         ENDIF
!
!
         N_PARAMETER=N_CLOUD_FIT_PARAMETER(I_FIT_TYPE, IP_ASYMMETRY)
!
         WRITE(IU_MONITOR, '(/A ,/A)')
     &      '------------------------------------'
     &      , 'FITTING THE ASYMMETRY.'
!
         CALL CONJUGATE_GRADIENT_CLOUD(IERR
     &      , N_DATA, I_FIT_TYPE, N_PARAMETER
     &      , VOL_FRACTION, RADIUS_SCALED
     &      , ASYMMETRY, PARAMETER_LOCAL(1, IP_ASYMMETRY)
     &      , ASYMMETRY_RESIDUAL, ASYMMETRY_GRADIENT
     &      )
         IF (IERR.NE.I_NORMAL) THEN
!           RECOVER FORM MILD ERRORS OR RETURN FROM FATAL ONES.
            IF (IERR.EQ.I_ABORT_CALCULATION) THEN
               IERR=I_NORMAL
            ELSE
               RETURN
            ENDIF
         ENDIF
!        RESCALE PARAMETERS
         IF (I_FIT_TYPE.EQ.IP_SLINGO_SCHRECKER) THEN
            PARAMETER_LOCAL(2, IP_ASYMMETRY)
     &         =PARAMETER_LOCAL(2, IP_ASYMMETRY)/RADIUS_MEAN
         ELSE IF (I_FIT_TYPE.EQ.IP_ACKERMAN_STEPHENS) THEN
            PARAMETER_LOCAL(2, IP_ASYMMETRY)
     &         =PARAMETER_LOCAL(2, IP_ASYMMETRY)
     &         *EXP(-PARAMETER_LOCAL(3, IP_ASYMMETRY)*LOG(RADIUS_MEAN))
         ENDIF
!
!        ASSIGN THE LOCAL ARRAYS OF PARAMETERS TO THE MAIN ARRAY.
         CALL BUILD_PARAMETER_LIST(I_FIT_TYPE
     &      , PARAMETER_LOCAL, CLOUD_PARAMETER_LIST)
      ELSE
         WRITE(IU_ERR, '(/A)')
     &      '***ERROR: ILLEGAL CLOUD PARAMETRIZATION SELECTED.'
         IERR=I_ERR_FATAL
         RETURN
      ENDIF
!
!
      RETURN
      END
