C **********************************************************************
C *                                                                    *
C *                                SOFTWARE LICENSING                  *
C *                                                                    *
C *                This program is free software; you can redistribute *
C *                it and/or modify it under the terms of the GNU      *
C *                General Public License as published by the Free     *
C *                Software Foundation, either Version 2 of the        *
C *                license, or (at your option) any later version.     *
C *                                                                    *
C *                This program is distributed in the hope that it     *
C *                will be useful, but without any warranty; without   *
C *                even the implied warranty of merchantability or     *
C *                fitness for a particular purpose. See the GNU       *
C *                General Public License for more details.            *
C *                                                                    *
C *                A copy of the GNU General Public License is         *
C *                available at http://www.gnu.org/copyleft/gpl.html   *
C *                or by writing to the Free Software Foundation, Inc.,*
C *                59 Temple Place - Suite 330, Boston, MA 02111, USA. *
C *                                                                    *
C **********************************************************************

      SUBROUTINE WETBULB(RH,TDRY,BP,DEWP,TWET)
C
CNG   RH IN PERCENT, TDRY IN DEG C, BP IN MB
CNG   COMPUTE DEW POINT AND WET BULB TEMPERATURE
CNG   USING THE BISECTION METHOD
CNG 
CNG   N.GEORGAS, 10/26/2011, AFTER TVA APPENDIX D
CNG 
      NITER=50

      A1=7.5 ! 9.5 for ice
      B1=237.3 ! 265.5 for ice
      C1=0.7858 ! same for ice

      TDRY1=AMAX1(TDRY,0.01) ! DO NOT ALLOW ICE/SNOW FOR NOW
      TDRY1=AMIN1(TDRY,40.0) ! NOR REALLY HIGH TEMPERATURES
      RH=AMAX1(RH,0.1)
      RH=AMIN1(RH,100.)

CNG   COMPUTE PARTIAL WATER VAPOR PRESSURE OF DRY AIR
      EA = RH/100.*EXP(2.303*(A1*TDRY1/(TDRY1+B1)+C1))
CNG   COMPUTE DEW POINT TEMPERATURE (ANCILLARY)
      DEWP=237.3*(ALOG10(EA)-0.7858)/(7.5+0.7858-ALOG10(EA))
CNG   COMPUTE LATENT HEAT OF VAPORIZATION AT DRY BULB TEMP
      HLTDRY=2500.-2.39*TDRY1
CNG   COMPUTE MIXING RATIO OF AIR
      RDRY=0.622*EA/(BP-EA)

CNG   ASSUME INITIAL WET BULB DEPRESSION BOUNDS
      WETDEP1=25.
      WETDEP2=0.

CNG   WET BULB TEMPERATURE BOUNDS
      TWET0=TDRY1-WETDEP1
      TWET2=TDRY1-WETDEP2

      DO ITER=1,NITER
       TWET1=0.5*(TWET0+TWET2)
       IF(ABS(WETBCHK(TWET1,BP,HLTDRY,RDRY)).LT.1.E-7) GOTO 10
       IF(WETBCHK(TWET1,BP,HLTDRY,RDRY).GT.0) THEN
        TWET0=TWET1
       ELSE
        TWET2=TWET1
       ENDIF
      ENDDO
 10   TWET=TWET1
      TWET=AMIN1(50.,AMAX1(0.,TWET)) ! Make sure it is reasonable and not ice
CNG      WRITE (*,*)'FOUND',TWET0,TWET,TWET2,WETBCHK(TWET,BP,HLTDRY,RDRY)
CNG      WRITE (*,*)'IN',ITER,'ITERATIONS'
      RETURN
      END

      REAL FUNCTION WETBCHK(TWET,BP,HLTDRY,RDRY)
      A1=7.5 ! 9.5 for ice
      B1=237.3 ! 265.5 for ice
      C1=0.7858 ! same for ice
      CVAPOR1=0.45 ! kcal/(kg*C)
      CWATER1=0.999 ! kcal/(kg*C)
      CDRYAR1=0.24 ! kcal/(kg*C)
CNG   COMPUTE PARTIAL WATER VAPOR PRESSURE OF AIR AT WET BULB
      EAW=EXP(2.3026*(A1*TWET/(TWET+B1)+C1))
CNG   COMPUTE MIXING RATIO OF AIR AT WET BULB
      RWET=0.622*EAW/(BP-EAW)
CNG   COMPUTE LATENT HEAT OF VAPORIZATION AT WET BULB TEMP
      HLTWET=2500.-2.39*TWET
CNG   COMPUTE LEFT AND RIGHT HAND SIDES OF THE THERMODYNAMIC WET BULB EQUATION
      ELEFT=ALOG10(HLTWET/HLTDRY)
      ERIGHT=(CVAPOR1-CWATER1)/CVAPOR1*
     .  ALOG10( (CDRYAR1+CVAPOR1*RDRY)/(CDRYAR1+CVAPOR1*RWET) )
      WETBCHK=ELEFT-ERIGHT ! RESIDUAL
      RETURN
      END