* Last processed by NICE on * Customized for : EEEI HPUX UNIX SUBROUTINE CHOPPE ( MEANCH, MEANCO, MEANAT, $FREQSI, FREQIM, ELEV, TCALR, TSYSR, TRECR, TAMB, PAMB, ALTI) C---------------------------------------------------------------------- C compute tcal, tsys, tref, from the chop and atm values observed C according to the usual chopper method (Ulich and Haas) C C Double side band operation C C If CAL_MODE is 'AU' and a cold load is available, or if it C is 'TR', derive the sky emission TEMI from the receiver C temperature TREC, and then CAL_WH2O from the best model fitting. C Else, if CAL_MODE is 'AU', the sky emission is computed from C the given CAL_WH2O and if CAL_MODE is 'MA', the sky emission is C computed from the given opacities and atmospheric temperatures. C C Uses a curved atmosphere, equivalent height 5.5 KM C C Input: C MEANCH R counts on chopper C MEANCO R counts on cold load C MEANAT R counts on sky C FREQSI R frequency in signal side band C FREQIM R frequency in image side band C ELEV R elevation (degrees) C TRECR R receiver dsb temperature (K) C TAMB C PAMB C Output: C tcalr R calibration temperature (K) C tsysr R system temperature (K) C trecr R receiver dsb temperature (K) C * Global --------------------------------------------------------------------- * INCLUDE 'calibr.inc' * Dummy --------------------------------------------------------------------- REAL*4 MEANCH, MEANCO, MEANAT, $ELEV, TCALR, TSYSR, TRECR, FREQSI, FREQIM, TAMB, PAMB, ALTI * Local ---------------------------------------------------------------------- REAL EM, DEM, EPS, GAMMA, HZ, AIRMAS, DELEV, $TATM, TAUT, DWATER, ERR, DELTAW, $ATTSIG, TEMISI, TEMIIM, TSKY, TEMI INTEGER K, IER REAL H0, R0 PARAMETER (H0= 5.5, R0= 6370.) CHARACTER*25 MESS(3) DATA MESS $/'Zero atmospheric opacity', $'No oxygen in atmosphere', $'No water in atmosphere'/ *------------------------------------------------------------------------ * * Check for bad atmosphere IF (MEANAT.GE.MEANCH) THEN CALL MESSAGE (6,3,'CHOPPER','Bad atmosphere') RETURN ENDIF IF (MEANCO.GE.MEANCH) THEN CALL MESSAGE (6,3,'CHOPPER', $ 'Signal stronger on COLD than on CHOPPER') RETURN ENDIF * DELEV = ELEV ! in radian EPS= ASIN(R0/(R0+H0)*COS(DELEV)) GAMMA= DELEV + EPS HZ = R0*R0 + (R0+H0)**2 - 2. * R0 * (R0+H0) * SIN(GAMMA) HZ = SQRT (MAX(H0**2,HZ)) AIRMAS = HZ/H0 AIRMAS = MIN (AIRMAS, 20.) * IF ((CAL_MODE.EQ.'AU'.AND. MEANCO.GT.0) $.OR. CAL_MODE.EQ.'TR') THEN * * AUTO mode with cold load or fixed TREC : minimize for water vapor * content CALL ATM_ATMOSP(TAMB, PAMB, ALTI) * * AUTO mode with cold load : compute TREC IF (MEANCO.NE.0) THEN TRECR = ( CAL_TCOL * MEANCH - $ CAL_TCHO * MEANCO ) / (MEANCO - MEANCH) ENDIF * * Compute TEMI *! TEMI = (CAL_TCHO + TRECR) * MEANAT/MEANCH - TRECR * Introducing the soft calibration RL 14-oct-1989 TEMI = ((CAL_CEFF*CAL_TCHO+TRECR)*MEANAT/MEANCH-TRECR) $ / (1.-MEANAT/MEANCH*(1.-CAL_CEFF) ) * * Strictly equivalent to : * TEMI = CAL_TCHO - (CAL_TCHO - CAL_TCOL) * * + (MEANCH-MEANAT) / (MEANCH-MEANCO) * * Should use a weighted average between OBsATM and OBsCHO since some * losses are within the cabin TSKY = (TEMI - (1.-CAL_FEFF) * TAMB) / CAL_FEFF ERR = 1E10 DWATER = .01 K = 0 10 CONTINUE IF (K.GE.10) THEN CALL MESSAGE (6,2,'CHOPPER', $ 'Did not converge after 10 iterations') RETURN ENDIF K = K + 1 CALL EMISS (DWATER, EM, DEM, FREQSI, $ FREQIM, AIRMAS) ERR = EM - TSKY IF (DEM.EQ.0.) THEN CALL MESSAGE (6,2,'CHOPPER', $ 'Blocked in a local minimum') RETURN ENDIF DELTAW = DWATER * ERR / DEM IF (ABS(DWATER).GE.ABS(DELTAW/10.)) DWATER = DELTAW/10. CAL_WH2O = CAL_WH2O - DELTAW ERR = ABS(ERR) IF (ERR.GE.0.1) GOTO 10 CALL ATM_TRANSM (CAL_WH2O,AIRMAS,FREQSI, $ TEMISI,CAL_ATMS,TAUOXS(IFE),TAUWS(IFE), $ TAUT,IER) CAL_TAUS = TAUOXS(IFE) + TAUWS(IFE) CALL ATM_TRANSM (CAL_WH2O,AIRMAS,FREQIM, $ TEMIIM,CAL_ATMI,TAUOXI(IFE),TAUWI(IFE), $ TAUT,IER) CAL_TAUI = TAUOXI(IFE) + TAUWI(IFE) IF (IER.NE.0) THEN CALL MESSAGE (6,2,'CHOPPER','Stupib calibration ' $ //MESS(IER)) ENDIF ELSEIF (CAL_MODE.EQ.'AU') THEN * * Auto Mode without cold load, or when the water vapor content is known CALL ATM_ATMOSP (TAMB, PAMB, 2.552) ! press at 2550m CALL ATM_TRANSM (CAL_WH2O,AIRMAS,FREQSI, $ TEMISI,CAL_ATMS,TAUOXS(IFE),TAUWS(IFE), $ TAUT,IER) CAL_TAUS = TAUOXS(IFE) + TAUWS(IFE) CALL ATM_TRANSM (CAL_WH2O,AIRMAS,FREQIM, $ TEMIIM,CAL_ATMI,TAUOXI(IFE),TAUWI(IFE), $ TAUT,IER) CAL_TAUI = TAUOXI(IFE) + TAUWI(IFE) IF (IER.NE.0) THEN CALL MESSAGE (6,2,'CHOPPER','Stupib calibration ' $ //MESS(IER)) ENDIF EM = (TEMISI + TEMIIM * CAL_GAIN) / (1.+CAL_GAIN) * * Should use a weighted average between TAMB and CAL_TCHO since some * losses are within the cabin TEMI = (1.-CAL_FEFF)*TAMB + CAL_FEFF * EM * * Manual Mode ELSE TEMISI = CAL_ATMS * (1.-EXP(-CAL_TAUS*AIRMAS)) TEMIIM = CAL_ATMI * (1.-EXP(-CAL_TAUI*AIRMAS)) EM = (TEMISI + TEMIIM * CAL_GAIN) / (1.+CAL_GAIN) * * Should use a weighted average between OBsATM and OBsCHO since some * losses are within the cabin TEMI = (1.-CAL_FEFF)*TAMB + CAL_FEFF * EM ENDIF * * Compute TCAL at zenith TCALR = (CAL_TCHO-TEMI) * (1.+CAL_GAIN) / CAL_FEFF * * Compute TSYS ATTSIG = EXP (-CAL_TAUS * AIRMAS) TSYSR = TCALR * MEANAT / (MEANCH - MEANAT) / ATTSIG * * Compute TREC TRECR = (TEMI * (MEANCH-MEANAT*(1.0-CAL_CEFF)) $- CAL_CEFF*CAL_TCHO * MEANAT) $/ (MEANAT - MEANCH) END * SUBROUTINE EMISS (DWATER,EM,DEM,FREQSI,FREQIM,AIRMAS) C---------------------------------------------------------------------- C computes the combined atmospheric emission of both sidebands C and its variation C C Input: C water1 R water vapor content (mm) C DWATER R increment of above C FREQSI R signal frequency C FREQIM R image frequency C AIRMAS R number of air masses C Output: C em R atmospheric radiation temperature (K) C DEM R variation of above for increment of C water content = DWATER C---------------------------------------------------------------------- INCLUDE 'inc:parameter.inc' INCLUDE 'calibr.inc' * INTEGER IER REAL FREQSI, FREQIM, AIRMAS, DWATER, EM, DEM, $TATM, TAUT, TEMISI, TEMIIM, T1, T2 * CALL ATM_TRANSM (CAL_WH2O, AIRMAS, FREQSI, TEMISI, $CAL_ATMS, T1, T2, TAUT, IER) CAL_TAUS = T1 + T2 CALL ATM_TRANSM (CAL_WH2O, AIRMAS, FREQIM, TEMIIM, $CAL_ATMI, T1, T2, TAUT, IER) CAL_TAUI = T1 + T2 EM = (TEMISI + TEMIIM*CAL_GAIN) / (1. + CAL_GAIN) CALL ATM_TRANSM (CAL_WH2O+DWATER, AIRMAS, FREQSI, $TEMISI, TATM, T1, T2, TAUT, IER) CALL ATM_TRANSM (CAL_WH2O+DWATER, AIRMAS, FREQIM, $TEMIIM, TATM, T1, T2, TAUT, IER) DEM = (TEMISI + TEMIIM*CAL_GAIN) / (1.+CAL_GAIN) - EM END