SUBROUTINE NIMS_CDS2OHMS(ICHNL,NDATA,OHMS,IR) C_TITLE NIMS_CDS2OHMS Compute resistance for given CDS channel and data number C_ARGS INTEGER ICHNL !I Channel number (see DESC below) INTEGER NDATA !I Data number for channel ICHNL (range 0-255) REAL OHMS !O Resistance corresponding to ICHNL and NDATA INTEGER IR !O Return code: 0 OK, -1 ICHNL invalid C_DESC This subroutine applies to the following CDS channels: C ICHNL CDS CHANNEL LOCATION SENSOR SE TERMINALS C ===== =========== ======== ====== ============ C 1 E-1911 Shield 15319 L - M C 2 E-1912 Telescope 15322 F - G C 3 E-1913 Gr. Mech. 15323 D - E C 4 E-1914 Chopper ? B - C C 5 E-1915 Electronics UA-75 N - d C Values for E-1910 and E-1916, both of which relate C to the FPA temperature, are handled in subroutine NIMS_FPATMP. C In this routine, we use the given CDS coefficients to C determine a resistance. In the first four cases, the C answer is found by solving a cubic equation. The C electronics temperature is linearly related to the DN C value. C Temperature may be computed from OHMS by routine NIMS_OHMS2TMP. C NOTE: We are solving the equation C DN = C3*OHMS^3 + C2*OHMS^2 + C1*OHMS + C0 C The C's are from GLL 625-308 C REFERENCES: GLL Telemetry Conversion Handbook, 625-308 C Appendices A,B of NIMS calibration document. C_CALLS U_FIND_ROOT3 C_HIST 16-AUG-88 R Carlson JPL ORIGINAL VERSION C 20-FEB-93 R Mehlman UCLA/IGPP Renamed from CDS2OHMS, renamed U_ROOT3 C and changed argument order, error parameter added, C internal documentation put in ISIS format. C_END DIMENSION C0(5), C1(5), C2(5), C3(5) ! CDS coefficients DATA C0/ -53.2668, -52.8356, -51.8332, -51.1541, & -209.424/ DATA C1/ 0.745103, 0.742733, 0.737517, 0.746692, & 0.719420/ DATA C2/-0.537878E-03, -0.525888E-03, -0.512928E-03, & -0.520592E-03, -0.000000E-03/ DATA C3/ 0.251498E-06, 0.237182E-06, 0.226003E-06, & 0.204892E-06, 0.000000E-06/ EPS = 0.001 ! Error parameter IF(ICHNL.LT.1) GO TO 750 ! Out of range IF(ICHNL.GT.5) GO TO 750 ! Out of range DN = FLOAT(NDATA) RTRY = (DN - C0(ICHNL))/C1(ICHNL) ! Linear approximation IF(ICHNL.EQ.5) GO TO 800 ! Linear is adequate for chnl 5 CALL U_FIND_ROOT3(C3(ICHNL),C2(ICHNL),C1(ICHNL),C0(ICHNL)-DN,EPS, 1 RTRY,OHMS) GO TO 900 750 IR = -1 C 750 TYPE *,'ICHNL out of range' RETURN 800 OHMS = RTRY ! For channel 5 900 IR = 0 RETURN END SUBROUTINE NIMS_DET_MODE(KHRS,LPTAB,KPTAB,KMF, 1 MRPT,MIROP,MGSTART,MGDELTA,MGSTEPS,KMODE) C_TITLE NIMS_DET_MODE Determine NIMS instrument mode, etc. from HRS hskpng, etc C_ARGS INPUT: BYTE KHRS(96,10) !I 10 packets of NIMS HRS data: PTABs C are in (3-6,p), p=2,7 BYTE LPTAB(4,2) !I LRS PTABs (acquired in mf 0-2) BYTE KPTAB !I Current PTAB: '00'X=1, '80'X=2 C from LRS. Good in mf 11-65 & 68-90 INTEGER KMF !I Current minor frame (0-90) C OUTPUT FROM CURRENT PTAB: INTEGER MRPT !O Mode repeat count LOGICAL MIROP !O Mirror scan (True or False) INTEGER MGSTART !O (Logical) grating start position INTEGER MGDELTA !O Grating delta INTEGER MGSTEPS !O Steps in grating cycle C OUTPUT MODE INTERPRETATION: INTEGER KMODE !O Mode determined by routine C 1,2: full map, spectrometer C 3,4: long map, spectrometer C 5,6: short map, spectrometer C 7,0: fixed map, spectrometer (safe) C 8,9: bandedge map, spectrometer C -1: undetermined C_DESC Extract current PTAB, determine NIMS mode and grating cycle parameters, C using HRS housekeeping and possibly LRS housekeeping values. C References: NIMS instrument paper, and GLL 3-280 C_LIMS LRS input not yet used in analysis. C_CALLS U_MOVE1 (ISISLIB) C_HIST 03nov92 R Mehlman, UCLA/IGPP Entirely rewritten from 31aug84 MODET. C 06dec92 RMX Added more output arguments: MRPT,MIROP,MGSTART,MGDELTA,MGSTEPS C 08feb92 RMX Corrected test for bandedge modes (PTAB(4)=1) C 09mar93 RMX Changed from function to subroutine, input override C mode removed, parameters renamed, LRS input added. C_END C Instrument mode tables: position n corresponds to mode n-1 C Actual # grating positions (to be used in cube label) C Nits per mirror cycle (fixed grating cycles use full mirror cycle) C Entries for modes 10-15 . C For modes 0-9, mirror moving iff mode is odd. BYTE IBYTE2 C Determine mode from HRS hskpng MIROP=.FALSE. MGSTEPS=KHRS(6,2) IF (MGSTEPS.EQ.0) THEN MGSTEPS=KHRS(6,7) IF (MGSTEPS.EQ.0) GO TO 50 MRPT=KHRS(3,7) IF (KHRS(4,7).LT.0) MIROP=.TRUE. CALL U_MOVE1(1,KHRS(4,7),IBYTE2) MGDELTA=KHRS(5,7) ELSE MRPT=KHRS(3,2) IF (KHRS(4,2).LT.0) MIROP=.TRUE. CALL U_MOVE1(1,KHRS(4,2),IBYTE2) MGDELTA=KHRS(5,2) END IF MGSTART=IBYTE2.AND.'3F'X IF (MGSTEPS.EQ.24) KMODE=3 IF (MGSTEPS.EQ.12) THEN KMODE=1 IF (MGDELTA.EQ.0) KMODE=7 END IF IF (MGSTEPS.EQ.6) KMODE=5 IF (MGSTEPS.EQ.1) KMODE=8 IF (.NOT.MIROP) THEN KMODE=KMODE+1 IF (MGDELTA.EQ.0) KMODE=0 END IF GO TO 80 C Undetermined 50 KMODE=-1 80 CONTINUE C WRITE (6,90) ((KHRS(I,J),I=3,6),J=2,7,5),LPTAB,KPTAB,KMF, C 1 MRPT,MIROP,MGSTART,MGDELTA,MGSTEPS,KMODE C 90 FORMAT (' * NIMS_MODE in: ',8Z3,2X,8Z3,Z5,I3/ C 1 ' out: ',I3,1X,L1,4I3) RETURN END SUBROUTINE NIMS_ENGR(ENGR,TDN,ODN,TEMP,RESIST) C_TITLE NIMS_ENGR Extract NIMS Engineering DNs, compute temperatures, etc. C_ARGS BYTE ENGR(2,91) !I Engineering items, 2/frame over RIM INTEGER TDN(7,6) !O Temperature DNs, 6 kinds, 7 per RIM C 1 Focal Plane Assembly C 2 Radiator Shield C 3 Telescope C 4 Grating C 5 Chopper C 6 Electronics INTEGER ODN(4) !O Other 4 items, DNs, 1 per RIM C 1 Flexprint temperature C 2 RCT-Pt temperature C 3 RCT-Ni temperature C 4 RCT-reference resistance REAL*4 TEMP(7,7) !O Temperatures (deg K) corresponding C to TDN, plus alternate TFPA (from C composite Pt + Flexprint (SS)), C 7 per RIM REAL*4 RESIST(7,4) !O Resistances (ohms), 7 per RIM C 1 RKA: R across terminals K,a C 2 RHPRK: RH + RK C 3 RHK: R across H,K (RH+RK+Rpt) C 4 RPT: R of Pt sensor C_DESC Extract NIMS Engineering DNs from RIM's worth of the 2 Engineering C bytes on the EDR; compute corresponding temperatures and resistances. C Item ID's, names and frequencies are: C E-1910 Focal Plane Temperature 7 per RIM C E-1911 Radiator Shield Temperature 7 per RIM C E-1912 Telescope Temperature 7 per RIM C E-1913 Grating Mechanism Temperature 7 per RIM C E-1914 Optical Chopper Temperature 7 per RIM C E-1915 Electronics Chassis Temperature 7 per RIM C E-1916 Flexprint Temperature 1 per RIM C E-1945 RCT-NIMS Temperature-Pt 1 per RIM C E-1946 RCT-NIMS Temperature-Ni 1 per RIM C E-1947 RCT-NIMS Reference R(esistance) 1 per RIM C Reference: GLL 3-280 Table A2.2.8 (Engineering Measurements) C See NIMS_* subroutines for details of temperature computations. C_CALLS U_FILL1, U_MOVE1 (isislib) C NIMS_FPATMP, NIMS_CDS2OHMS, NIMS_OHMS2TMP (nimslib) C_HIST 20feb93, R Mehlman, UCLA/IGPP ORIGINAL VERSION C_END INTEGER TFRAME(6) /0,0,1,1,2,2/ ! Temperature DN frames (mod 13, 7/RIM) INTEGER TBYTE(6) /1,2,1,2,1,2/ ! Temperature DN bytes INTEGER OFRAME(4) /56,29,58,85/ ! Other item DN frames (of 0-90) INTEGER OBYTE(4) / 1, 1, 1, 1/ ! Other item DN bytes REAL EPS /.001/ ! Epsilon for NIMS_FPATMP DO J=1,6 ! Temp. DNs IFRAME = TFRAME(J) + 1 DO I=1,7 CALL U_FILL1( 0, 4, TDN(I,J)) CALL U_MOVE1( 1, ENGR( TBYTE(J), IFRAME), TDN(I,J) ) IFRAME = IFRAME + 13 ENDDO ENDDO DO I=1,4 ! Other DNs CALL U_FILL1( 0, 4, ODN(I)) CALL U_MOVE1( 1, ENGR( OBYTE(I), OFRAME(I) + 1), ODN(I) ) ENDDO DO I=1,7 ! Compute: CALL NIMS_FPATMP( TDN(I,1), ODN(1), EPS, RESIST(I,1), ! FPA temp, R's 1 RESIST(I,2), RESIST(I,3), RESIST(I,4), TEMP(I,7), TEMP(I,1) ) DO J=2,6 ! Other temps CALL NIMS_CDS2OHMS( J-1, TDN(I,J), TOHM, IR) CALL NIMS_OHMS2TMP( J-1, TOHM, TCEN, TEMP(I,J), IR) ENDDO ENDDO RETURN END SUBROUTINE NIMS_FPATMP(N1910,N1916,EPS,RKA,RHPRK,RHK,RPT, 1 TPTSS,TPT) C_TITLE NIMS_FPATMP Compute FPA temperature from CDS engineering data C_ARGS INTEGER N1910 !I Data number for CDS channel E-1910 INTEGER N1916 !I Same for E-1915 (flexprint) REAL EPS !I Error parameter for cubic root routine REAL RKA !O Resistance across terminals K,a REAL RHPRK !O RH + RK REAL RHK !O Resistance across H,K (RH + RK _ RPt) REAL RPT !O Resistance of Pt sensor (#15604, FPA 002) REAL TPTSS !O FPA temperature from composite Pt + C flexprint (SS), in K REAL TPT !O FPA temperature from Pt, in K (USE THIS) C_DESC Computes FPA temperature from CDS engineering data for P/F unit. C The principal output FPA temperature is TPT, which uses just the C platinum sensor after subtracting out the simultaneously measured C flexprint resistance. The other output FPA temperature, TPTSS, C assumes the flexprint resistance is also temperature sensitive C (slightly). The combination of the two is a "resistance thermometer". C TPT and TPTSS should agree, providing the optics are cold. C References: Appendix A, B OF calibration document, C GLL Telemetry Conversion Handbook, 625-308 C_CALLS U_FIND_ROOT3 C_HIST 04-SEP-85 R Carlson JPL ORIGINAL VERSION C 20-FEB-93 R Mehlman UCLA/IGPP Renamed from FPATMP, renamed U_ROOT3 C and changed argument order, internal documentation put C in ISIS format. C_END C01910 = -0.100047619E+03 ! CDS coeff, E-1910 C11910 = 0.680503115E+00 C21910 = -0.342859039E-03 C31910 = 0.925943766E-07 C01916 = -0.246681101E+03 ! For E-1916 C11916 = 0.721500002E+00 RF = 413.91 ! Ohms, fixed R in FPA ETA = 0.98698 ! (RH + RK)/(RK + Ra) A0PT = -76.64005 ! For Pt R, R=A0 + A1*T + A2*T**2 A1PT = 2.265692 A2PT = -5.6604762E-04 A0SS = 70.83 ! For SS flexprint resistance A1SS = 0.41932 A2SS = -1.8061E-03 B0 = 34.20093 ! For PT, T = B0 + B1*R + B2*R**2 B1 = 0.4473015 B2 = 6.119418E-05 C A0 = A0PT + A0SS ! For composite Pt + SS A1 = A1PT + A1SS A2 = A2PT + A2SS DN1910 = FLOAT(N1910) DN1916 = FLOAT(N1916) RKA = (DN1916 - C01916)/C11916 ! E-1916 is linear RHPRK = ETA*(RKA - RF) ! from least squares linear fit RTRY = (DN1910 - C01910)/C11910 ! Initial guess at RHK CALL U_FIND_ROOT3(C31910,C21910,C11910,C01910-DN1910,EPS,RTRY,RHK) C = A0 - RHK ! For quadratic equation B = A1 ! For composite RSS + RPt A = A2 TPTSS = (-B + SQRT(B**2-4.*A*C))/(2.*A) RPT = RHK - RHPRK TPT = B0 + B1*RPT + B2*RPT**2 GO TO 800 C DATA IOLU /6/ C WRITE(IOLU,200) C WRITE(IOLU,201) RKA C WRITE(IOLU,202) RHPRK C WRITE(IOLU,203) RHK C WRITE(IOLU,204) RPT C WRITE(IOLU,205) TPTSS C WRITE(IOLU,206) TPT C WRITE(IOLU,207) N1910 C WRITE(IOLU,208) N1916 C WRITE(IOLU,209) C 200 FORMAT(1H1,5X,'NIMS P/F FPA Temperature From CDS Values') C 201 FORMAT(1H0,5X,'R(K-a)...............................',F7.2) C 202 FORMAT(1H ,5X,'R(H) + R(K)..........................',F7.2) C 203 FORMAT(1H ,5X,'R(H-K)...............................',F7.2) C 204 FORMAT(1H ,5X,'R(Pt)................................',F7.2) C 205 FORMAT(1H ,5X,'Temperature(Pt+SS)...................',F7.2) C 206 FORMAT(1H ,5X,'Temperature(Pt)......................',F7.2) C 207 FORMAT(1H ,5X,'DN (E-1910)..........................',I4) C 208 FORMAT(1H ,5X,'DN (E-1916)..........................',I4) C 209 FORMAT(1H) 800 CONTINUE RETURN END SUBROUTINE NIMS_GET_MODE_INFO(KMODE,KNGP,KNITS) C_TITLE NIMS_GET_MODE_INFO Get grating & mirror cycle info for NIMS mode C_ARGS INTEGER KMODE !I NIMS instrument mode C 1,2: full map, spectrometer C 3,4: long map, spectrometer C 5,6: short map, spectrometer C 7,0: fixed map, spectrometer (safe) C 8,9: bandedge map, spectrometer C 10-15: INTEGER KNGP !O Actual # grating positions in cycle INTEGER KNITS !O Number of nits in "mirror" cycle C_HIST 10mar93 R Mehlman UCLA/IGPP ORIGINAL VERSION C_END C Instrument mode tables: position n corresponds to mode n-1 C Actual # grating positions (to be used in cube label) C Nits per mirror cycle (fixed grating cycles use full mirror cycle) C Entries for modes 10-15 . C For modes 0-9, mirror moving iff mode is odd. INTEGER NGPOS(16) / 1,12,12,24,24,6,6, 1,2,2,0,0,0,0,0,0/ ! GP/CYCLE INTEGER NITSPC(16) /13,13,13,26,26,7,7,13,4,4,0,0,0,0,0,0/ ! NITS/CYCLE IF (KMODE.GE.0.AND.KMODE.LE.15) THEN KNGP=NGPOS(KMODE+1) KNITS=NITSPC(KMODE+1) ELSE KNGP=0 KNITS=0 END IF RETURN END SUBROUTINE NIMS_HSKPNG(LHSK,HHSK,CMD,HSKSTAT,DCHAN,IGS,ICHM) C_TITLE NIMS_HSKPNG Extract NIMS housekeeping items from LRS and HRS for RIM C_ARGS BYTE LHSK(3,91) !I LRS bytes, 3/frame over RIM BYTE HHSK(6,10,91) !I HRS hskpng items, 6/mod10 over RIM BYTE CMD(7,2) !O NIMS command (7 bytes: LRS, HRS) BYTE HSKSTAT(10,2) !O Housekeeping status bytes (LRS, HRS) C 1 Hardware status C 2 Software status C 3 Transaction bus parity errors C 4 All bus parity errors C 5 Power supply input C 6 Ave grating drive current C 7 Ave mirror drive current C 8 Reference voltage C 9 Optics cal source C 10 ROM checksum INTEGER*2 DCHAN(2) !O Si & InSb channels (10-bit LRS only) INTEGER IGS(2) !O Gain state (LRS, HRS) 1-4 INTEGER ICHM(2) !O Chopper mode (LRS, HRS) 1 reference C 2 63-hz 3 free-run 4 off C_DESC Extract NIMS housekeeping from both LRS and HRS hskpng streams. C Names and locations are (frames 0-90, packets 0-9, bytes 1-*): INTEGER LCFRAME(7) / 3, 3, 3, 4, 4, 4, 5/ ! LRS command buffer INTEGER LCBYTE(7) / 1, 2, 3, 1, 2, 3, 1/ INTEGER HCFRAME(7) / 0, 0, 0, 0, 0, 0, 0/ ! HRS command buffer INTEGER HCPACKET(7)/ 8, 8, 8, 8, 8, 8, 9/ INTEGER HCBYTE(7) / 1, 2, 3, 4, 5, 6, 6/ INTEGER LSFRAME(10) /10,10, 5, 6, 7, 8, 7, 8,16,90/ ! LRS status INTEGER LSBYTE(10) / 1, 2, 2, 2, 1, 1, 2, 2, 1, 1/ INTEGER HSFRAME(10) / 1, 1, 2, 4, 8,16,16,32,48,64/ ! HRS status INTEGER HSPACKET(10)/ 8, 8, 8, 8, 8, 8, 8, 8, 8, 8/ INTEGER HSBYTE(10) / 5, 6, 6, 6, 6, 5, 6, 6, 6, 6/ INTEGER LDFRAME(2) /66,67/ ! Data channels, LS 10-bits of bytes 1 & 2 INTEGER IGAIN(4) /2,4,3,1/ ! Gain state of (status bits + 1) INTEGER ICHOP(4) /2,1,4,3/ ! Chopper mode of (status bits + 1) C Reference: GLL 3-280 Tables A2.8.1-4 C_CALLS U_MOVE1 (isislib), NIMS_EXPAND4 (nimslib) C_HIST 26feb93 R Mehlman, UCLA/IGPP ORIGINAL VERSION C 08mar93 RMX Added gain state and chopper mode parameters, C corrected some LRS status item locations. C_END BYTE ITEMP1(5), IHW INTEGER*2 ITEMP2(4) DO I=1,7 ! Command buffers CMD(I,1) = LHSK( LCBYTE(I), LCFRAME(I)+1) CMD(I,2) = HHSK( HCBYTE(I), HCPACKET(I)+1, HCFRAME(I)+1) ENDDO DO I=1,10 ! Status items HSKSTAT(I,1) = LHSK( LSBYTE(I), LSFRAME(I)+1) HSKSTAT(I,2) = HHSK( HSBYTE(I), HSPACKET(I)+1, HSFRAME(I)+1) ENDDO DO I=1,2 ! LRS data channels CALL U_MOVE1( 2, LHSK(1, LDFRAME(I)+1), ITEMP1(4)) ! expand CALL NIMS_EXPAND4(ITEMP1, ITEMP2) ! 10-bit DCHAN(I) = ITEMP2(4) ! values ENDDO DO J=1,2 IHW = HSKSTAT(1,J) ! Hardware status I = (IHW.AND.'30'X)/16 IGS(J) = IGAIN(I+1) ! Gain state I = IHW.AND.'03'X ICHM(J) = ICHOP(I+1) ! Chopper mode END DO RETURN END SUBROUTINE NIMS_OHMS2TMP(ICHNL,OHMS,TCEN,TKEL,IR) C_TITLE NIMS_OHMS2TMP Compute NIMS P/F temperatures from measured resistances C_ARGS INTEGER ICHNL !I Channel number (see DESC below) REAL OHMS !I Resistance for channel ICHNL REAL TCEN !O Temperature (degrees Centigrade) REAL TKEL !O Temperature (degrees Kelvin) INTEGER IR !O Return code: 0 OK, -1 invalid channel C_DESC This routine is for the shield, telescope, grating, chopper, and C electronics. (Use NIMS_FPATMP for focal plane temperatures.) C Uses measured resistance and corresponding Callander - van Dusen C coefficients to figure out temperatures. Channel identification C given below. The identity of the chopper sensor was not recorded C by Bulova, so we use parameters for 15323 (grating mechanism). C NO. LOCATION TERMINALS CDS CHANNEL SENSOR C === ======== ========= =========== ====== C 1 Shield L - M E-1911 15319 C 2 Telescope F - G E-1912 15322 C 3 Grating Mech. D - E E-1913 15323 C 4 Opt. Chopper B - C E-1914 ? C 5 Electronics N - d E-1915 UA75 C_CALLS NIMS_PT_RES_THERM C_HIST 25-AUG-88 R Carlson JPL Adapted from INSTMP C 01-FEB-93 R Mehlman UCLA/IGPP Renamed from OHMS2TMP, renamed C subroutine, error parameter added, internal C documentation put in ISIS format. C_END DIMENSION RZ(5), ALPHA(5), DELTA(5), BETA(5) ! COEFFICIENTS DATA RZ/ 497.793, & 498.235, & 498.434, & 498.434, & 500.089/ DATA ALPHA/0.0038896874, & 0.0038918771, & 0.0038922303, & 0.0038922303, & 0.0038830074/ DATA DELTA/1.386378, & 1.403585, & 1.396086, & 1.396086, & 1.386402/ DATA BETA/ 0.2247264, & 0.2131825, & 0.2143358, & 0.2143358, & 0.2327845/ 100 I = ICHNL IF(I.LT.(1)) GO TO 900 ! Out of range IF(I.GT.(5)) GO TO 900 ! Out of range R = OHMS CALL NIMS_PT_RES_THERM(R,RZ(I),ALPHA(I),DELTA(I),BETA(I), 1 TCEN,TKEL) IR = 0 RETURN 900 IR = -1 RETURN END SUBROUTINE NIMS_PT_RES_THERM(R,RZ,ALPHA,DELTA,BETA,TCEN,TKEL) C_TITLE NIMS_PT_RES_THERM Compute Pt resistance thermometer temperatures C_ARGS REAL R !I Measured resistance REAL RZ !I Resistance at 0 degrees C REAL ALPHA !I | Constants of REAL DELTA !I | equation described REAL BETA !I | in DESC below REAL TCEN !O Temperature in degrees Centigrade REAL TKEL !O Temperature in degrees Kelvin C_DESC Computes NIMS Platinum resistance thermometer temperatures from C measured resistance. Uses Callandar - van Dusen equation, measured C value of resistance (R), and previously determined constants C (RZ - resistance at 0 deg C, ALPHA, DELTA, BETA) to determine C temperature in Centigrade and Kelvin (TCEN,TKEL). For TAU = TCEN/100, C the empirical Callandar - van Dusen equation is: C R/RZ = 1 + ALPHA[100*TAU - DELTA*TAU*(TAU-1) C - BETA*(TAU-1)*TAU**3] C with BETA = 0 if TCEN > 0. (Note that the value of BETA supplied C to this routine is finite and for TCEN < 0; it is ignored as C appropriate.) The program first solves the quadratic equation C (BETA = 0). For RHO = R/RZ > 1, the resulting TAU gives the C temperatures. For RHO < 1, the BETA term gives a quartic equation, C and the above-derived TAU is used as a starting point to solve C this quartic. C_CALLS U_FIND_ROOT4 C_HIST 24-APR-87 R Carlson JPL ORIGINAL VERSION C 20-FEB-93 R Mehlman UCLA/IGPP Renamed from PRT, renamed U_ROOT4 and C changed argument order, internal documentation put in C ISIS format. C_END RHO = R/RZ A = -ALPHA*DELTA ! Coeff in A*X**2 + B*X + C B = ALPHA*(100.+DELTA) ! Linear coefficient C = 1. - RHO ! C TAU = (-B + SQRT(B**2-4.*A*C))/(2.*A) ! Quadratic sol'n IF(RHO.GE.(1.00)) GO TO 800 ! for TCEN > 0 TAUS = TAU ! Initial approximation for TCEN < 0 D = ALPHA*BETA ! Cubic coefficient E = -D ! Coefficient of TAU**4 EPS = 0.0001 ! Iterate to 0.01 degree CALL U_FIND_ROOT4(E,D,A,B,C,EPS,TAUS,TAU) 800 TCEN = 100.*TAU ! In Centigrade TKEL = TCEN + 273.15 ! Ice point = 0 deg C = 273.15 K RETURN END SUBROUTINE NIMS_UNPK_ERT(IERT,KERT) C_TITLE NIMS_UNPK_ERT Unpack ERT from NIMS EDR into component parts C_ARGS INTEGER*4 IERT !I Packed ERT: bits 0-6 year-1900 C 7-15 day-of-year C 16-31 minute-of-day INTEGER*4 KERT(4) !O Unpacked ERT: (1) year-1900 C (2) day-of-year C (3) hour-of-day C (4) minute-of-hour C_HIST 28aug91 R. Mehlman UCLA/IGPP ORIGINAL VERSION C 02nov92 RMX Renamed from U_UNPK_ERT. C_END INTEGER*4 JERT INTEGER*2 JMIN,JYD(2) BYTE LERT(4) EQUIVALENCE (JERT,LERT,JYD),(JMIN,LERT(3)) JERT = IERT KERT(1) = IAND(JERT,'007F'X) ! year KERT(2) = ISHFT(JYD(2),-7) ! day KERT(3) = JYD(1)/60 KERT(4) = JYD(1) - 60*KERT(3) ! minute RETURN END SUBROUTINE U_FILL1(IVAL,N,IA) C_TITLE U_FILL1 Fill byte array with constant value C_ARGS BYTE IVAL ! I Value to be inserted INTEGER*4 N ! I Length of array in bytes BYTE IA(*) ! O Array to be filled C_KEYS UTILITY, BYTE, L1, ARRAY C_HIST 09feb93 RMX PORTABLE VERSION of ISIS routine C_END IF (N.LE.0) RETURN DO I=1,N IA(I)=IVAL END DO RETURN END SUBROUTINE U_FIND_ROOT3(A3,A2,A1,A0,EPS,XTRY,X) C_TITLE U_FIND_ROOT3 Find root of cubic equation by Newton's method C_ARGS REAL A3 !I Coefficients of cubic equation REAL A2 !I REAL A1 !I REAL A0 !I REAL EPS !I Rough accuracy REAL XTRY !I First guess as to where root is REAL X !O Result C_DESC Find root of cubic equation: 0 = A3*X**3 + A2*X**2 + A1*X + A0 C using Newton's method. C_HIST ??-JUL-85 R Carlson JPL ORIGINAL VERSION C 01-FEB-93 R Mehlman UCLA/IGPP Renamed from ROOT3, argument order C and internal documentation to ISIS standard. C_END XNEW = XTRY 100 XOLD = XNEW Y = A3*XOLD**3 + A2*XOLD**2 + A1*XOLD + A0 YP = 3.*A3*XOLD**2 + 2.*A2*XOLD + A1 XNEW = XOLD - (Y/YP) IF(ABS(XNEW-XOLD).GT.EPS) GO TO 100 X = XNEW RETURN END SUBROUTINE U_FIND_ROOT4(A4,A3,A2,A1,A0,EPS,XTRY,X) C_TITLE U_FIND_ROOT4 Find root of quartic equation by Newton's method C_ARGS REAL A4 !I Coefficients of quartic equation REAL A3 !I REAL A2 !I REAL A1 !I REAL A0 !I REAL EPS !I Rough accuracy REAL XTRY !I First guess as to where root is REAL X !O Result C_DESC Find root of quartic equation: C 0 = A4*X**4 + A3*X**3 + A2*X**2 + A1*X + A0 C using Newton's method. C_HIST 22-APR-87 R Carlson JPL Adapted from ROOT3. C 01-FEB-93 R Mehlman UCLA/IGPP Renamed from ROOT4, argument order C and internal documentation to ISIS standard. C_END XNEW = XTRY 100 XOLD = XNEW Y = A4*XOLD**4 + A3*XOLD**3 + A2*XOLD**2 + A1*XOLD + A0 YP = 4.*A4*XOLD**3 + 3.*A3*XOLD**2 + 2.*A2*XOLD + A1 XNEW = XOLD - (Y/YP) IF(ABS(XNEW-XOLD).GT.EPS) GO TO 100 X = XNEW RETURN END SUBROUTINE U_MOVE1(N,IA,JA) C_TITLE U_MOVE1 Move byte array C_ARGS INTEGER*4 N ! I Length of array in bytes BYTE IA(*) ! I Source array BYTE JA(*) ! O Destination array C_DESC Move bytes from source array to destination array. Source C and destination must not overlap C_KEYS UTILITY, BYTE, L1, ARRAY, SYSTEM, C_HIST 16oct92 RMX PORTABLE VERSION of ISIS routine C_END DO I=1,N JA(I)=IA(I) END DO RETURN END