C
C   ~crb/ing/signal.f
C
C   Chris Benn  La Palma  crb@ing.iac.es
C
C   Predicts photon counts from target and from sky, and resulting 
C   signal-to-noise, for spectroscopic and direct-imaging observations.
C   The program is self-contained - quantitative information about
C   the wavelength bands, sky, telescopes, instruments, gratings and 
C   detectors is provided in the /data statements.  
C   There are no input data files.
C
C   This version is for the instruments of the Isaac Newton Group, La Palma,
C   but only the data statements are site- or instrument-specific.
C
C   If you are viewing this code with a web browser, you may
C   find chunks of concatenated lines near the end of the program.
C   This is triggered by html 'pre' tags embedded in write statements (used
C   when SIGNAL is run via the web interface), but if the code is saved
C   to a file, it should look normal.
C
C   Wavelength ranges of UBVRIZ bandpasses are given starting ~ line 800
C   Numbered list of instruments starts ~ line 1580 - search for '((true'
C
C   Prominent reminders for editor (crb): 
C   (1) Special flag in code IP = 6 for standalone version, 8 for 
C   version to run with web interface.
C   (2) IW = 0 usually, otherwise writes summary information about parameters
C   to file outsum.
C
C==============================================================================
C
C   COMPILING AND LINKING
C   ---------------------
C   To compile and link under unix, use:
C      f77 signal.f -o signal.exe
C   (option -e needed if modifications extend any lines beyond 80 characters).
C   The program is self-contained, no libraries or data files needed.
C   No graphics (those provided by the web interface are made at a higher
C   level, after repeated calls to signal.f to get the data for plotting).
C   If you have any problems running this code on your machine, please let me 
C   know.
C   The program is written in Fortran 77, with the following widely-used
C   extensions (as flagged by f77 -ansi): IMPLICIT NONE, DO...ENDDO,
C   subroutine names > 6 characters or including '_',
C   non-standard 'name=' (why?), non-standard comment (separating data
C   statement and data?).  Code checked with ftnchek (for date of last
C   check, see history section of comments).
C
C   If you obtained the executable by ftp, you may need 
C      chmod +x signal
C   to reset the protection on the file, before running it.
C
C
C   WHAT THE PROGRAM DOES
C   ---------------------
C   The program calculates the number of photons/A from the object
C   (SOBJ) and photons/A/arcsec**2 from the sky (SSKY) as =
C      10**(-mag/2.5)                 where mag is OBJM or SKYMAG(IBAND)
C                                     OBJM is magnitude for a point object,
C                                     mag/arcsec**2 for extended object
C     *TIME                           integration time in sec
C     *PHOTONS(IBAND)*10000           conversion mag -> 
C                                     photons/m**2/sec/A
C                                     (factor 10000 is because PHOTONS
C                                     defined /cm**2 not /m**2)
C     *EXTIN(IBAND)**AIRMASS          zenith transmission by atmosphere
C                                     (implemented from v11.0 onwards)
C     *TELAREA(ITEL(IINSTR))          telescope area (m**2)
C     *EFFMIRROR**IMIRROR(IINSTR)     efficiency of main mirrors
C     *VIGN                           fraction of light passed by slit,
C                                     fibre or lesnlet (spectroscopic obs
C                                     of point sources only)
C     *EFFINSTR(IBAND,IINSTR)         throughput of instrument (including
C                                     derotator e.g Nasmyth and corrector
C                                     e.g. PF, if any, and excluding other
C                                     components listed here)
C     *EFFGRAT(IBAND,IGRAT,IINSTR)    efficiency of grating
C     *EFFDET(IBAND,IDET)             q. efficiency of detector
C     *TRUETHRU(IBAND,IINSTR)         ratio empirical/theoretical throughput,
C                                     now determined for most ING instruments,
C                                     typically 0.7 - 1.1, other values 
C                                     extrapolated or set to 1.0 (all were
C                                     1.0 for pre-1995 versions of SIGNAL).
C
C   If IPX=0 (point-source calculation), OBJM is the apparent mag,
C   if IPX=1 (surface-brightness calc.), OBJM is apparent mag/arcsec**2
C   ITYPE = 1 for an imaging calculation, 2 for spectroscopic,
C   3 for integral-field or fibre spectroscopy (in which case IPX is
C   forced to be 0).
C
C   (1) IMAGING - POINT SOURCE (ITYPE=1, IPX=0)
C   Number of photons from object (SOBJ2) and sky/pixel (SSKY2) are 
C   calculated as:
C      SOBJ2=SOBJ*BAND(IBAND)
C      SSKY2=SSKY*BAND(IBAND)*ARCPIX*ARCPIX
C   And S:N ratio for an aperture of diameter 2 FWHM is:
C      SNR=SOBJ2/SQRT(SOBJ2+SEE*(SSKY2+READ**2)) 
C   where SEE is pi*(FWHM/ARCPIX)**2 pixels (number of pixels over 
C   which integrating)
C
C   (2) IMAGING - EXTENDED SOURCE (ITYPE=1, IPX=1)
C   Number of photons from object/pixel (SOBJ2) and sky/pixel (SSKY2) are 
C   calculated as:
C      SOBJ2=SOBJ*BAND(IBAND)*ARCPIX*ARCPIX
C      SSKY2=SSKY*BAND(IBAND)*ARCPIX*ARCPIX
C   And S:N ratio per pixel is:
C      SNR=SOBJ2/SQRT(SOBJ2+SSKY2+READ**2) 
C   i.e. SEE is assumed to be 1 pixel
C
C   (3) SPECTROSCOPY - POINT SOURCE (ITYPE=2, IPX=0)
C   Number of photons per pixel-step in wavelength (SOBJ2) and 
C   number of photons/pixel from the sky (SSKY2) are calculated as:
C      SOBJ2=SOBJ*DISPPIX
C      SSKY2=SSKY*SLIT*ARCPIX*DISPPIX
C   And S:N ratio per pixel-step in wavelength (i.e. per pixel of
C   extracted spectrum) is:
C      SNR=SOBJ2/SQRT(SOBJ2+SEE*(SSKY2+READ**2))
C   where SEE is assumed to be 2*(FWHM/ARCPIX)
C
C   (4) SPECTROSCOPY - EXTENDED SOURCE (ITYPE=2, IPX=1)
C   Number of photons/pixel from object from the object (SOBJ2) and 
C   photons/pixel from the sky are calculated as:
C      SOBJ2=SOBJ*SLIT*ARCPIX*DISPPIX
C      SSKY2=SSKY*SLIT*ARCPIX*DISPPIX
C   And S:N ratio per pixel is:
C      SNR=SOBJ2/SQRT(SOBJ2+SSKY2+READ**2)
C   i.e. SEE is assumed to be 1 pixel
C
C   (5) INTEGRAL-FIELD / FIBRE SPECTROSCOPY - POINT SOURCE 
C                                               (ITYPE=3, IPX=0)
C   Number of photons per pixel-step in wavelength (SOBJ2) and 
C   number of photons/pixel from the sky (SSKY2) are calculated as:
C      SOBJ2=SOBJ*DISPPIX
C      SSKY2=SSKY*SLIT*SLIT*3.14159/4.*DISPPIX 
C                 (i.e. phot/pixel-step, not per pixel)
C   And S:N ratio per pixel-step in wavelength (i.e. per pixel of
C   extracted spectrum, for one lenslet or fibre) is:
C      SNR=SOBJ2/SQRT(SOBJ2+SSKY2+SEE*READ**2)
C   where SEE is assumed to be SPREAD pixels (spread of light orthogonal
C   to dispersion direction)
C
C   (6) INTEGRAL-FIELD / FIBRE SPECTROSCOPY - EXTENDED SOURCE 
C                                               (ITYPE=3, IPX=1)
C   Calculation made as above, except that object mag is assumed
C   to be mag within solid angle subtended by circular aperture
C   of diameter SLIT.
C
C   The program also calculates the scales at the detector in arcsec/pixel
C   (ARCPIX) and A/pixel (DISPPIX) as:
C      ARCPIX=ARCMM(IINSTR)     *PIXEL(IDET)
C     DISPPIX=DISP(IGRAT,IINSTR)*PIXEL(IDET)
C
C   Where no empirical data are available, SIGNAL uses interpolated,
C   extrapolated or (as a last resort) theoretical throughputs.
C   Empirical throughputs are typically 80 - 90% of theoretical.
C   NB some of the instruments and detectors documented by SIGNAL are not
C   currently offered by ING (e.g. WHT FOS, JKT RBS, IPCS).
C
C   The program loops, at each loop either obtaining the parameters
C   interactively from the user, via the menu, or reading them 
C   non-interactively from the data file SIGNAL.INPUT (one loop only
C   if this file has been generated from the web interface).
C
C==============================================================================
C
C   LIMITATIONS (some of this is summarised in subroutine helpsub)
C   -----------
C   The program does not take into account the following:
C
C   - Variation of atmospheric (non-grey) and dust (grey) extinction.
C     The default values are from La Palma technical note 31. 
C     La Palma suffers from Saharan dust in the summer, extinction typically
C     0.2 - 0.4 mag, up to 1 mag.  The dust is grey but the extinction
C     varies from night to night (see the ING web page for
C     a nightly record from the Carlsberg Meridian Circle).
C
C   - Degradation of the main telescope mirrors between aluminisings.  This may
C     be up to 10%, probably more in the UV.
C
C   - Vignetting and optical aberrations (broadening the image) at large
C     field radius (see the ING Observers' Guide for more information).  
C
C   - Wavelength-dependence of vignetting by the spectrograph slit, 
C     if slit is not at parallactic angle.
C     (Vignetting due to slit, fibre or lenslet aperture was included
C     from v13.1 onwards).  
C
C   - Loss of light to colour, ND and polarisation filters in the light path.
C     But SIGNAL assumes that order-sorting filters HAVE been inserted where
C     necessary (e.g. to exclude second-order light in ISIS red arm!).
C
C   - Loss of grating efficiency at large angles, or when the grating is
C     reversed or overfilled.
C
C   - Conversion between recorded counts from CCD and number of detected
C     photons (i.e. the effects of gain = e-/ADU, bias level, dark current,
C     non-linearity of response, cosmic rays etc.)
C     Quoted CCD readout noise and e-/ADU are for the commonest readout speed
C     for the instrument in question.
C
C   - S:N depends on the brightness of the sky (specified by the user),
C     and this can less easily be predicted than e.g. the throughput of the
C     system from atmosphere to CCD.  Some of the contributions to the
C     brightness of the night sky are:
C     - twilight (if sun < 18 deg below horizon);
C     - moonlight (if moon not below horizon);
C     - solar cycle (sky brightens up to 0.4 mag at solar maximum);
C     - zodiacal light (sky is 0.4 mag brighter on ecliptic than at poles);
C     - Milky Way;
C     - airmass (sky is 0.2 mag brighter at airmass 1.5 than at zenith);
C     - light pollution - e.g. the sky is significantly brighter 
C       when very dusty (> 0.3 mag extinction)
C     See ING technical note 115 for further information.
C
C   - Signal-to-noise for spectroscopy of unresolved objects is given
C     per pixel step in wavelength, which may not be the same as the
C     wavelength resolution (typically 2 pixels).
C
C   Notes specific to ING instruments, filters and detectors:
C
C   - ISIS polarisation optics have throughput 0.9 approx, ISIS dichroics 
C     about 0.85 
C
C   - For FOS, U,B are for second order, V,R,I are for first order
C
C   - The BVRI filters are assumed to be Harris.  Note that the throughput
C     of the KPNO B filter is about 25% greater than Harris
C     The U filter is the 50-mm type.  The 125-mm WHT PF U filter has a 
C     transmission about a factor of two lower.
C
C==============================================================================
C
C   HISTORY
C   -------
C   1990 Nov 26  Version 1.0 in VAX Fortran
C                Includes non-interactive option 'X'
C   1991 Jun 11  Version 2.0
C   1992 Feb  8  Corrected erroneous I, Z CCD sensitivities
C   1992 Mar 12  Added detector QE for Tek CCD (version 3.0)
C   1992 Mar 25  Updated GEC/EEV CCD responses with PRJ's latest data
C                (version 3.1)
C   1992 Apr  3  Added LDSS spectroscopic and imaging 
C                Generalised code to be instrument-independent (version 4.0)
C   1993 Oct 14  Updated UES efficiencies using Max Pettini's data
C   1993 Dec  9  Added INT FORD and GEC7 QE 
C   1994 Sep  6  Version 6.0: added WHT prime (instr 19), removed JKT QUBES
C                and grism options; added WHT WYFFOS and hitchhiker;
C                updated readout noises and default CCD;
C                corrected IPCS2 pixel size to 10 microns;
C                added derotator eff 0.75 for UES;
C                added R1200B grating (assume eff as R600B)
C   1995 Jan 23  Added TRUETHRU empirical/predicted data (+ extrapolations)
C                WHT FOS throughput based on 1/1/92 measurement
C                Version 7.0 (-> UK/NL/ES site managers 6/3/95)
C   1995 Mar 24  Removed incorrect warning that no throughput data for WYFFOS
C                Added correct DISP for WYFFOS
C   1995 May 10  Removed tabs.  Version 7.1
C   1995 May 16  Corrected WHT PF TRUETHRU(B) from 0.4 to 0.6.  Version 7.2
C   1995 Nov 14  Revised (slightly upwards) aux-port TRUETHRU on basis of
C                throughput measurements of night 13/11/95.  Version 7.3
C   1995 Nov 24  Added table showing status of empirical data
C   1996 Aug 14  Broke up data statements to avoid more than 19 
C                continuation lines (Fortran 77 standard, adhered to by
C                some UNIX machines)
C                Version 8.0
C   1996 Sep 20  Removed carriagecontrol='list' which makes Linux choke
C                Added Loral QE
C                Version 9.0 (emailed to ING/RGO webmans + UK/NL/ES
C                site managers)
C   1996 Oct  1  LORAL readnoise=8, version 9.1
C   1996 Oct 10  Added help comment about B(KP) vs B(Harris)
C   1996 Nov 11  Updated sky-brightness values (version 9.2)
C   1996 Dec 17  Added ISIS R grating efficiencies in Z (version 9.3)
C   1997 Mar  3  Corrected WHT PF U effinstr to 1.0->0.6, truethru 0.7->1.1 
C                (version 9.4)
C        Mar 14  Added WYFFOS empirical throughput data from 11/96 comm run,
C                provisonal result that 0.5 mag poorer than ISIS (version 10)
C        Dec  2  Added EEV10 CCD
C                Removed dependence of ADU on INSTRUMENT
C                Fine-tuned sky-brightness values
C   1998 Feb 13  Structural changes and enhancements not affecting 
C                object and sky count rates for main instruments:
C                - Added calculations of S:N for spectrsocopy and for
C                  extended images both imaging and spectroscopy
C                - Option to set effective bandwidth, useful if using
C                  a filter with very different transmission from
C                  default set (e.g. WHT PF U) or narrow band
C                - Added JHK bands, WHIRCAM InSb detector, WHIRCAM instrument
C                - Added STJ detector, human eye, scrapped IPCS and redundant
C                  CCDs
C                - FWHM now refers to fwhm of star only, IPX takes over
C                  function of flagging whether calculation is for star
C                  or per arcsec**2 for extended object.
C                - Added calculation of detector sizes in both directions.
C                - Reorganised to make addition of new telescopes, 
C                  instruments etc. easier.  Some variable names changed.
C                  Extended arrays to allow 10 bands instead of 6 
C                  and 50 instruments instead of 20.
C                - Made instrument-menu 2-column
C                - Help simplified, most additional info now in comments
C                  in code
C                - Updated readout-noise and gain values, now reports both
C                - Corrected error in WHT FOS V throughput
C                - Added extinction-correction code, but set to 0.0
C                - Made program check command issued by user
C                - Added option to change telescope area e.g. for JOSE,
C                  Shack-H, which use sub-apertures
C                (version 10.2, checked ISIS + INT-imaging calcs give same
C                counts as for previous versions)
C   1998 Feb 17  Changes which DO affect object and sky count rates:
C                - Implemented extinction correction
C                - Updated sky-brightness values (small changes)
C                - Updated effective filter bandwidths BAND(I) for Harris
C                  filters and new U filters (still waiting for Harris I scan)
C                  NB big change in U bandwidth, prob. offset by inclusion of
C                  atmospheric ext. in U, so TRUETHRU not much different
C                  for imaging.
C                - Changes to theoretical throughput values:
C                - UES throughput changed 0.36 to 0.45 (new derotator)
C                  ISIS throughput changed 0.45 to 0.51 or 0.42 (two arms)
C                  WYFFOS throughput changed 0.36 to 0.51 (error)
C                  WHT PF B throughput changed 1.0 to 0.8 (PF corrector)
C                  IDS500 throughput changed 0.41 to 0.44 (not same as 235)
C                - Multiplied INT PF U, B th. throughputs by 0.6, 0.8 because
C                  of corrector.
C                - Updated table of empirical throughputs and adjusted
C                  TRUETHRU values accordingly.
C                - Changed area over which sky flux counted from pi/4 FWHM**2
C                  to pi FWHM**2 for imaging calculations (Hans Deeg IAC 
C                  pointed out an inconsistency in previous versions, gave
C                  correct counts from star and from sky, but an overoptimistic
C                  estimate of resulting S:N for imaging, by integrating
C                  sky over aperture of diameter FWHM, which doesn't include 
C                  all the flux)
C                (version 11.0)
C   1998 Feb 20  - Changed to display band with wavelength when 
C                  a spectroscopic instr. (version 11.1)
C   1998 Feb 21  - Changed to accept Oke AB mags too, altered some of 
C                  TRUETHRU to reflect difference between apparent mag
C                  and Oke AB mag (version 11.2)
C   1998 Apr 14  - Calculated object counts, for imaging of extended source,
C                  corrected to be per pixel instead of per arcsec^2
C                  (error noted by Nigel Metcalfe)
C                - Two coding errors, accepted by SUN f77, rejected by
C                  Richard McMahon's gnu g77 compiler:
C                  (1) if(itype(iinstr.eq.1)) print *,... corrected to
C                      if(itype(iinstr).eq.1) print *,...
C                  (2) lines ending in comma, followed by lines starting
C                      with comma, modified so only one comma
C                      (even f77 -ansi doesn't catch this)
C                - IROT data statement accidentally commented out, 
C                  90-deg rotation of EEV10a not being used
C                  (error noted by Sally Oey)
C                - Minor changes to help + output format
C                (version 11.3)                  
C   1998 Apr 20  Cleared up a few minor inconsistencies flagged by ftnchek
C                running on lpss20.  Version 11.4
C   1998 Jun  3  Added CIRSI (INT or WHT).  Version 11.5
C   1998 Jul  3  Added default detector numbers (IDEFDET) for each
C                instrument.  Added INGRID.  Version 11.6.
C   1998 Jul 19  Tweaked sky-brightness values to match those given in 
C                LP tech note 115 (Benn + Ellison).  Version 11.7.
C   1998 Aug  5  Adding sky-brightness calculator - estimates sky brightness
C                in each band with or without moonlight.  Version 11.8.
C                Checked SIGNAL's output for ISIS and INT WFC agree previous
C                versions.  Checked with ftnchek.
C   1998 Sep  9  Adding lines at start of code to allow it to be used
C                in one of two modes after minor editing:
C                IP = 6 - usual interactive mode, writes to screen (unit = 6)
C                IP = 7 - reads parameters from signal.input,
C                         writes to signal.output (unit=7),
C                         which allows it be called from a web browser,
C                         via the cgi (perl) script written by
C                         Ashley James (summer student 1998)
C                Moved all parameter-setting code into two subroutines,
C                one for interactive, one for non-interactive.
C                Non interactive mode can also be triggered from menu, with
C                undocumented option X (untested).
C                Replaced 'print *' common to both modes by writes 
C                (variable put free format in BUFFA, BUFFB)
C                Ruthlessly eliminated all data type logical
C                Version 12.0.
C   1998 Sep 16  Minor changes to some format statements (error picked
C                up by Paul O'Brien's compiler).  FTNCHEK OK.  Version 12.1.
C   1999 Mar  2  Updating data on Loral CCD so that it refers to the one in 
C                INT FOS rather than to the defunct WHT ones.  Version 12.2.
C   1999 Mar 29  Correcting large-EEV read noise 11 -> 4 e-.  Version 12.3.
C   2000 Jun 21  Added some IDS gratings for WYFFOS.  Version 12.4.
C   2000 Sep  4  Removed print * formatting of output values from subroutine
C                free.
C                Changed K SKYMAG (i.e. sky + thermal) to 12.5 (for INGRID),
C                corrected INGRID gain to 5.6 e/ADU, QE changed from
C                estimate of 0.6 to .4, .55, .55, corrected ZJHK bandpasses
C                to be for INGRID.
C                Version 12.5
C   2000 Sep 13  New INGRID and UES throughput measurements, version 12.6.
C   2001 Feb 12  New INGRID Ks skybr 11.8 mag/arcsec^2 from Johan Knapen's
C                measurements 4-12/00.
C                Rearranged instrument menu.
C                Added NAOMI/INGRID.
C                Added comment line to appear after o/p table for each instr.
C                Version 12.7
C   2001 Jul 19  New INGRID gain 5.3 (was 5.5), readnoise 32 (was 10),
c                scale 0.238 "/pixel,
C                J,H,K skybr values (fainter) and zeropoints (factor 1.6
C                fainter, largely due to factor 1.4 increased throughput
C                for the new aspheric)
C                Version 12.8
C   2001 Oct 19  Commented K bgr high for NAOMI/INGRID
C                Version 12.9
C   2002 Feb  5  Added ITYPE(IINSTR) = 3 option for ifu/fibre spectroscopy,
C                  so that now get sky counts and S:N for e.g. AF2/WYFFOS.
C                Added OASIS and MIT/LL CCD.
C                Added variable DEFSLIT - default slit width (or fibre
C                  diam etc) for given instruments.
C                Corrected mismatches between instr names in .pl script
C                  and html (including removal of multiple spaces, which
C                  are not treated in a straightforward way).
C                Option to allow setting sky-br as D G or B.
C                Now outputs sky photons/A/sq.arcsec as well (McMahon's
C                  suggestion 10/98).
C                Updated wyffos copy of data for ISIS gratings
C                  (redward extrapolation had been done only for ISIS).
C                Version 13.0
C   2002 Feb  7  Added correction for vignetting by slit, fibre or
C                lenslet (VIGN).
C                Corrected parameter-order errors in sky-brightness
C                  .pl script (noted by Ian Skillen 6/01).
C                Converted CBAND* and QBAND to 2-character, for future
C                  implementation of additional filters.
C                Added variable SPREAD - number of pixels over which light
C                  distributed orthognal to dispersion direction, for
C                  fibre and lenslet spectrographs (2 for wyffos, 3 for 
C                  oasis).
C                Version 13.1
C   2002 Mar  2  Small changes to default dark, grey and bright sky 
C                intensities.
C                Version 13.2
C   2003 Jun  4  Added NOT ALFOSC imaging + spectroscopy test
C   2003 Jun  5  Changed number of detectors from 10 to 30, version 13.4
C                Removed detector options: INT FOS Loral, old GEC
C                Added: Marconi 2148*4700, NOT CCD7 (Loral)
C   2003 Jun 10  Adding NOT NOTCAM, MOSCA, version 13.5
C                Removed superfluous print of readnoise to output text
C                Changed lengths and sizes of some arrays for
C                  greater flexibility
C                Version 13.5
C   2003 Aug 25  Robert Greimel added mode ip=8, for use with his
C                new .php web interface.
C                Version 13.6
C   2003 Dec  3  Added calc of peak counts per pixel for imaging of
C                point sources (sugg by Almudena Zurita)
C                Version 13.7
C   2003 Dec 30  Adding LIRIS
C                Version 13.8
C   2009 Nov 18  Version 13.9 (which may be same as 13.8), checked that
C                it works when copied into 
C                /export/localhome/god/signal/current/cgi-bin
C   2009 Nov 19  Version 14.0:
C                Added ACAM imaging and spectroscopy (instr = 35, 36)
C                Marked aux-port as N/A
C                Added detectors AUXCAM, RED+, EEV12, LIRIS, EEV10a, WFC (INT) 
C                (det = 14 - 19)
C                Converted detector 2 from 'EEV' to 'WHTWFC',
C                  renamed 'Rockw' as 'INGRID'
C                Made default det on startup = RED+
C                Updated AF2/WYFFOS throughput values to match 2004
C                  small-fibres values on AF2 web page, changed 
C                  dispmm from short-camera to long-camera values,
C                  made default detector WHTWFC.
C                Added LIRIS grisms
C                Added section to write instrument/detector
C                  data summary to outsum (unit 11),
C                  defined IW, IWD, IWI, IWB
C                Defined IDEFDET2 as second default detector (basically
C                  for above, to allow for fact that 2 detectors
C                  used with ISIS)
C                Current text block output to user is 79 chars wide.
C
C   2009 Dec 16  Version 14.1:
C                Added OASIS imaging (instr 37).
C                Added calculation of *predicted* zeropoint ZP
C                  for each instrument,
C                  and recorded measured zeropoints (previously only
C                  in comment lines) in ZMEAS, ZMEASC, ZMEASD 
C                  (these last 3 vars not used in any calcs)
C                Addressed problem that IPX always set to 0 if ITYPE=3,
C                  with a fudge: at the start of the main loop,
C                  the mag/arcsec^2 is converted to a mag within 
C                  the circular aperture of diameter SLIT, and then
C                  treated as a single object, but the numbers displayed 
C                  to the user are IPX = 1 and mag/arcsec^2 (NB VIGN
C                  is forced to be 1).
C                
C   2009 Dec 17  Version 14.2:
C                Added default gratings IDEFGRAT, IDEFGRAT2, to
C                  be used with IDEFDET and IDEFDET2
C                Changed LIRIS observed zeropoints (the values
C                  quoted on the web pages are for 1 ph/sec, not
C                  1 ct/sec).  LIRIS now publicly available, but
C                  only with warning 'Approx...'
C
C   2010 Jan 21  Version 14.3:
C                Extended Atm/tel/instr throughput line beyond 80 chars
C                  with calc and reporting of thru(4)
C                Moved reporting section to subroutine WRITEOUT,
C                  put subroutines in alphabetic order
C                Corrected ACAM zeropoints, esp U
C                Corrected IDEFDET for INT WFC to 19 not 2.
C
C   2013 Aug 13  Version 14.4:
C                Added RED+4 detector (sl diff QE from other RED+), 
C                ndet now = 21.
C                Revised (by factors ~ 2 in B V R I) the estimated
C                throughputs effinstr for AF2/WYFFOS, to be consistent
C                with Dave King's analysis of 12/2011.
C                Updated the AF2/WYFFOS truethru to be consistent with
C                the new effinstr and with the new on-sky zeropoints 
C                measured by Cecilia and Lilian (logged in ZMEASD)
C                during 2013.  Now ~ 1 in R I Z; ~0.5 in B, V; ~ 0.04 in U
C                (but throughput changes fast with wavelength in UV,
C                so averaging over a broad band has limited meaning).
C                Added extra decimal place to output value of truethru.
C                Reviewed and updated comment lines.
C                Changed ZMEASD(19) from '(4/2009 Carter' to '(9/1994 Carter)' 
C                
C   2013 May 6   Version 14.5:
C                Updated AF2 zeropoints, see details in zeropoints section.
C
C   Old values of numerical parameters are often preserved as comment
C   lines (in the data block).
C
C
C==============================================================================
C
C   SUBROUTINES
C   -----------
C   display
C   free
C   getpar
C   help
C   helpsky
C   pause
C   pause2
C   slitvign
C   skybrcalc
C   skybrcalc2
C   welcome
C   writeout
C
C   
C==============================================================================
C
C   VARIABLES DEFINING INSTRUMENTS, DETECTORS, BANDS ETC.
C   -----------------------------------------------------
C   Data for new telescopes, instruments, gratings, detectors and bands
C   can be added by modifying the DATA statements (only).
C   In the table below are listed the principal data arrays, with 
C   the array dimensions currently being:
C                                            B (bands)       = 10
C                                            G (gratings)    = 20
C                                            I (instruments) = 50
C                                            D (detectors)   = 30
C   (e.g. 'G I' = one number required for each grating for each instrument).
C
C    Variable   Indices    Meaning (see full list above for more detail)
C
C    BAND       B          effective bandwidths of filters i.e. integral 
C                          over T(l)dl, where T(l) is the transmission and 
C                          l is wavelength = area under filter transmission 
C                          curve
C    CBAND      B          band name (uppercase)
C    CBANDLOWER B          band name (lowercase)
C    EXTIN      B          extinction in mag, at airmass = 1.0
C    FLUXJY     B          Jy flux (1e-26 W/Hz/m**2) for mag=0 in each band 
C                          (Z is approx)
C    SKYMAG     B          sky brightness in mag/arcsec^2
C    WAVE       B          wavelength in microns
C                          (used in v11.0 and after, when calculating 
C                          PHOTONS from FLUXJY)
C
C    ARCMM          I      "/mm at detector
C    CINSTR         I      instrument name
C    COMMENT        I      1-line comment for instrument (printed by signal)
C    DEFSLIT        I      default slit width arcsec
C    IAXIS          I      normal direction on CCD of spectral dispersion 
C                          (1 = x, 2 = y)
C    IDEFDET        I      default detector number
C    IDEFDET2       I      second default detector number (usu 0)
C                          (used when writing to unit 11)
C    IDEFGRAT       I      default grating number (usu 1)
C    IDEFGRAT2      I      second default detector number (usu 0)
C                          (used when writing to unit 11)
C    IMIRROR        I      number of aluminium tel mirrors between sky and 
C                          instrument (the throughput of any non-aluminium 
C                          tel mirrors, e.g. the silvered A&G flat, is 
C                          included in EFFINSTR)
C    ITEL           I      telescope number for the given IINSTR (to get tel 
C                          geom area)
C                          1 = JKT     2 = INT     3 = WHT    4 = EYE   5 = NOT
C                          6 = GTC
C    ITRANS         I      order in which instrument menu written
C    ITYPE          I      1 = imager, 2 = spec, 3 = fibre or IFU spec
C    NGRATING       I      number of gratings 
C    SPREAD         I      number of pixels over which spectrum spread perp
C                          to wavelength, for fibre and lenslet spec
C                          (typically 2 - 3).
C                          Note that this assumes always same detector.
C
C    ADU              D    gain e-/count
C    CDET             D    detector name
C    IROT             D    rotation of each detector: 0 means spectral 
C                          direction would be in X if mounted on ISIS, 
C                          1 means 90 deg rotated
C                          (for ING instr only)
C    NPIXELS          D    no. active pixels in x and y on detector
C                          (second dimension of arrays is 2, for x and y)
C    PIXEL            D    pixel size on detector (mm)
C    READ             D    read noise (e-)
C
C    CGRATING     G I      grating name (G,I)
C    DISPMM       G I      dispersion A/mm (G,I)
C    EFFDET     B     D    detector QE (B,D)
C    EFFGRAT    B G I      grating throughput (B,G,I)
C    EFFINSTR   B   I      instrument throughput (B,I)
C    TRUETHRU   B   I      empirical/predicted instrument throughput -
C                          i.e. a fudge factor to make prediction agree with
C                          observed zeropoint (B,I)
C    ZMEAS      B   I      Measured zeropoint for each instrument
C                          = mag for 1 det ph/sec for an imaging instrument
C                          = mag for 1 det ph/sec/A for a spec instrument
C                          Third index of array is default det number (usu 1)
C                          ZMEAS is not currently used in any calculation,
C                          just for reporting to file outsum
C    ZMEASC         I      Configuration (grating/det etc.) used for 
C                          above ZMEAS
C    ZMEASD         I      Date and name of person who made measurement
C
C
C    CTEL                  names of telescopes
C    TELAREA               geometric collecting area of telescope 
C
C    NBAND, NINSTR, NGRATING(IINSTR), NDET 
C                          numbers of bands, 
C                          instruments, gratings and detectors
C
C  NBAND, NGRATING, NINSTR and NDET may need changing.
C  If the detector arrays are changed, IDEFDET may need changing.
C  If the grating arrays are changed, IDEFGRAT may need changing.
C  To add data involving a new telescope, may need to modify some or all
C  of IMIRROR, ITEL, TELAREA.
C  ITRANS will probably need changing to allow correct ordering of options.
C
C
C   OTHER VARIABLES
C   ---------------
C     AIRMASS   airmass
C     ANS       1-char answer
C     ARCPIX    scale "/pixel at detector (derived)
C     AREA      geometrical collecting area * mirror, instrument, grating
C               and detector efficiencies
C     BUFFA, BUFFB, BUFFC, BUFFD  
C               15-char buffers to which reals are written free format,
C               prior to output
C               CGRATING is (IGRAT,IINSTR)
C     CBANDLOWER
C     CINSTRB   buffer to hold instr menu 
C     DISPPIX   dispersion A/pixel at detector (derived)
C     DX, DY    sizes of detector in X and Y in " or A (whichever is
C               appropriate)
C     LOOP      main loop counter
C     EFFMIRROR reflective efficiency assumed for aluminised
C               telescope mirror 
C     FWHM      FWHM arcsec (combined effects of seeing and intinsic
C               object FWHM)
C     I,J,K,KK  loop counters
C     IBAND, IINSTR, IGRAT, IDET, ITEL - band, instrument, grating, detector 
C               and telescope numbers
C     IDX,IDY   axis labels: 0 = arcsec, 1 = A
C     IID,IIG   default det and grating numbers, used as temporary
C               variables in loop calculating assuming zeropoints
C     IIN,IOUT
C     IP        6  => interactive (writes to unit 6)
C               7  => non-interactive reads from SIGNAL.DAT, writes to 
C                     SIGNAL.OUT (unit 7)
C               8  => mode added by Robert Greimel, to allow signal.f to 
C                     be used with web interface
C     IPX       calculation for stellar object (0) or per sq arcsec for
C               extended image (1)
C     IPX0      Remembers IPX during main calculation loop
C     IW        = 1 to write out instrument or detector summary
C               (Invoked by editing the code)
C     IWB       Number of bands to write out e.g. 9
C     IWD       Detector numbers for which summary required
C     IWI       Instrument numbers for which summary required
C     LOOP      Loop counter for main calculation loop
C     MODE      mode of use after parameters obtained:
C               0 = S:N calc; 1 = sky br calc; 2 = help; 3 = skybr help;
C               9 = end
C     OBJM      magnitude of object
C               (per sq arcsec for spectroscopy with IPX=1)
C     OBJM0     Remembers OBJM during main calculation loop
C     PEAK      peak counts / pixel for imaging = SOBJ/1.14/FWHM^2
C               (i.e. assumes gaussian profile)
C     PHOTONS   photons/s/sq.cm/A for mag = 0 for each band (derived)
C     SEE       3.14159 * (FWHM/PIXSIZE)**2 i.e. number of pixels over which
C               2D image is assumed to be spread (was 3.14159/4. * ... before
C               version 11.0), affects S:N ratio only
C     SKYMAG0   set to SKYMAG at start, but SKYMAG0 cannot be changed
C               by the user (the values are passed to subr SKYBRCALC)
C     SLIT, SLIT2  slit-width (arcsec)
C     SNR       signal-to-noise for imaging calculation
C     SOBJ      number of photons/A from object
C               (per sq arcsec for spectroscopy with IPX=1)
C     SOBJ2     For spectroscopy, photons/pixel (sky or ext object with IPX=1)
C               or photons/pixel step in wavelength i.e. after integrating
C               in the spatial direction (IPX = 0)
C     SSKY      number of photons/arcsec**2/A from sky
C     SSKY2     SSKY per pixel for imaging or for normal spectroscopy
C               SSKY per pixel step in wavelength for IFU/fibre spec.
C     THRU(4)   throughputs of atm, tel, instr, grating, used in displaying
C               results to user
C     TIME      integration time (sec)
C     TRUE0     value of TRUETHRU passed to DISPLAY
C     UNIT(2)   " or A (for writing out pixels scales in the two directions)
C     VERSION   Version number and date
C     VIGN
C     ZP        Predicted zeropoint for each instrument, calculated from
C               EFFINSTR, EFFGRAT, EFFDET etc.
C               Same array size as ZMEAS, above.
C               Not currently used to calculate anything else, it's
C               just reported to output file outsum.
C
C     Variables for passing to DISPLAY:
C     QINST,QBAND,QGRAT,QGRAT2,QDET,EXT,
C     GEOMAREA,THROUGHPUT,PIXSIZE,DETEFF,
C     PHOTONS2,READOUT,BWIDTH,SKY
C     (the ones starting with Q are char strings)
C
C   Variables used by sky-brightness calculator SKYBRCALC only:
C     SOL,ECLAT,GLAT,AIR,IM,ALPHA,R,ZM,AK
C
C   Last checked above list against variable declarations below: 23/11/09
C  
C==============================================================================
C
      implicit none
      integer i,j,k,kk,ninstr,ndet,iband,iinstr,igrat,ip,iin,iout,mode
      integer itrans(50),itel(50),itype(50),imirror(50),ngrating(50)
      integer iaxis(50),npixels(30,2),idefdet(50),idefdet2(50)
      integer idet,nband,loop,ipx,idx,idy,iw,iwd(30),iwi(50),irot(30)
      integer iwb,ipx0,idefgrat(50),idefgrat2(50),iid,iig
      real defslit(50),spread(50),vign,thru(4)
      real effinstr(10,50),effgrat(10,20,50),effdet(10,30)
      real zmeas(10,50,5),zp(10,50,5)
      real truethru(10,50),dispmm(20,50),arcmm(50),telarea(50)
      real objm,time,snr,slit,arcpix,disppix,effmirror,area,sobj,ssky
      real sky,sobj2,ssky2,see,slit2,true0,fwhm,geomarea,ext,airmass
      real throughput,pixsize,deteff,photons2,bwidth,readout,dx,dy
      real pixel(30),read(30),adu(30),peak,objm0
      real band(10),fluxjy(10),photons(10),skymag(10),skymag0(10),
     1 extin(10),wave(10)
      character*1 ans
      character*2 cband(10),cbandlower(10),qband
      character*1 unit(2)
      character*6  qdet
      character*10 qgrat,qgrat2,ctel(100)
      character*15 buffa,buffb,buffc,buffd
      character*20 cgrating(20,50)
      character*25 cinstr(50),qinst
      character*28 version
      character*40 cdet(30)
      character*80 comment(50),zmeasc(50),zmeasd(50)

C
C   Variables used by sky-brightness calculator (only)
C
      integer im
      real sol,eclat,glat,air,alpha,r,zm,ak

C==============================================================================

      data unit/'"','A'/
C
C   Version number and date
C
      data version/'Version 14.5      6 May 2014'/

C  
C   Lists of instruments and detectors to be summarised
C   on outsum (unit 11), if IW = 1
C  
      data iwi/1,35,36,19,16,33,34,26,27,37,9,10,38*0/
      data iwd/16,15,4,20,14,2,3,17,10,19,18,1,9,17*0/
      data iwb/9/
C
C
C##############################################################################
C   
C   DATA SECTION
C
C   The data statements below define the properties of the sky 
C   (La Palma only), telescopes, instruments and detectors.
C   They are organised in 4 sections, in the same order as the variables
C   list above, i.e.:
C     - Those depending on band only
C     - Those depending on instrument only
C     - Those depending on detector only
C     - Those depending on some combination of band, instrument, 
C       grating and detector.
C   Most of the instrument/detector data (and the measured zeropoints,
C   documented in the comments) are taken from the relevant
C   ING web pages.  Many of the predicted instrument throughputs 
C   EFFINSTR are a product of the throughputs of indivdual 
C   components in the optical train, and where possible, these are
C   documented.
C
C
C   Number of instruments, detectors and bands on offer
C
      data ninstr/37/
      data ndet  /21/
      data nband / 9/
      data ctel/'JKT','INT','WHT','Eye','NOT','GTC',94*' '/

      data effmirror/0.85/

      data telarea/0.66,4.41,12.47,8.e-5,4.6,74.8,44*0.0/
C     (approx est of area for GTC, itel=6, as circle diam 10.4 m,
C     minus hexagon of side 1.4 m)

C
C==============================================================================
C
C   VARIABLES DEPENDING ON BAND ONLY 
C   --------------------------------
C
C   Central wavelengths of bands used (from output of ~crb/ing/filter.f).
C   Typical wavelength ranges are given below.  Actual values of effective
C   bandwidth (integral over T(l)dl, where T is transmission, l is wavelength)
C   are given in array BAND
C                                         BAND values are for:
C      U  3600 A  (3300 - 3800 A)         50-mm square and 139-mm liq-CuSO4 U
C      B  4300 A  (4000 - 5000 A)         Harris
C      V  5500 A  (5000 - 6000 A)         Harris
C      R  6500 A  (6000 - 7200 A)         Harris
C      I  8200 A  (7200 - 9200 A)         Harris
C      Z  9500 A  (8500 - 9700 A)         RGO glass
C      J 12500 A  (1.10 - 1.40 mic)       WHIRCAM
C      H 16700 A  (1.50 - 1.82 mic)       WHIRCAM 
C      K 21600 A  (2.00 - 2.30 mic)       WHIRCAM Kshort
C
C   Note that the central wavelengths and wavelength range are different
C   for different types of filter.  E.g. the NOT ALFOSC filters have:
C      Centre   FWHM Thruput
C    U  3620 A   600  60%
C    B  4400 A  1000  68%
C    V  5300 A   800  85%
C    R  6500 A  1300  81%
C    i  7970 A  1570  91%  
C

      data cband     
     1/'U ','B ','V ','R ','I ','Z ','J ','H ','Ks','? '/
      data cbandlower
     1/'u ','b ','v ','r ','i ','z ','j ','h ','ks','? '/

      data wave/3600.,4300.,5500.,6500.,8200.,9500.,
     1  12500.,16700.,21600.,0.0/
C
C   Effective bandwidths UBVRI taken from INT prime-focus guide for
C   versions up to 10.2, and by integrating over UBVR filter curves with
C   ~crb/ing/filter.f for signal version 11.0 onward.
C   NB still waiting for reliable scans of the I and Z filters
C
      data band
C    1 / 260., 780., 675.,1000.,1400.,1200.,2070.,2300. ,2800. ,0./ ! <  10.2
C    1 / 500., 720., 860.,1330.,1400.,1200.,2070.,2300. ,2800. ,0./ ! >= 11.0
     1 / 500., 720., 860.,1330.,1400., 700.,1600.,2900. ,3200. ,0./ ! >= 12.4
C
C   Mag -> Jy conversions from Bessell (1979 PASP 91 589)
C   and Bessell & Brett (1988 PASP 100 1134)
C   (prior to v11.0, were from Johnson H.L., 1966, ARAA, 4, 193)
C
      data fluxjy
C    1 /1880.,4440.,3810.,3010.,2430.,2200.,1570.,1020. , 640. ,0./ ! < 10.2
     1 /1810.,4260.,3640.,3080.,2550.,2200.,1570.,1020. , 640. ,0./ ! >=11.0
C
C
C   Photon rate for mag=0 (not used from v11.0 onward, since it can be
C   calculated from FLUXJY) taken from INT prime-focus manual
C     data photons
C    1 / 783.,1463.,1050., 873., 376., 329., 194.,  94.4,  44.1,0./
C   New values are similar +-10% except R now 701, I now 449)

C
C   Sky-brightness values are median dark-of-moon ING 1994-6 (solar
C   minimum) from Sara Ellison's summer 1996 work (for UBVRI;
C   Z is an interpolation between these and the known JHK values),
C   as reported in LP technical note 115.
C   K is for INGRID (Thomas Augusteijn email 7/8/00, previous values
C   were for WHIRCAM)
C   Sky brightness is similar (+- 0.2 mag) at other dark sites.
C
      data skymag/22.0,22.7,21.9,21.0,20.0,18.8,16.1,14.7,12.5,0.0/
C                 22.0,22.7,21.9,21.0,20.0,18.8,16.6,14.4,11.8,0.0 to v12.7
C                 22.0,22.7,21.9,21.0,20.0,18.8,16.6,14.4,12.5,0.0 to v12.6 
C                 22.0,22.7,21.9,21.0,20.0,18.8,16.6,14.4,12.0,0.0 to v12.4
C                 22.0,22.7,21.7,20.9,20.0,18.8,16.6,14.4,12.0,0.0 to v11.6
C                 21.4,22.3,21.4,20.4,19.3,19.2 used until version 9.1
C                 22.0 23.0 21.9 21.2 20.2      (La Silla, ESO manual)

C
C   Zenith extinction for La Palma, taken from the ING Observers' Guide
C   = constant (wavelength-dependent) gas + variable (grey) dust (= 0.05)
C   Used for version 11.0 onwards.
C   The values for other sites are similar.  For CTIO, Frogel 1998,
C   PASP 110 200 gives .11, .06 and .09 in J, H and K.
C
      data extin/0.55,0.25,0.15,0.09,0.06,0.05,0.1,0.11,0.07,0.0/

C
C==============================================================================
C
C   VARIABLES DEPENDING ON INSTRUMENT ONLY 
C   --------------------------------------
      data arcmm/
     1  14.9,35.0,23.25,12.4,23.25, 12.4,15.1,4.51,24.7,29.4, ! 1 - 10
     1  13.8,58.0,13.8,49.47,0.,    0.,24.6,24.6,17.5,13.6,   ! 11 - 20
     1  8.0,20000.,17.5,24.7,12.87, 2.16,0.0,12.6,12.6,12.76, ! 21 - 30
     1  4.32,7.33,13.51,13.51,16.7, 16.7,1.33, ! 31-37
     1     13*0.0/

      data (cinstr(i),i=1,10)/'WHT ISIS                 ',
     1  'WHT FOS             (N/A)',
     1  'WHT TAURUS imag f/2 (N/A)',
     1  'WHT TAURUS imag f/4 (N/A)',
     1  'WHT TAURUS etal f/2 (N/A)',
     1  'WHT TAURUS etal f/4 (N/A)',
     1  'WHT UES             (N/A)',
     1  'WHT aux-port imaging(N/A)',
     1  'INT wide-field camera    ',
     1  'INT IDS 235 camera       '/
      data (cinstr(i),i=11,50)/
     1  'INT IDS 500 camera  (N/A)',
     1  'INT FOS             (N/A)',
     1  'JKT imaging         (N/A)',
     1  'JKT RBS             (N/A)',
     1  'JKT QUBES           (N/A)',
     1  'WHT AF2/WYFFOS           ',
     1  'WHT LDSS imaging    (N/A)',
     1  'WHT LDSS spectrosc  (N/A)',
     1  'WHT prime-focus imaging  ',
     1  'WHT Hitchhiker      (N/A)',
     1  'WHT WHIRCAM         (N/A)',
     1  'Human eye (approx.)      ',
     1  'WHT CIRSI           (N/A)',
     1  'INT CIRSI           (N/A)',
     1  'WHT INGRID               ',
     1  'WHT NAOMI/INGRID         ',
     1  'WHT NAOMI/OASIS spec.    ',
     1  'NOT ALFOSC imaging       ',
     1  'NOT ALFOSC spectroscopy  ',
     1  'NOT NOTCAM IR imaging WF ',
     1  'NOT NOTCAM IR imaging HR ',
     1  'NOT MOSCA imaging        ',
     1  'WHT LIRIS imaging        ',
     1  'WHT LIRIS spectroscopy   ',
     1  'WHT ACAM imaging         ',
     1  'WHT ACAM spectroscopy    ',
     1  'WHT NAOMI/OASIS imaging ',
     1     13*' '/

      data comment/24*'                                  ',
     1 '   (NB assumes read noise due to one IR-chip readout only)', ! 25
     1 '   (NB K bgr in NAOMI/INGRID is about 10 mag/arcsec^2)', ! 26
     1 '   (NB test version only - use official OASIS S:N calculator)', ! 27
     1 '   (NB for ALFOSC ignore effective bandwidth quoted above)', ! 28
     1 '   ', ! 29
     1 '   (NB for NOTCAM ignore effective bandwidth quoted above)', ! 30
     1 '   (NB for NOTCAM ignore effective bandwidth quoted above)', ! 31
     1 '   (NB for MOSCA ignore effective bandwidth quoted above)', ! 32
     1 '   (NB approx only - use official LIRIS S:N calculator)', ! 33
     1 '   (NB approx only - use official LIRIS S:N calculator)', ! 34
     1 '   ', ! 35
     1 '   ', ! 36
     1 '   (NB test version only - use official OASIS S:N calculator)', ! 37
     1     13*'                                     '/

C       1 2 3 4 5 6 7 8 9  0  1 2 3 4 5 6 7 8 9 0 1 2 3 4 5 6 7
      data iaxis  
     1 /1,2,0,0,0, 0,2,0,0,1, 1,2,0,1,1, 2,0,1,0,0, 
     1  0,0,0,0,0, 0,2,0,1,0, 0,0,0,0,0, 2,0,13*0/ 

      data idefdet
     1 /15,6,1,1,1, 1,5,1,19,18, 18,4,1,1,1, 2,5,5,2,6, 7,9,10,10,10, 
     1  10,3,6,6,11, 11,13,17,17,14, 14,3,13*0/

      data idefdet2/16,14*0,2,34*0/

      data idefgrat/15*1,3,34*1/

      data idefgrat2/5,14*0,7,34*0/

C                   1   2   3   4   5   6   7   8   9   0

      data defslit/1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,0.0,1.0,
     1             1.0,1.0,1.0,1.0,1.0,1.6,1.0,1.0,1.0,0.0,
     1             1.0,1.0,1.0,1.0,1.0,0.0,0.2,1.0,1.0,1.0,
     1             1.0,1.0,0.0,1.0,0.0,1.0,0.0,1.0,1.0,1.0,
     1     10*1.0/             

C
C   IMIRROR - NB WHT Nasmyth flat is aluminium, A&G for ACAM is silvered
C
      data imirror
     1 /2,2,2,2,2, 2,3,3,1,2, 2,2,2,2,2, 1,2,2,1,3,
     1  3,0,1,1,3, 3,3,2,2,2, 2,2,2,2,2, 2,3,13*0/

      data itel   
     1 /3,3,3,3,3, 3,3,3,2,2, 2,2,1,1,1, 3,3,3,3,3, 
     1  3,4,3,2,3, 3,3,5,5,5, 5,5,3,3,3, 3,3,13*0/

      data itrans/
     1     19,35,25,26,33,37,9,
     1     1,36,16,27,34,8,2,20,23,24,
     1     17,18,3,4,5,6,7,21,12,10,11,
     1     13,14,15,22,28,29,30,31,
     1     32,13*0/

      data itype  
     1 /2,2,1,1,2, 2,2,1,1,2, 2,2,1,2,2, 3,1,2,1,1,
     1  1,1,1,1,1, 1,3,1,2,1, 1,1,1,2,1, 2,1,13*0/

      data ngrating
     1 /9,2,0,0,0, 0,2,0,0,16, 16,2,0,4,0, 19,0,4,0,0,
     1  0,0,0,0,0, 0,3,0,13,0, 0,0,0,5,0, 1,0,13*0/

      data spread /0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,
     1             0.0,0.0,0.0,0.0,0.0,2.0,0.0,0.0,0.0,0.0,
     1             0.0,0.0,0.0,0.0,0.0,0.0,3.0,0.0,0.0,0.0,
     1     20*0.0/             

C
C==============================================================================
C
C   VARIABLES DEPENDING ON DETECTOR ONLY 
C   ------------------------------------
      data adu  /1.7,0.9,1.0,1.2,1.7, 1.04,6.,1.0,1.0,5.3,
     1   3.0,2.8,1.3,0.9,
     1   1.0,1.2,3.6,1.2,1.5,0.0,0.69,9*0.0/
C   old EEV had adu = 1.1

      data cdet(1) /'TEK    ING 1024*1024, 24.0-mic'/
      data cdet(2) /'WHTWFC ING 2048*4096,13.5-mic'/
      data cdet(3) /'MITLL  ING 2048*4096, 15.0-mic'/
      data cdet(4) /'Marcon ING 2148*4700, 13.5-mic'/
      data cdet(5) /'SITE   ING 2048*2048, 24.0-mic'/
      data cdet(6) /'CCD7   NOT 2048*2048  15.0-mic'/
      data cdet(7) /'InSb   ING 256*256,   30.0-mic'/
      data cdet(8) /'STJ    ESA 6*6,       25.0-mic'/
      data cdet(9) /'Eye                          '/
      data cdet(10)/'INGRID ING 1024*1024, 18.5-mic'/
      data cdet(11)/'Rockw  NOT 1024*1024, 18.5-mic, 24 e- rms'/
      data cdet(12)/'Rockw as (11) but read-up-ramps, 14 e- rms'/
      data cdet(13)/'Loral  NOT 4096*4096, 15.0-mic'/
      data cdet(14)/'AUXCAM ING 2048*2048, 15.0-mic'/
      data cdet(15)/'RED+   ING 2048*4096, 15.0-mic'/
      data cdet(16)/'EEV12  ING 2048*4100, 13.5-mic'/
      data cdet(17)/'LIRIS  ING 1024*1024, 18.5-mic'/
      data cdet(18)/'EEV10a ING 2048*4096, 13.5-mic'/
      data cdet(19)/'INTWFC ING 6144*6144, 13.5-mic'/
      data cdet(20)/'QUCAM2/3 ING 1024*1024, 13-mic'/
      data cdet(21)/'RED+4  ING 4096*4096, 15.0-mic'/

C     data cdet(3) /'EEV   ING 1242*1152, 22.5-mic'/

      data irot/0,1,1,0,0, 0,0,0,0,0, 0,0,0,0,1, 1,0,0,0,0,0,9*0/

      data npixels/1024,4096,2048,2148,2048,2048,256,6,6000,1024,
     1   1024,1024,4096,2048,2048,2048,1024,2048,6144,1024,4096,9*0,
     1             1024,4096,4096,4700,2048,2048,256,6,6000,1024,
     1   1024,1024,4096,2048,4096,4100,1024,4096,6144,1024,4112,9*0/
C   old EEV was 1242x1152 pixels, old FOS INT Loral was 512*1024

      data pixel/0.024,0.0135,0.015,0.0135,0.024,
     1           0.015,0.030,0.025,0.003,0.0185,
     1           0.0185,0.0185,0.015,0.015,0.015,0.0135,
     1           0.0185,0.0135,0.0135,0.013,0.015,9*0.0/

      data read /6.,4.,4.,4.,6.,6.,40.,0.,0.,30.,24.0,14.0,8.5,3.2,
     1     3.9,4.8,10.0,4.2,3.5,0.0,4.2,9*0.0/
C   read noise for old EEV was 6

C
C==============================================================================
C
C   VARIABLES DEPENDING ON SOME COMB BAND, INSTRUMENT, GRATING, DETECTOR
C   --------------------------------------------------------------------
C
C
C   *** CGRATING ***
C
      data ((cgrating(i,j),i=1,20),j=1,10)/
     1  '   R158R ','   R316R ','   R600R ',  ! WHT ISIS
     1  '   R1200R',  
     1  '   R158B ','   R300B ','   R600B ','   R1200B',
     1  '   H2400B',
     1  11*'      ',
     1  '1stORD','2ndORD',18*'      ',        ! WHT FOS
     1  20*'      ',                          ! WHT TAU
     1  20*'      ',                          ! WHT TAU
     1  20*'      ',                          ! WHT TAU
     1  20*'      ',                          ! WHT TAU
     1  '       E31 ','       E79 ',18*' ',   ! WHT UES
     1  20*'      ',                          ! WHT aux-port
     1  20*'      ',                          ! INT prime
     1  'R150V ','R300V ','R400B ',
     1  'R400V ','R400R ','R600R ',           !IDS
     1  'R600IR','R632V ','R831R ','R900V ','R1200U',
     1  'R1200B','R1200Y','R1200R','H1800V','H2400B',
     1  4*'      '/                           ! IDS
      data ((cgrating(i,j),i=1,20),j=11,25)/
     1  'R150V ','R300V ','R400B ',
     1  'R400V ','R400R ','R600R ',           
     1  'R600IR','R632V ','R831R ','R900V ','R1200U',
     1  'R1200B','R1200Y','R1200R','H1800V','H2400B',
     1  4*'      ',                           ! IDS
     1  '1stORD','2ndORD',18*'      ',        ! INT FOS
     1  20*'      ',                          ! JKT imaging
     1  ' 300  ',' 600  ','1200  ','2400  ',16*'      ',  ! RBS
     1  20*'      ',                          ! JKT QUBES
     1  '    R158R ','    R316R ','    R600R ','    R1200R',  ! WHT WYFFOS (i.e.
     1  '    R158B ','    R300B ','    R600B ','    R1200B',
     1  '    H2400B', ! copy ISIS)
     1  '    R600IR','    R632V ','    R831R ','    R900V ',
     1  '    R1200U',
     1  '    R1200B','    R1200Y','    R1200R','    H1800V',
     1  '    H2400B',
     1  1*'      ',
     1  20*'      ',                          ! LDSS imaging
     1  'Low   ','MedBlu','MedRed','High  ',16*'      ',  ! LDSS
     1  20*'      ',                          ! 
     1  20*'      ',100*' '/                          !
      data (cgrating(i,26),i=1,20)/20*'      '/
      data (cgrating(i,27),i=1,20)/'LR','MR','HR',17*'      '/
      data (cgrating(i,28),i=1,20)/20*'      '/
      data (cgrating(i,29),i=1,20)/
     1 '       g3','       g4','       g5','       g6',
     1 '       g7','       g8',
     1 '      g10','      g11','      g12','      g14',
     1 '   g9+g10  echelle','  g13+g10  echelle',
     1 'g13+OrdSt  echelle',
     1 7*'      '/
      data ((cgrating(i,j),i=1,20),j=30,33)/80*'      '/
      data  (cgrating(i,34),i=1,20)/'lr_zj','lr_hk',
     1 'hr_j','hr_h','hr_k',15*'      '/
      data  (cgrating(i,35),i=1,20)/20*'      '/
      data  (cgrating(i,36),i=1,20)/'V400',19*'      '/
      data ((cgrating(i,j),i=1,20),j=37,50)/280*'      '/

C
C   *** DISPMM ***
C
      data ((dispmm(i,j),i=1,20),j=1,10)/
     1  121.,62.,33.,17.,120.,64.,33.,17.,8.,11*0.0, ! WHT ISIS
     1  395.,195.,18*0.0,                            ! WHT FOS
     1  20*0.0,                                      ! WHT TAURUS
     1  20*0.0,                                      ! WHT TAURUS
     1  20*0.0,                                      ! WHT TAURUS
     1  20*0.0,                                      ! WHT TAURUS
     1  2.5,2.5,18*0.0,                              ! WHT UES
     1  20*0.0,                                      ! WHT aux-port
     1  20*0.0,                                      ! INT imaging
     1  271.3,138.5,104.5,104.5,104.4,69.8,70.2,66.5, ! IDS235
     1  50.7, 46.4, 35.3, 35.3, 35.2,34.8,23.2,17.5,4*0.0/
      data ((dispmm(i,j),i=1,20),j=11,25)/
     1  132.2, 66.1, 49.6, 49.7, 49.9,33.3,33.4,31.6, ! IDS500
     1  23.9, 22.1, 16.7, 16.7, 16.4,15.5,10.5, 8.1,4*0.0,
     1  486.,245.,18*0.0,                            ! INT FOS
     1  20*0.0,                                      ! JKT imaging
     1  160.,80.,40.,19.,16*0.0,                     ! JKT RBS
     1  1.8,1.5,1.28,1.12,1.0,0.9,0.82,0.75,0.7,11*0.0,! JKT QUBES
C     1  479.,245.,131.,67.,479.,254.,131.,67.,33.,
C     1  141.,134.,102.,94.,71.,71.,71.,71.,47.,35.,
C     1  1*0.0,                                       ! WHT WYFFOS 
C    Above values are for short camera, pre 2005??
     1  229.,114.,61.,31.,223.,118.,61.,31.,15.3,
     1  66.,62.,47.,44.,33.,33.,33.,33.,21.9,16.3,
     1  1*0.0,                                       ! WHT WYFFOS
C    Above numbers from Romano Corradi's email 2/3/05  long camera
     1  20*0.0,                                      ! WHT LDSS im.
     1  470.,220.,220.,100.,16*0.0,                  ! WHT LDSS
     1  20*0.0,                                      ! WHT prime
     1  20*0.0,                                      ! 20 Hitchhiker
     1  100*0.0/
      data (dispmm(i,26),i=1,20)/20*0.0/
      data (dispmm(i,27),i=1,20)/316.,147.,67.,17*0.0/ ! 27 WHT NAOMI/OASIS
      data (dispmm(i,28),i=1,20)/20*0.0/               ! 28 NOT ALFOSC im
      data (dispmm(i,29),i=1,20)/153.,200.,207.,100.,100.,87.,
     1 433.,320.,867.,93.,36.7,33.3,33.3,7*0.0/ ! 29 NOT ALFOSC spec
      data ((dispmm(i,j),i=1,20),j=30,33)/80*0.0/
      data (dispmm(i,34),i=1,20)
     1       /330.,524.,99.,132.,182.,15*0.0/ ! 34 WHT LIRIS spec
      data (dispmm(i,35),i=1,20)/0.0,0.0,18*0.0/ ! 35 WHT ACAM imag
      data (dispmm(i,36),i=1,20)/220.0,0.0,18*0.0/ ! 36 WHT ACAM spec
      data (dispmm(i,36),i=1,20)/20*0.0/ ! 37 OASIS imag
      data ((dispmm(i,j),i=1,20),j=38,50)/260*0.0/

C
C   *** EFFDET ***
C
C   Detector quantum efficiences etc are from ING web page or articles
C   in SPECTRUM.  Readout noise and gain are for typical readout
C   speeds, usually FAST.  For WHIRCAM, readout noise is for ND-STARE.
C
      data effdet /
     1  0.30,0.60,0.65,0.68,0.45,0.25,0.00,0.00,0.00,0.00, !  1 TEK
     1  0.55,0.82,0.80,0.74,0.45,0.20,0.00,0.00,0.00,0.00, !  2 WHTWFC mosaic
     1  0.20,0.65,0.80,0.93,0.95,0.56,0.00,0.00,0.00,0.00, !  3 MIT/LL (BIV)
     1  0.48,0.67,0.90,0.90,0.70,0.30,0.00,0.00,0.00,0.00, !  4 Marconi
     1  0.30,0.60,0.65,0.68,0.45,0.25,0.00,0.00,0.00,0.00, !  5 SITE
     1  0.60,0.90,0.93,0.85,0.50,0.00,0.00,0.00,0.00,0.00, !  6 Loral (NOT)
     1  0.00,0.00,0.00,0.00,0.00,0.00,0.80,0.80,0.80,0.00, !  7 InSb (WHIRCAM)
     1  0.60,0.60,0.60,0.60,0.60,0.60,0.60,0.60,0.60,0.00, !  8 STJ
     1  0.00,1.00,1.00,1.00,0.00,0.00,0.00,0.00,0.00,0.00, !  9 Eye
     1  0.00,0.00,0.00,0.00,0.00,0.40,0.40,0.55,0.55,0.00, ! 10 INGRID (also
C                                                                OK for CIRSI)
     1  0.00,0.00,0.00,0.00,0.00,0.40,0.40,0.55,0.55,0.00, ! 11 NOT Hawaii
     1  0.00,0.00,0.00,0.00,0.00,0.40,0.40,0.55,0.55,0.00, ! 12 NOT Hawaii RUR
     1  0.60,0.90,0.93,0.85,0.50,0.00,0.00,0.00,0.00,0.00, ! 13 Lor (MOSCA NOT)
     1  0.35,0.55,0.85,0.90,0.80,0.30,0.00,0.00,0.00,0.00, ! 14 AUXCAM
     1  0.35,0.55,0.85,0.90,0.80,0.30,0.00,0.00,0.00,0.00, ! 15 RED+
     1  0.52,0.82,0.80,0.76,0.45,0.12,0.00,0.00,0.00,0.00, ! 16 EEV12
     1  0.00,0.00,0.00,0.00,0.00,0.80,0.80,0.85,0.90,0.00, ! 17 LIRIS
C Z QE is a guess, based on extrapolation of J, K values
     1  0.55,0.82,0.80,0.74,0.45,0.20,0.00,0.00,0.00,0.00, ! 18 EEV10a
     1  0.55,0.82,0.80,0.74,0.45,0.20,0.00,0.00,0.00,0.00, ! 19 INTWFC
     1  0.10,0.65,0.90,0.89,0.60,0.25,0.00,0.00,0.00,0.00, ! 20 QUCAM2/3 (same)
     1  0.50,0.70,0.90,0.96,0.73,0.56,0.00,0.00,0.00,0.00, ! 21 RED+4
     1  90*0.0
     1 /
C   Old detector QEs:
C     1  0.60,0.81,0.81,0.78,0.58,0.23,0.00,0.00,0.00,0.00, ! Loral (INT FOS)
C     1  0.13,0.14,0.32,0.45,0.35,0.12,0.0,0.0,0.0,0.0,    ! EEVold
C     1  0.81,0.80,0.93,0.93,0.50,0.15,0.0,0.0,0.0,0.0,    ! Loral (WHT)
C     1  0.10,0.20,0.33,0.36,0.28,0.13,0.0,0.0,0.0,0.0,    ! FORD
C     1  0.28,0.40,0.67,0.76,0.61,0.33,0.0,0.0,0.0,0.0,    ! GEC7 (coated)
C     1  0.15,0.60,0.70,0.70,0.55,0.20,0.0,0.0,0.0,0.0,    ! RCA
C     1  0.12,0.11,0.06,0.03,0.00,0.00,0.0,0.0,0.0,0.0,    ! IPCS1
C     1  0.12,0.11,0.06,0.03,0.00,0.00,0.0,0.0,0.0,0.0,    ! IPCS2
C     1  0.13,0.14,0.32,0.45,0.35,0.12,0.00,0.00,0.00,0.00, ! GEC

C
C   *** EFFGRAT ***
C
      data (((effgrat(i,j,k),i=1,10),j=1,20),k=1,9)/
     1  0.0,0.0,0.65,0.7,0.6,0.60,0.0,0.0,0.0,0.0,   ! WHT ISIS 158R
     1  0.0,0.0,0.75,0.7,0.45,0.50,0.0,0.0,0.0,0.0,  ! 316R
     1  0.0,0.0,0.55,0.60,0.5,0.55,0.0,0.0,0.0,0.0,  ! 600R
     1  0.0,0.0,0.45,0.5,0.5,0.50,0.0,0.0,0.0,0.0,   ! 1200R
     1  0.7,0.55,0.4,0.35,0.0,0.0,0.0,0.0,0.0,0.0,   ! 158B
     1  0.7,0.65,0.55,0.42,0.0,0.0,0.0,0.0,0.0,0.0,  ! 300B
     1  0.65,0.6,0.5,0.38,0.0,0.0,0.0,0.0,0.0,0.0,   ! 600B
     1  0.5,0.55,0.6,0.38,0.0,0.0,0.0,0.0,0.0,0.0,   ! 1200B added 13/10/94
     1  0.2,0.4,0.3,0.0,0.0,0.0,0.0,0.0,0.0,0.0,110*0.0, ! 2400B
     1  200*1.0,                    ! WHT FOS
     1  200*1.0,                    ! WHT TAURUS imaging
     1  200*1.0,                    ! WHT TAURUS imaging
     1  200*1.0,                    ! WHT TAURUS + etalon
     1  200*1.0,                    ! WHT TAURUS + etalon
     1  200*0.55,                   ! WHT UES
     1  200*1.0,                    ! WHT aux-port imaging
     1  200*1.0/                    ! INT prime
      data (((effgrat(i,j,k),i=1,10),j=1,20),k=10,10)/
     1  0.4,0.55,0.65,0.55,0.40,0.0,0.0,0.0,0.0,0.0,  !IDS 235
     1  0.4,0.55,0.7,0.65,0.4,0.0,0.0,0.0,0.0,0.0,
     1  0.55,0.45,0.3,0.0,0.0,0.0,0.0,0.0,0.0,0.0,
     1  0.4,0.65,0.65,0.5,0.35,0.0,0.0,0.0,0.0,0.0,
     1  0.0,0.0,0.4,0.55,0.55,0.55,0.0,0.0,0.0,0.0,
     1  0.0,0.0,0.55,0.6,0.55,0.5,0.0,0.0,0.0,0.0,
     1  0.0,0.0,0.0,0.3,0.65,0.7,0.0,0.0,0.0,0.0,
     1  0.4,0.6,0.65,0.6,0.45,0.4,0.0,0.0,0.0,0.0,
     1  0.0,0.0,0.4,0.55,0.65,0.65,0.0,0.0,0.0,0.0,
     1  0.0,0.65,0.65,0.55,0.0,0.0,0.0,0.0,0.0,0.0,
     1  0.7,0.6,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,
     1  0.65,0.6,0.5,0.4,0.0,0.0,0.0,0.0,0.0,0.0,
     1  0.35,0.55,0.65,0.65,0.4,0.4,0.0,0.0,0.0,0.0,
     1  0.0,0.0,0.4,0.6,0.65,0.0,0.0,0.0,0.0,0.0,
     1  0.0 ,0.3 ,0.45,0.5,0.0,0.0,0.0,0.0,0.0,0.0,
     1  0.25,0.55,0.5 ,0.0,0.0,0.0,0.0,0.0,0.0,0.0,40*0.0/
      data (((effgrat(i,j,k),i=1,10),j=1,20),k=11,11)/
     1  0.4,0.55,0.65,0.55,0.40,0.0,0.0,0.0,0.0,0.0,  ! IDS 500
     1  0.4,0.55,0.7,0.65,0.4,0.0,0.0,0.0,0.0,0.0,
     1  0.55,0.45,0.3,0.0,0.0,0.0,0.0,0.0,0.0,0.0,
     1  0.4,0.65,0.65,0.5,0.35,0.0,0.0,0.0,0.0,0.0,
     1  0.0,0.0,0.4,0.55,0.55,0.55,0.0,0.0,0.0,0.0,
     1  0.0,0.0,0.55,0.6,0.55,0.5,0.0,0.0,0.0,0.0,
     1  0.0,0.0,0.0,0.3,0.65,0.7,0.0,0.0,0.0,0.0,
     1  0.4,0.6,0.65,0.6,0.45,0.4,0.0,0.0,0.0,0.0,
     1  0.0,0.0,0.4,0.55,0.65,0.65,0.0,0.0,0.0,0.0,
     1  0.0,0.65,0.65,0.55,0.0,0.0,0.0,0.0,0.0,0.0,
     1  0.7,0.6,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,
     1  0.65,0.6,0.5,0.4,0.0,0.0,0.0,0.0,0.0,0.0,
     1  0.35,0.55,0.65,0.65,0.4,0.4,0.0,0.0,0.0,0.0,
     1  0.0,0.0,0.4,0.6,0.65,0.0,0.0,0.0,0.0,0.0,
     1  0.0,0.3,0.45,0.5,0.0,0.0,0.0,0.0,0.0,0.0,
     1  0.25,0.55,0.5,0.0,0.0,0.0,0.0,0.0,0.0,0.0,40*0.0/
      data (((effgrat(i,j,k),i=1,10),j=1,20),k=12,16)/
     1  200*1.0,                    ! INT FOS
     1  200*1.0,                    ! JKT imaging
     1  1.,1.,1.,1.,1.,1.,0.0,0.0,0.0,0.0,           ! WHTRBS
     1  1.,1.,1.,1.,1.,1.,0.0,0.0,0.0,0.0,
     1  1.,1.,1.,1.,1.,1.,0.0,0.0,0.0,0.0,
     1  1.,1.,1.,1.,1.,1.,0.0,0.0,0.0,0.0,160*0.0,
     1  200*1.0,                    ! JKT QUBES
     1  0.0,0.0,0.65,0.7,0.6,0.6,0.0,0.0,0.0,0.0,   ! WYFFOS (copy ISIS data)
     1  0.0,0.0,0.75,0.7,0.45,0.5,0.0,0.0,0.0,0.0,
     1  0.0,0.0,0.55,0.60,0.5,0.5,0.0,0.0,0.0,0.0,
     1  0.0,0.0,0.45,0.5,0.5,0.5,0.0,0.0,0.0,0.0,
     1  0.7,0.55,0.4,0.35,0.0,0.0,0.0,0.0,0.0,0.0,
     1  0.7,0.65,0.55,0.42,0.0,0.0,0.0,0.0,0.0,0.0,
     1  0.65,0.6,0.5,0.38,0.0,0.0,0.0,0.0,0.0,0.0,
     1  0.65,0.6,0.5,0.38,0.0,0.0,0.0,0.0,0.0,0.0,
     1  0.2,0.4,0.3,0.0,0.0,0.0,0.0,0.0,0.0,0.0,
     1  0.0,0.0,0.0,0.3,0.65,0.7,0.0,0.0,0.0,0.0, ! WYFFOS (copy IDS data)
     1  0.4,0.6,0.65,0.6,0.45,0.4,0.0,0.0,0.0,0.0,
     1  0.0,0.0,0.4,0.55,0.65,0.65,0.0,0.0,0.0,0.0,
     1  0.0,0.65,0.65,0.55,0.0,0.0,0.0,0.0,0.0,0.0,
     1  0.7,0.6,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,
     1  0.65,0.6,0.5,0.4,0.0,0.0,0.0,0.0,0.0,0.0,
     1  0.35,0.55,0.65,0.65,0.4,0.4,0.0,0.0,0.0,0.0,
     1  0.0,0.0,0.4,0.6,0.65,0.0,0.0,0.0,0.0,0.0,
     1  0.0 ,0.3 ,0.45,0.5,0.0,0.0,0.0,0.0,0.0,0.0,
     1  0.25,0.55,0.5 ,0.0,0.0,0.0,0.0,0.0,0.0,0.0,10*0.0/
      data (((effgrat(i,j,k),i=1,10),j=1,20),k=17,20)/
     1  200*1.0,                          ! WHT imaging
     1  0.2,0.6,0.8,0.7,0.6,0.5,0.0,0.0,0.0,0.0,          ! WHT LDSS low
     1  0.4,0.7,0.75,0.65,0.4,0.3,0.0,0.0,0.0,0.0,        ! WHT LDSS blue 
     1  0.0,0.25,0.7,0.75,0.6,0.5,0.0,0.0,0.0,0.0,        ! WHT LDSS red
     1  0.35,0.55,0.6,0.3,0.0,0.0,0.0,0.0,0.0,0.0,160*0.0, ! WHT LDSS high
     1  200*1.0,
     1  200*1.0/
      data (((effgrat(i,j,k),i=1,10),j=1,20),k=21,26)/1200*1.0/
      data ((effgrat(i,j,27),i=1,10),j=1,20)/200*1.0/    ! WHT NAOMI/OASIS
      data ((effgrat(i,j,28),i=1,10),j=1,20)/200*1.0/    ! NOT ALFOSC imaging
      data ((effgrat(i,j,29),i=1,10),j=1,20)/
     1 0.75,0.70,0.51,0.36,0.00,5*0.0,                   ! NOT ALFOSC g3
     1 0.35,0.67,0.73,0.61,0.42,5*0.0,                   ! NOT ALFOSC g4
     1 0.00,0.00,0.62,0.72,0.59,5*0.0,                   ! NOT ALFOSC g5
     1 0.80,0.70,0.49,0.00,0.00,5*0.0,                   ! NOT ALFOSC g6
     1 0.00,0.57,0.68,0.55,0.00,5*0.0,                   ! NOT ALFOSC g7
     1 0.00,0.00,0.00,0.61,0.45,5*0.0,                   ! NOT ALFOSC g8
     1 0.70,0.80,0.59,0.42,0.21,5*0.0,                   ! NOT ALFOSC g10
     1 0.00,0.67,0.81,0.67,0.42,5*0.0,                   ! NOT ALFOSC g11
     1 0.00,0.00,0.43,0.69,0.63,5*0.0,                   ! NOT ALFOSC g12
     1 0.60,0.70,0.49,0.00,0.00,5*0.0,                   ! NOT ALFOSC g14
     1 0.00,0.30,0.26,0.18,0.10,5*0.0,                   ! NOT ALFOSC g3+g10
     1 0.23,0.34,0.30,0.00,0.13,5*0.0,                   ! NOT ALFOSC g13+g10
     1 0.00,0.00,0.38,0.00,0.00,5*0.0,                   ! NOT ALFOSC g13+Ordsort
     1 70*1.0/    ! NOT ALFOSC unused
      data (((effgrat(i,j,k),i=1,10),j=1,20),k=30,33)/800*1.0/
      data ((effgrat(i,j,34),i=1,10),j=1,20)/200*1.0/    ! 34 LIRIS spec
      data ((effgrat(i,j,35),i=1,10),j=1,20)/200*1.0/    ! 35 ACAM imag
      data ((effgrat(i,j,36),i=1,10),j=1,20)
     1 /0.55,0.58,0.80,0.70,0.50,0.00,194*1.0/    ! 36 LIRIS spec
      data ((effgrat(i,j,37),i=1,10),j=1,20)/200*1.0/    ! 37 OASIS imag
      data (((effgrat(i,j,k),i=1,10),j=1,20),k=38,50)/2600*1.0/

C
C   *** EFFINSTR ***
C
C   Predicted instrumental throughputs, taken from
C   the August 1995 Observers' Guide unless otherwise stated.
C
      data ((effinstr(i,j),i=1,10),j=1,10)/
C        U    B    V    R    I    Z    J   H   K
C   ISIS throughput 51% on blue arm, 42% on red arm:
     1  0.51,0.51,0.51,0.42,0.42,0.42,0.0,0.0,0.0,0.0,  ! WHT ISIS
     1  0.30,0.50,0.45,0.70,0.60,0.00,0.0,0.0,0.0,0.0,  ! WHT FOS
C   TAURUS focal-reducer throughput:
     1  0.60,0.60,0.60,0.60,0.60,0.60,0.0,0.0,0.0,0.0,  ! WHT TAURUS IMAGING
     1  0.60,0.60,0.60,0.60,0.60,0.60,0.0,0.0,0.0,0.0,  ! WHT TAURUS IMAGING
C   * 0.8 for etalon:
     1  0.48,0.48,0.48,0.48,0.48,0.48,0.0,0.0,0.0,0.0,  ! WHT TAURUS ETALON
     1  0.48,0.48,0.48,0.48,0.48,0.48,0.0,0.0,0.0,0.0,  ! WHT TAURUS ETALON
C   UES efficiency is 0.86 (cross disperser) * 0.7 (camera) * 0.8 (grating) 
C   = 0.48 * 0.93 (new 1997 short-slit derotator) = 0.45
     1  0.45,0.45,0.45,0.45,0.45,0.45,0.0,0.0,0.0,0.0,  ! WHT UES
C   No optics:
     1  1.00,1.00,1.00,1.00,1.00,1.00,0.0,0.0,0.0,0.0,  ! WHT aux-port imaging
C   INT prime efficiency is 0.98**6 for 6 surfaces of corrector
C   * 0.5 in U, * 0.8 in B for transmission
     1  0.45,0.71,0.89,0.89,0.89,0.89,0.0,0.0,0.0,0.0,  ! INT prime imaging
C   IDS 235 throughput = 0.45 (camera) * 0.9 (collimator):
     1  0.41,0.41,0.41,0.41,0.41,0.41,0.0,0.0,0.0,0.0/  ! INT IDS 235 camera
      data ((effinstr(i,j),i=1,10),j=11,50)/
C   IDS 500 throughput = 0.49 (camera) * 0.9 (collimator):
     1  0.44,0.44,0.44,0.44,0.44,0.44,0.0,0.0,0.0,0.0,  ! INT IDS 500 camera
C   
     1  0.30,0.60,0.30,0.60,0.70,0.00,0.0,0.0,0.0,0.0,  ! INT FOS
C   No optics:
     1  1.00,1.00,1.00,1.00,1.00,1.00,0.0,0.0,0.0,0.0,  ! JKT imaging
C   Unknown:
     1  1.00,1.00,1.00,1.00,1.00,1.00,0.0,0.0,0.0,0.0,  ! JKT RBS
C   Unknown:
     1  1.00,1.00,1.00,1.00,1.00,1.00,0.0,0.0,0.0,0.0,  ! JKT QUBES
C   WYFFOS = 0.7 fibres * 0.92 coll. * 0.95 relay mirror *
C   0.9 camera * 0.94 other components = 0.51,
C   (*0.6, 0.8 in U, B for prime-focus corrector)
C   - these numbers for large fibres, not changed for small fibres
C   The thruput of the fibres is particularly uncertain - 0.7
C   for a 26-m fibre run sounds optimistic?, implies 0.96 per 3-m.
C
C    1  0.31,0.41,0.51,0.51,0.51,0.51,0.0,0.0,0.0,0.0,  ! WHT WYFFOS up to v14.3
C
C   New estimates below for AF2/WYFFOS, excluding atm, tel, grating and CCD,
C   are from Dave King's 12/2011 report. E.g. in V band,
C   0.22 = 0.82 (PF corr) x 0.82 (fibre feed) x 0.41 wyffos x 0.75 geom losses.
C   The 0.41 for wyffos = .87 (collimator) x .87 (decollimator) 
C   x .85 (relay mirror) x .96 (corrector) x .85 (camera fold mirror) 
C   x .85 (camera primary) x .96 (field flattener) x .96 (cryostat window).
C
     1  0.08,0.17,0.22,0.21,0.20,0.20,0.0,0.0,0.0,0.0,  ! WHT WYFFOS from v14.4
C
     1  0.20,0.60,0.75,0.70,0.40,0.10,0.0,0.0,0.0,0.0,  ! LDSS imaging
C
     1  0.20,0.60,0.75,0.70,0.40,0.10,0.0,0.0,0.0,0.0,  ! LDSS + grating
C   WHT PF 4-element field corrector + 2-element ADC:
     1  0.60,0.80,1.00,1.00,1.00,1.00,0.0,0.0,0.0,0.0,  ! WHT prime imaging
C
     1  0.00,1.00,1.00,1.00,1.00,0.00,0.0,0.0,0.0,0.0,  ! WHT hitchhiker
C   From WHIRCAM users' guide v1.0 = derotator+dichroic+optics
C   (throughput is a factor 2 lower for MARTINI):
     1  0.00,0.00,0.00,0.00,0.00,0.00,0.6,0.6,0.6,0.0,  ! WHT WHIRCAM
     1  0.00,0.06,0.06,0.06,0.00,0.00,0.0,0.0,0.0,0.0,  ! Human eye
C   Assuming 6/98 blocking filter 0.5 * PF corrector 0.5 = 0.25
     1  0.00,0.00,0.00,0.00,0.00,0.00,0.25,0.25,0.0,0.0,! WHT CIRSI
     1  0.00,0.00,0.00,0.00,0.00,0.00,0.25,0.25,0.0,0.0,! INT CIRSI
C   From ING web page, giving total QE / det QE = 0.5
     1     5*0.0,0.7,0.7,0.7,0.7,0.0,                   ! WHT INGRID
     1     5*0.0,0.5,0.5,0.5,0.5,0.0,                   ! WHT NAOMI/INGRID
C   NAOMI thruput to INGRID = 0.51 = .85 derot (2004) * .87 OAP1
C   to OAP2 including DM * .99 dichroic * 0.7 INGRID optics
     1  0.0,0.0,0.25,0.28,0.27,0.16,0.00,0.00,0.00,0.00,! WHT NAOMI/OASIS
C    1  0.0,0.0,0.38,0.38,0.38,0.25,0.00,0.00,0.00,0.00,! WHT NAOMI/OASIS
C    (where is this from??)
C   NAOMI peak thruput = 0.28 = .85 derot * .87 OAP1 to OAP2
C   including DM * .99 dichroic * .98 flat
C   to OASIS * .94 focus lens * 0.41 OASIS optics (incl gratings?)
C
     1  0.65,0.71,0.70,0.69,0.67,5*0.0,                 ! NOT ALFOSC imaging
     1  0.65,0.71,0.70,0.69,0.67,5*0.0,                 ! NOT ALFOSC spec.
     1  5*0.0,5*1.0,             ! NOT NOTCAM WF
     1  5*0.0,5*1.0, ! NOT NOTCAM HR
     1  10*1.0, ! NOT MOSCA imaging
     1  10*0.79, ! 33 WHT LIRIS imaging
     1  10*0.79, ! 34 WHT LIRIS spec
     1  0.52,0.64,0.71,0.71,0.70,0.70,0.0,0.0,0.0,0.0,  ! 35 WHT ACAM imag
     1  0.52,0.64,0.71,0.71,0.70,0.70,0.0,0.0,0.0,0.0,  ! 36 WHT ACAM imag
C   The above are Tibor's 20/1/10 predicted throughput values for ACAM
C   as-built (including lack of AR coating lens 5), multiplied by the
C   reflectivity of the A&G flat.  See ACAM imaging web page for values. 
C    1  0.32,0.72,0.78,0.78,0.78,0.78,0.0,0.0,0.0,0.0,  ! ACAM im/spec values
C      used before 20/1/10
     1  0.0,0.0,0.38,0.38,0.38,0.25,0.00,0.00,0.00,0.00,! 37 WHT NAOMI/OASIS imag
C  new est, as above for spec
     1  130*0.0/  

C
C   *** TRUETHRU ***
C
C   The TRUETHRU(BAND,INSTR) values are derived from the instrumental
C   zeropoints (apparent mag for 1 detected photon/sec for imaging,
C   or Oke AB mag for 1 photon/sec/A for spectroscopy) given below.  
C   Where two sets of values given, newest is used.
C   These measurements are likely to be
C   affected at the 10% level by changes in the primary-mirror
C   reflectivity between aluminisings.
C   Unless otherwise stated, the imaging data are through the 50-mm 
C   CuSO4/glass U filters and the BVRI Harris set.
C
C   Dates last checked against instrument web pages:
C     20091128: acam, af2, ids (nothing > 2000), intwfc (nothing on 
C               page), isis red/blue, naomi/ingrid, whtpf
C
C
C                      Z E R O P O I N T S
C  Instrument        U    B    V    R    I    Date  Source of information
C  ----------       ---- ---- ---- ---- ----  ----  ---------------------
C  ISIS
C  ----
C                   18.1 18.2 17.6 18.5 17.4  9708  Smartt   B/158B/LOR UBV
C                                                            R/158R/TEK RI
C                        18.0 17.2            0101  Garcia   B/158B/EEV12
C                             18.0~18.0                      R/158R/TEK4
C                                                            (no dichroic)
      data (zmeas(i,1,1),i=1,10)
     1            /18.3,18.6,17.8,18.3,17.5,5*0.0/        
      data zmeasc(1)/'B/R300B/EEV12(UB) R/R158R/RED+ (VRI):'/
      data zmeasd(1)/'(''for new 5300 dichroic'', date unknown)'/
C
C  WHT ACAM imag.
C  --------------
      data (zmeas(i,35,1),i=1,10)
     1           /24.3,26.0,26.7,26.8,27.1,5*0.0/
      data zmeasc(35)/'AUXCAM/Bessell filters, conv to std bwidth:'/
      data zmeasd(35)/'(6/2009, Benn)'/
C
C  WHT ACAM spec.
C  --------------
      data (zmeas(i,36,1),i=1,10)/3*0.0,18.5,6*0.0/
      data zmeasc(36)/'V400/AUXCAM:'/
      data zmeasd(36)/'(6/2009, Benn)'/
C
C  WHT PF imaging
C  --------------
      data (zmeas(i,19,1),i=1,10)/24.1,26.2,26.4,26.5,26.0,5*0.0/
      data zmeasc(19)/'TEK:'/
      data zmeasd(19)/'(9/1994 Carter)'/
C
C   WHT AF2/WYFFOS (refln mode)
C   --------------
C  Instrument        U    B    V    R    I    Date  Source of information
C  ----------       ---- ---- ---- ---- ----  ----  ---------------------
C                        16.5 16.5 16.5       9610? Bridges lower-limits
C                                                   on throughput with
C                                                   various gratings.
C                                                   H2400B is factor 2 poorer.
C                        16.5 16.8 16.6       0107  Corradi small-fibres R600B
C                   13.1 16.1 16.5 16.6 16.5  0410  Corradi small-fibres,
C   R600B (U B V) and R600R (R I), seeing variable 1 - 1.3 arcsec, numbers
C   scaled to the median fibre.
C
C   New data from peak-of-raster throughput measurements by Cecilia/Lilian,
C   and I've assumed seeing = 1.1" to take into account degradation of 
C   throughput by acquisition errors:   
C   
C     data (zmeas(i,16,1),i=1,10) /12.1,16.1,16.9,17.0,16.5,16.1,4*0.0/ ! 8/13
C     data zmeasc(16)/
C    1 'Small fib / long-cam / WHTWFC R300B (UB) / R316R (VRIZ):'/
C
C   Tweaked 5/14 for consistency with re-evaluated zenith zero points 
C   from Cecilia (new correction for airmass) 
C   + using actual seeing = 1.3" for UBV, 0.9" for VRI:
      data (zmeas(i,16,1),i=1,10) /11.0,16.1,16.9,17.0,16.5,16.1,4*0.0/ ! 5/14
      data zmeasc(16)/
     1 'Small fib / long-cam / WHTWFC R300B (UBV) / R316R (RIZ):'/
      data zmeasd(16)/'(5+2/2013, Farina/Dominguez)'/
C
C  WHT LIRIS imag.
C  ---------------
      data (zmeas(i,33,1),i=1,10) /6*0.0,24.83,25.17,24.55,0.0/
      data zmeasc(33)/'LIRIS:'/
      data zmeasd(33)/'(2003, LIRIS commissioning)'/
C      
C  WHT LIRIS spec.
C  ---------------
C
C  WHT NAOMI/INGRID
C  ----------------
      data (zmeas(i,26,1),i=1,10) /6*0.0,23.5,23.6,22.8,0.0/
      data zmeasc(26)/'INGRID (but gain uncertain): '/
      data zmeasd(26)/'(1/2001, Benn)'/
C
C  WHT OASIS imag.  No zeropoints yet, only predicted thruput
C  WHT OASIS spec.    "
C
C   INT IDS 235
C   -----------
C                   13.8 16.3 16.4 16.2 15.1  97?   WEB      post-96-cleaning
C                                                         R150V, TEK for BVRI
C                                             9411  Harlaftis  R150V TEK for U
C                   13.5 16.4 16.5 16.4 15.2  9810  Garcia  EEV10
C                   14.5 16.4 16.4 16.3 14.8  0010  Garcia  TEK5
C 
      data (zmeas(i,10,1),i=1,10) 
     1             /14.5,16.4,16.4,16.3,14.8,5*0.0/
      data zmeasc(10) /'TEK5:'/
      data zmeasd(10) /'(10/2000, Garcia)'/
C
C   INT PF imaging
C   --------------
C                   23.8 25.5 25.5 25.5 25.0  97?   WEB      Loral no. 4
      data (zmeas(i,9,1),i=1,10)
     1             /23.6,25.7,25.6,25.6,24.9,5*0.0/
      data zmeasc(9)/'INTWFC:'/
      data zmeasd(9)/'(3/2005)'/
C  (above values passed on by Robert Greimel, email 9/3/05)
C
C
C  Old, or non-ING, instruments
C
C  WHT INGRID                                 0107  Almudena Zurita
C                                                   25.20 J, 25.28 H, 24.54 Ks
C                                                   assuming 5.3 e-/count
C  WHT aux-port          26.2 26.1 26.3 25.4  9511  Benn     TEK
C  WHT TAURUS imag.      24.5 25.1 25.5 24.9  9203  GEMINI   f/2, GEC
C  JKT imaging      21.0 23.1 23.3 23.5 23.0  9407  Peletier
C          (checked 21.1 23.1 23.2 23.4 23.0  9707  Benn     TEK4)
C  NOT ALFOSC       23.6 25.7 25.7 25.5 24.9  0306  NOT page Loral
C  NOT MOSCA        24.6 26.2 26.0 25.7 25.2  0306  Telting  Loral mosaic
C  NOTCAM HR + WF                             0306  Telting  Rockwell
C                                                            J 23.8, H 23.8, K 23.3
C
C  WHT FOS                    17.0 18.1 17.2  9203  GEMINI
C  WHT LDSS spec    15.9 18.4 18.6 18.2       9404  Carter   U (high), 
C                                                            B, V (blue),
C                                                            R (red), TEK
C  WHT UES          16.3 17.0 17.0 16.7 15.5  9311  Pettini  old derot, TEK
C                        16.7 16.9 16.9 16.0  0009  Telting  new derot, SIT1
C                   17.0 17.4 17.3 17.0       0101  Telting  new derot, 2EEV
C  WHT Hitchhiker        24.5 25.4 25.5 24.7  9206  GEMINI   GEC
C  WHT WHIRCAM                                9506  WHIRCAM users' guide
C  WHT CIRSI                                  9806  Richard McMahon
C  INT CIRSI                                  9806  Richard McMahon
C  INT FOS          13.7 15.7 15.7 17.0 15.9  9402  Harlaftis  GEC
C  NOT ALFOSC                                 0306  Telting
C
C   No empirical data for TAURUS FP, LDSS imaging, IDS500,
C   RBS or QUBES
C 
C  To obtain apparent mag, add the following numbers to an Oke AB mag:
C
C                     U    B    V    R    I
C                   -0.8  0.2  0.0 -0.2 -0.4     
C  
C   The values of TRUETHRU below have been derived by judicious smoothing,
C   interpolation and extrapolation of the data above.  Changes in the
C   instrument since the measurements were made (e.g. the new derotator 
C   for UES) have been taken into account.
C   In 3 cases, the derived TRUETHRU were significantly LARGER than 1.0.
C   ISIS R 1.7, ISIS I 1.4 are at risk of second-order contamination, 
C   and have somewhat arbitrarily been set to 1.0.  
C   The INGRID JHK values have been left as measured.
C
      data ((truethru(i,j),i=1,10),j=1,10)/
C       U   B   V   R   I   Z   J   H   K   ?
     1 0.8,1.0,0.8,1.0,1.0,1.0,1.0,1.0,1.0,1.0, !  1 ISIS
     1 0.7,0.7,0.7,0.6,0.8,0.8,1.0,1.0,1.0,1.0, !  2 WHT FOS
     1 0.8,0.8,0.8,0.7,0.8,0.8,1.0,1.0,1.0,1.0, !  3 TAURUS imaging f/2
     1 0.8,0.8,0.8,0.7,0.8,0.8,1.0,1.0,1.0,1.0, !  4 TAURUS imaging f/4
     1 0.8,0.8,0.8,0.7,0.8,0.8,1.0,1.0,1.0,1.0, !  5 TAURUS etalon f/2
     1 0.8,0.8,0.8,0.7,0.8,0.8,1.0,1.0,1.0,1.0, !  6 TAURUS etalon f/4
C    1 0.8,0.9,0.7,0.6,0.4,0.5,1.0,1.0,1.0,1.0, !  7 UES to v12.4
     1 0.4,0.4,0.6,0.7,0.7,0.7,1.0,1.0,1.0,1.0, !  7 UES > v12.5
     1 0.6,0.6,0.7,0.7,0.8,0.8,1.0,1.0,1.0,1.0, !  8 aux-port
     1 1.0,0.9,0.7,0.7,0.9,0.7,1.0,1.0,1.0,1.0, !  9 INT prime
     1 0.3,0.7,0.7,0.7,0.7,1.0,1.0,1.0,1.0,1.0/ ! 10 IDS 235
      data ((truethru(i,j),i=1,10),j=11,50)/
     1 0.3,0.7,0.7,0.7,0.7,1.0,1.0,1.0,1.0,1.0, ! 11 IDS 500
     1 0.3,0.7,0.8,0.8,0.7,0.8,1.0,1.0,1.0,1.0, ! 12 INT FOS
     1 0.7,0.7,0.9,0.9,1.0,1.1,1.0,1.0,1.0,1.0, ! 13 JKT imaging
     1 1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0, ! 14 JKT RBS
     1 1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0, ! 15 JKT QUBES
C    1 0.3,0.3,0.3,0.3,0.3,0.3,0.0,0.0,0.0,0.0, ! 16 AF2/WYFFOS large fibres
C    1 0.03,0.11,0.21,0.3,0.6,0.3,0.0,0.0,0.0,0.0, ! 16 AF2/WYFFOS small fibres
C                                                  up to v14.3
C    1 0.04,0.18,0.24,0.29,0.46,0.43,4*0.0,     ! 16 AF2/WYFFOS small fibres
C   (revised values 14/8/13 to be consistent with new on-sky
C   measurements by Cecilia and Lilian (R300B/RED+4, R316R/RED+4)
C   and v14.3 estimates of effinstr for af2/wyffos.  These values
C   not used in any released version of signal.
C                                                    
C    1 0.04,0.43,0.55,0.70,1.20,1.08,4*0.0,     ! 16 AF2/WYFFOS small fibres
C                                                 from v14.4,
C   for consistency with new zeropoints from Cecilia/Lilian 8/2013
C   and with new effinstr from Dave King 12/2011
C
C   Revised 6/5/14, as noted above under zero-points:
     1 0.014,0.37,0.61,0.57,0.96,0.81,4*0.0,     ! 16 AF2/WYFFOS small fibres
C                                                  from v14.5
     1 1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0, ! 17 LDSS imaging
     1 1.2,1.1,0.9,0.7,0.7,0.7,1.0,1.0,1.0,1.0, ! 18 LDSS spectroscopy
     1 0.9,0.9,0.7,0.7,0.7,0.7,1.0,1.0,1.0,1.0, ! 19 WHT prime
     1 0.7,0.7,0.8,0.5,0.5,1.0,1.0,1.0,1.0,1.0, ! 20 Hitchhiker
     1 0.0,0.0,0.0,0.0,0.0,0.0,0.7,0.7,0.7,0.0, ! 21 WHT WHIRCAM
     1 0.0,1.0,1.0,1.0,0.0,0.0,0.0,0.0,0.0,0.0, ! 22 Human eye
     1 0.0,0.0,0.0,0.0,0.0,0.7,0.7,0.7,0.7,0.0, ! 23 WHT CIRSI
     1 0.0,0.0,0.0,0.0,0.0,0.7,0.7,0.7,0.7,0.0, ! 24 INT CIRSI
C    1 0.0,0.0,0.0,0.0,0.0,1.6,1.6,1.6,1.6,0.0, ! 25 WHT INGRID to v12.7
     1 0.0,0.0,0.0,0.0,0.0,2.0,2.0,1.8,1.6,0.0, ! 25 WHT INGRID
     1 10*1.0,                                  ! 26 WHT NAOMI/INGRID
     1 10*1.0,                                  ! 27 WHT NAOMI/OASIS spec
     1 0.95,1.03,1.14,0.93,1.28,0.00,0.00,0.00,0.00,0.00, ! 28 NOT ALFOSC imag
     1 0.73,0.89,0.98,1.01,1.14,0.00,0.00,0.00,0.00,0.00, ! 29 NOT ALFOSC spec
     1 0.00,0.00,0.00,0.00,0.00,0.00,0.90,0.74,0.85,0.00, ! 30 NOT NOTCAM WF
     1 0.00,0.00,0.00,0.00,0.00,0.00,0.90,0.74,0.85,0.00, ! 31 NOT NOTCAM HR
     1 1.52,1.18,1.09,0.77,1.16,5*1.0,                    ! 32 NOT MOSCA imag
     1 5*0.0,0.29,0.53,3*0.8,                             ! 33 WHT LIRIS imag
     1 5*0.0,0.29,0.53,3*0.8,                             ! 34 WHT LIRIS spec
     1 6*1.0,0.00,0.00,0.00,0.00, ! 35 WHT ACAM imag  
     1 6*1.0,0.00,0.00,0.00,0.00, ! 36 WHT ACAM spec  
     1 10*1.0,                                        ! 37 WHT NAOMI/OASIS imag
     1 130*0.0/

C
C   Set default values of parameters
C
      data iinstr/1/
      data iband/3/
      data idet/15/
      data igrat/1/
      data qgrat/'     R158R'/
      data fwhm/1/
      data objm/20/
      data slit/1/
      data time/100/
      data ipx/0/
      data airmass/1.0/
C
C##############################################################################
C
C   Calculate PHOTONS from FLUXJY (v11.0 onward), set SKYMAG0 = SKYMAG
C
      do i=1,nband
         photons(i)=fluxjy(i)*10000/6.6/wave(i)
         skymag0(i)=skymag(i)
C         print *,'Band, photons ',i,photons(i)
      enddo
C
C   IP = 6 - interactive (normal mode)
C        7 - non-interactive, reads data from a file signal.dat
C            writes output to signal.out
C            (e.g. for web interface)
C        8 - non-interactive, reads data from a file standard input
C            writes output to standard output
C            (e.g. for web interface)
C   Non-interactive mode can also be triggered by typing X at the prompt.
C
C      print *,'IP (6, 7 or 8)? '
C      read(*,*) ip
       ip=6
C
C   If interactive, type welcome lines
C
      if(ip.eq.6) then
        iin = 1
        iout = 6
        call welcome(version)
      else if (ip.eq.7) then
        iin = 1
        iout = 7
        open(unit=iin,name='/tmp/signal.dat')
        open(unit=iout,name='/tmp/signal.out')
      else
        iin = 5
        iout = 6
      endif

C   
C   Calculate predicted instrument zeropoints, for default detector
C   and first grating in list for that instrument
C
      do j=1,ninstr

      do k=1,2 
C   K loop over default instrument/grating combinations
      if(k.eq.1) then
        iid=idefdet(j)
        iig=idefgrat(j)
      else if(k.eq.2) then
        iid=idefdet2(j)
        iig=idefgrat2(j)
      endif

      do i=1,nband
      zp(i,j,k)=photons(i)*10000*telarea(itel(j))*
     1  effmirror**float(imirror(j))*effinstr(i,j)*
     1  truethru(i,j)*effdet(i,iid)

      if(itype(j).eq.1) then
         zp(i,j,k)=zp(i,j,k)*band(i)
      else
         zp(i,j,k)=zp(i,j,k)*effgrat(i,iig,j)
      endif
      if(zp(i,j,k).ne.0) zp(i,j,k)=2.5*log10(zp(i,j,k))-extin(i)
      enddo
      enddo
      enddo

C========================================================================

C
C   Write out summary of instrument / detector parameters
C
C     
      iw=0
      if(iw.gt.0) call writeout(
     1 iw,iwb,iwd,iwi,
     1 ninstr,ndet,
     1 itel,imirror,itype,ngrating,iaxis,itrans,
     1 idefdet,idefdet2,idefgrat,idefgrat2,
     1 npixels,irot,
     1 telarea,arcmm,defslit,spread,dispmm,
     1 read,adu,pixel,
     1 effmirror,wave,band,fluxjy,photons,skymag,extin,
     1 effinstr,effgrat,effdet,truethru,zp,zmeas,
     1 zmeasc,zmeasd,
     1 version,cband,cinstr,comment,cgrating,cdet
     1 )


C===========================================================================

C
C   Main menu / calculation loop
C
      do loop=1,100000
C
C   Get parameters (from user or from SIGNAL.INPUT)
C
      call getpar(ip,iin,iout,ninstr,ngrating,ndet,nband,loop,
     1                  cinstr,cgrating,cdet,cband,cbandlower,
     1                         effgrat, effdet,
     1 wave,itel,idefdet,defslit,itrans,
     1                  iinstr,igrat,   idet,iband,
     1 objm,ipx,time,fwhm,slit,airmass,telarea,skymag,skymag0,
     1 extin,band,
     1 sol,eclat,glat,air,im,alpha,r,zm,ak,
     1 mode)

C
C   Program flow now determined by MODE
C
C      print *,'MODE = ',mode
      if(mode.eq.1) then
        call skybrcalc(ip,iout,skymag0,
     1  sol,eclat,glat,air,im,alpha,r,zm,ak)
        goto 500
      else if(mode.eq.2) then
        call help(ip,iout)
        goto 500
      else if(mode.eq.3) then
        call helpsky(ip,iout)
        goto 500
      else if(mode.eq.9) then
        goto 500
      endif

C
C   ELSE IF MODE=0, continue with S:N calculation
C
C
C
C   If IP=7, write html tags
C
      if(ip.eq.7) write(iout,492)
492   format('',/,'',/,
     1' SIGNAL ',/,'',/,
     1'',/,
     1'

Predicted throughput and signal-to-noise', 1'

',/,'

')

C
C   Convert input mag to required quantity
C
        if(objm.gt.100) then
          objm=objm-100
          objm=-1.e-6*10.**((23.9-objm)/2.5)
          write(iout,812) -objm
812       format('Brightness in Jy:',e10.3)
        endif
        if(objm.lt.0) then
          objm=-2.5*alog10(-objm/fluxjy(iband))
          write(iout,811) objm
811       format('Apparent mag: ',f7.2)      
        endif

C
C   Set IPX = 0 for IFU/fibre spectroscopy
C
      ipx0=ipx
      objm0=objm

      if(itype(iinstr).eq.3.and.ipx.eq.1) then
        ipx=0
        objm=objm-2.5*log10(slit*slit*3.14/4.)
      endif

C
C   Reset QGRAT
C
        if(igrat.ne.0) then
          qgrat=cgrating(igrat,iinstr)(1:10)
        else
          qgrat=' '
        endif


C
C   Calculate effective collecting area for TEL * INSTR * DET
C
      geomarea=telarea(itel(iinstr))
      area=geomarea

      thru(1)=10.**(-extin(iband)/2.5*airmass)
      thru(2)=effmirror**(float(imirror(iinstr)))
      thru(3)=effinstr(iband,iinstr)
      thru(4)=effgrat(iband,igrat,iinstr)
      area=area*thru(1)*thru(2)*thru(3)*thru(4)

C      area=geomarea
C    1 *10.**(-extin(iband)/2.5*airmass)
C    1 *effmirror**(float(imirror(iinstr)))
C    1 *effinstr(iband,iinstr)*effgrat(iband,igrat,iinstr)

     
c      print *,'area, geomarea = ',area,geomarea
c      print *,'extin, airmass, effmirror ',extin(iband),airmass,
c     1 effmirror
c      print *,'imirror,effinstr,grat ',imirror(iinstr),
c     1 effinstr(iband,iinstr),effgrat(iband,igrat,iinstr)

      throughput=area/geomarea
C      print *,'*** thruput ***', throughput

      area=area*effdet(iband,idet)      

      if(area.eq.0) then
        write(iout,813) 
813   format(
     1  '*** WARNING - throughput = 0  Grating OK?',
     1  '  Detector OK?  Band OK?',/)
      endif


C
C   For spectroscopic obs of point sources, calculate fraction
C   of light VIGN passing through slit, fibre or lenslet.
C   NB VIGN = 1 for combination ITYPE = 3 and user selects
C   IPX = 1 (i.e. IPX = 0, IPX0 = 1 at this point)
C
      vign=1
      if(itype(iinstr).ne.1.and.ipx.eq.0.and.ipx0.eq.0)
     1  call slitvign(iinstr,itype,slit,fwhm,vign)

C
C   Calculate scales at detector in arcsec/pixel and Angstrom/pixel
C
      arcpix=arcmm(iinstr)*pixel(idet)
      disppix=dispmm(igrat,iinstr)*pixel(idet)

C
C   Calculate number of photons from object and from sky/arcsec**2
C   per Angstrom
C
      sobj=photons(iband)*10000*area*time*10.**(-objm/2.5)*vign
      ssky=
     1     photons(iband)*10000*area*time*10.**(-skymag(iband)/2.5)

C
C   Correct for empirical/predicted throughput ratio
C
      sobj=sobj*truethru(iband,iinstr)
      ssky=ssky*truethru(iband,iinstr)


C
C   Display parameter settings to user
C
      qinst=cinstr(iinstr)
      qband=cband(iband)
      qdet=cdet(idet)(1:6)
      qgrat2=qgrat
      sky=skymag(iband)
      deteff=effdet(iband,idet)
      pixsize=pixel(idet)
      photons2=photons(iband)
      bwidth=band(iband)
      readout=read(idet)
      slit2=slit
      true0=truethru(iband,iinstr)
      ext=extin(iband)

      if(ngrating(iinstr).gt.0) then      
        qgrat2=qgrat
      else
        qgrat2=' '
      endif

      if(itype(iinstr).eq.1) then
        slit2=0
      else
        bwidth=0
      endif

      call display(ip,iout,ans,iinstr,itype,ipx,ipx0,ext,
     1  iband,cband,wave,
     1  qinst,qband,qgrat2,qdet,fwhm,objm,objm0,slit2,airmass,
     1  sky,time,geomarea,thru,throughput,deteff,pixsize,photons2,
     1  bwidth,readout,true0,comment)
C
C   Direct imaging calculation and output
C
      if(itype(iinstr).eq.1) then
        write(iout,'('' '')')
        if(ipx.eq.0) then
          sobj2=sobj*band(iband)
          peak=sobj2/1.14/fwhm/fwhm*arcpix*arcpix
        else
          sobj2=sobj*band(iband)*arcpix*arcpix
        endif
        ssky2=ssky*band(iband)*arcpix*arcpix
        call free(sobj2,buffa)
        call free(ssky2,buffb)
        call free(ssky,buffc)
        call free(peak,buffd)
        if(ipx.eq.0) then
          write(iout,821) buffa(1:10),buffd(1:10),buffc(1:10),
     1    buffb(1:10)
821       format('Detected photons            from object =',a10,
     1           ' (max',a10,' /pixel)',/,
     1           'Detected photons/A/arcsec^2 from sky    =',a10,/,
     1           'Detected photons/pixel      from sky    =',a10)
        else if(ipx.eq.1) then
          write(iout,822) buffa(1:10),buffc(1:10),buffb(1:10)
822       format('Detected photons/pixel      from object =',a10,/,
     1           'Detected photons/A/arcsec^2 from sky    =',a10,/,
     1           'Detected photons/pixel      from sky    =',a10)
        endif
      else

C
C   Spectroscopic calculation and output
C
        call free(sobj,buffa)
        if(ipx.eq.0) then
          write(iout,823) buffa(1:10),vign
823       format('Detected photons/A from object =',a10,
     1    '  (from in-slit fraction ',f4.2,')')
        else if(ipx.eq.1) then
          write(iout,824) buffa(1:10)
824       format(
     1 'Detected photons per A from each sq arcsec of object =',a10)
        endif
      
        if(ipx.eq.0) then
          sobj2=sobj*disppix
        else
          sobj2=sobj*slit*arcpix*disppix
        endif
        if(sobj2.gt.0) then
          call free(sobj2,buffa)
          if(ipx.eq.0) then
            write(iout,825) buffa(1:10)
825         format(
     1      'Detected photons/pixel-step-in-wavelength from object =',
     1      a10)
          else if(ipx.eq.1) then
            write(iout,826) buffa(1:10)
826         format('Detected photons/pixel from object  =',a10)
          endif
        endif
        ssky2=ssky*slit*arcpix*disppix
        if(itype(iinstr).eq.3) ssky2=ssky*slit*slit*3.14159/4.*disppix
        call free(ssky2,buffa)
        call free(ssky,buffc)
        if(ssky2.ne.0) then
          write(iout,734) buffc(1:10)
734       format('Detected photons/A/sq.arcsec from sky =',a10)
          if(itype(iinstr).eq.2) then
            write(iout,827) buffa(1:10)
827         format('Detected photons/pixel from sky =',a10)
          else if(itype(iinstr).eq.3) then
            write(iout,927) buffa(1:10)
927         format('Detected photons/pixel-step-in-wavelength',
     1             ' from sky =',a10)
          endif
        endif
      endif
C
C   Calculate signal-to-noise
C

C
C   Until version 11.0, used SEE = 3.14159/4.* for imaging
C
      if(ipx.eq.0.and.arcpix.ne.0) then
        if(itype(iinstr).eq.1) see=3.14159*(fwhm/arcpix)**2
        if(itype(iinstr).eq.2) see=2.*fwhm/arcpix
      else if(ipx.eq.1) then
        see=1.
      endif
      if(itype(iinstr).eq.3) see=spread(iinstr)
      if(itype(iinstr).le.2) then
        snr=sobj2/sqrt(sobj2+see*(ssky2+read(idet)*read(idet)))
      else if(itype(iinstr).eq.3) then
        snr=sobj2/sqrt(sobj2+ssky2+see*read(idet)*read(idet))
      endif
      if(snr.ne.0) then
      if(arcpix.ne.0.or.itype(iinstr).eq.3) then
        call free(snr,buffa)
        if(ipx.eq.0) then
          if(itype(iinstr).eq.1) then
              write(iout,828) buffa(1:10)
828           format('Signal-to-noise                         =',a10)
          else if(itype(iinstr).eq.2) then
            write(iout,829) buffa(1:10)
829         format(
     1      'Signal-to-noise/pixel of extracted spectrum =',a10)
          else if(itype(iinstr).eq.3) then
            write(iout,929) buffa(1:10)
929         format(
     1      'Signal-to-noise/pixel of extracted spectrum',
     1      ' (fibre or lenslet) =',a10)
          endif
        else if(ipx.eq.1) then
          write(iout,830) buffa(1:10)
830       format('Signal-to-noise/pixel                   = ',a10)
        endif
      endif
      endif

C
C   Scales at the detector
C   For spectroscopy:
C   IAXIS(IINSTR)  IROT(IDET)  dispersion direction
C      1             0            X
C      2             0            Y
C      1             1            Y
C      2             1            X
C
      write(iout,'('' '')')
      if(itype(iinstr).eq.1) then
        dx=npixels(idet,1)*arcpix
        dy=npixels(idet,2)*arcpix
        idx=1
        idy=1
      else
        if(iaxis(iinstr)+irot(idet).ne.2) then
          dx=npixels(idet,1)*disppix
          dy=npixels(idet,2)*arcpix
          idx=2
          idy=1
        else
          dx=npixels(idet,1)*arcpix
          dy=npixels(idet,2)*disppix
          idx=1
          idy=2
        endif
      endif 
C     print *,'*** ',iinstr,iaxis(iinstr),idet,irot(idet),idx,idy
      if(dx.ne.0) write(iout,259) 
     1 dx/float(npixels(idet,1)),unit(idx),
     1 npixels(idet,1),dx,unit(idx)
259   format('Scale in X direction: ',f7.2,1x,a1,'/pixel; ',
     1 i4,' pixels = ',f9.2,1x,a1)
      if(dy.ne.0) write(iout,260) 
     1 dy/float(npixels(idet,2)),unit(idy),
     1 npixels(idet,2),dy,unit(idy)
260   format('Scale in Y direction: ',f7.2,1x,a1,'/pixel; ',
     1 i4,' pixels = ',f9.2,1x,a1)

C
C   Photons/ADU conversion
C
      write(iout,151) cdet(idet)(1:6),adu(idet)
C   ,int(read(idet))
151   format(a6,' gain =',f4.1,
     1       ' e-/ADU approx.')
C   ; readnoise =',i3,' electrons')

C
C   If IP=7, write closing html tags
C
      if(ip.eq.7) write(iout,493)
493   format('
',/,'',/,'') C C Reset IPX to IPX0 and OBJM to OBJM (changed above if this loop C was for IPX = 1, ITYPE = 3) C ipx=ipx0 objm=objm0 enddo 500 continue end C============================================================================== C C S U B R O U T I N E S C C============================================================================== C****************************************************************************** C***** D I S P L A Y ****************************************************** C****************************************************************************** C subroutine display(ip,iout,ans,iinstr,itype, 1 ipx,ipx0,ext, 1 iband,cband,wave, 1 qinst,qband,qgrat,qdet,fwhm,objm,objm0,slit,airmass, 1 sky,time,geomarea,thru,throughput,deteff,pixsize,photons, 1 bwidth,readout,true0,comment) C C Displays parameter settings to user C C NB as of version 14.1, writes IPX0 rather than IPX, and C OBJM0 rather than OBJM, C the combination ITYPE = 3, IPX = 0 C implicit none integer i,iinstr,itype(50),ipx,ipx0,iband,ip,iout real fwhm,objm,sky,time,geomarea,throughput,deteff,pixsize, 1 photons,bwidth,readout,slit,true0,ext,airmass,wave(10) real objm0,thru(4) character*6 qdet character*10 qgrat,qgratbuff character*37 text(20) character*25 qinst character*1 ans character*2 qband,cband(10) character*80 cthru character*80 comment(50) do i=1,20 text(i)=' ' enddo qgratbuff=qgrat if(qgratbuff.eq.' ') qgratbuff=' none ' print *,' ' text(1)='G Grating '//qgratbuff text(2)='D Detector '//qdet C C Text strings which depend on whether imaging or spectroscopic C if(itype(iinstr).eq.1) then C C Imaging C text(3)='B Band '//qband write(text(4),'( 1 ''M Apparent magnitude '',f6.1)') objm0 write(text(15),'( 1 ''W Effective bandwidth (A) '',f9.2)') bwidth text(8)=' ' write(text(12),'( 1 '' Atm * tel * instr throughput '',f6.2)') 1 throughput cthru=' ' write(cthru,200) (thru(i),i=1,3) 200 format(' =',f6.2,' *',f6.2,' *',f6.2) else C C Spectroscopic C if(itype(iinstr).eq.2) then write(text(8),'( 1 ''S Slit width (arcsec) '',f6.1)') slit else if(itype(iinstr).eq.3) then write(text(8),'( 1 ''S Fibre/lenslet diam (arcsec) '',f6.1)') slit endif if(itype(iinstr).eq.2.or.itype(iinstr).eq.3) then write(text(12),'( 1 '' Atm*tel*inst*grat throughput '',f6.2)') 1 throughput cthru=' ' write(cthru,201) (thru(i),i=1,4) 201 format(' =',f6.2,' *',f6.2,' *',f6.2,' *',f6.2) endif if(ipx0.eq.0) then write(text(3),'( 1 ''B Band '',a2, 1 '' = wavelength (A) '',i5)') 1 qband,int(wave(iband)) write(text(4),'( 1 ''M Apparent magnitude '',f6.1)') objm0 else write(text(4),'( 1 ''M Apparent magnitude/sq.arcsec'',f6.1)') objm0 endif text(15)=' ' endif C C Remaining text strings C write(text(5),'( 1 ''P Mag (0) or mag/sq.arcsec (1)'',i4)') ipx0 write(text(7),'( 1 ''F FWHM (object*seeing, arcsec)'',f6.1)') fwhm write(text(6),'( 1 ''T Integration time (sec) '',f9.1)') time write(text(9),'( 1 ''A Airmass '',f6.1)') airmass write(text(11),'( 1 ''U Tel. geometric area (sq.m) '',f6.2)') geomarea write(text(13),'( 1 '' Detector efficiency '',f6.2)') deteff write(text(19),'( 1 '' Pixel size (microns) '',f6.1)') pixsize*1000 write(text(18),'( 1 '' Readout noise (electrons) '',f6.1)') readout write(text(16),'( 1 '' Photons/s/A/sq.cm mag=0 '',i5)') int(photons) write(text(17),'( 1 ''K Sky brightness mag/sq.arcsec '',f6.2)') sky write(text(10),'( 1 ''E Extinction per airmass (mag) '',f6.2)') ext write(text(14),'( 1 '' Empirical/theoretical '',f6.2)') true0 if(ip.eq.6) then write(iout,101) qinst,text(10), 1 text(1),text(11), 1 text(2),text(12),cthru(1:32), 1(text(i),text(i+10),i=3,9), 1 comment(iinstr) 101 format(' I Instrument = ',a25,1x,a37,/, 1 1x,a37,4x,a37,/, 1 1x,a37,4x,a37,a32,/, 1 6(1x,a37,4x,a37,/), 1 1x,a37,4x,a37,/,a80) else write(iout,102) qinst,text(10)(3:37), 1 text(1)(3:37),text(11)(3:37), 1 text(2)(3:37),text(12)(3:37),cthru(1:32), 1 (text(i)(3:37),text(i+10)(3:37),i=3,9), 1 comment(iinstr) 102 format('Instrument = ',a25,1x,a35,/, 1 a35,4x,a35,/, 1 a35,4x,a35,a32,/, 1 7(a35,4x,a35,/),a80) endif end C****************************************************************************** C***** F R E E ************************************************************ C****************************************************************************** C subroutine free(x,buff15) C C Writes real X into buffer BUFF15 free-format, as in PRINT *,X C (unix fortran print * not convenient because renders numbers C outside 0.1049 - 1012863 (why?) in exponential format, and may C not yield same formats as other fortran) C implicit none real x character*15 buff15 C print *,'okfree 1',x,buff15 if(x.lt.1000000.and.x.gt.0.1) then write(buff15,'(f10.2,5x)') x else write(buff15,'(e10.3,5x)') x endif C print *,'okfree 2',x,buff15 end C****************************************************************************** C***** G E T P A R ******************************************************** C****************************************************************************** C subroutine getpar(ip,iin,iout,ninstr,ngrating,ndet,nband,loop, 1 cinstr,cgrating,cdet,cband,cbandlower, 1 effgrat, effdet, 1 wave,itel,idefdet,defslit,itrans, 1 iinstr,igrat, idet,iband, 1 objm,ipx,time,fwhm,slit,airmass,telarea,skymag,skymag0, 1 extin,band, 1 sol,eclat,glat,air,im,alpha,r,zm,ak, 1 mode) C C CMDOK string of acceptable user command characters C implicit none integer ninstr,ngrating(50),ndet,nband,itel(50),loop,iunit integer iinstr,igrat,idet,iband,ipx,mode,i,j,k,idefdet(50) integer itrans(50) real defslit(50) real effgrat(10,20,50),effdet(10,30),wave(10),objm,time real fwhm,slit,airmass,telarea(50),skymag(10),extin(10), 1 band(10),dummy character*1 ans character*2 cband(10),cbandlower(10),ans2 character*20 cgrating(20,50) character*25 cinstr(50) character*30 cinstrb(40) character*40 cdet(30) character*40 cmdok character*80 cline C C Variables used only in call to skybrcalc C integer ip,iin,iout,im real skymag0(10),sol,eclat,glat,air,alpha,r,zm,ak data cmdok/'IGDBMPTFSAKEWigdbmptfsakewQHqhUuXx '/ if(ip.ne.6) then C C ******************************** C C Get parameters from SIGNAL.INPUT C C ******************************** C mode=9 read(iin,*,end=500) mode if(mode.eq.0) then read(iin,*) iinstr read(iin,*) igrat if(igrat.eq.0) igrat=1 read(iin,*) idet read(iin,*) iband read(iin,*) objm read(iin,*) ipx read(iin,*) time read(iin,*) fwhm read(iin,*) slit read(iin,*) airmass read(iin,'(a80)') cline if(cline.eq.'D'.or.cline.eq.'d') then if(iband.le.6) skymag(iband)=skymag0(iband)-0.4 else if(cline.eq.'G'.or.cline.eq.'g') then if(iband.le.3) skymag(iband)=skymag0(iband)-2.15 if(iband.ge.4.and.iband.le.6) 1 skymag(iband)=skymag0(iband)-1.3 else if(cline.eq.'B'.or.cline.eq.'b') then if(iband.le.3) skymag(iband)=skymag0(iband)-3.4 if(iband.ge.4.and.iband.le.6) 1 skymag(iband)=skymag0(iband)-2.7 else read(cline,*) dummy if(dummy.ne.0) skymag(iband)=dummy endif read(iin,*) dummy if(dummy.ne.0) extin(iband)=dummy read(iin,*) dummy if(dummy.ne.0) band(iband)=dummy else if(mode.eq.1) then im=1 read(iin,*) sol read(iin,*) eclat read(iin,*) glat read(iin,*) air read(iin,*) alpha if(alpha.eq.-1) im=0 read(iin,*) zm read(iin,*) r read(iin,*) ak endif else C C ************************ C C Get parameters from user C C ************************ C iunit=5 print *,' ' 173 continue if(loop.gt.1) print *, 1 'Type I G D B M P T F S A U W K or E to change a parameter' print *, 1 'Type return to continue,', 1 ' Q to quit, H for more information' C print *,' ' read(iunit,'(a1)') ans C C Check command is one of set C do k=1,40 if(ans.eq.cmdok(k:k)) goto 174 enddo print *,'Invalid command: ',ans goto 173 174 continue C C Stop C if(ans.eq.'Q'.or.ans.eq.'q') mode=9 C C Type rudimentary help C if(ans.eq.'H'.or.ans.eq.'h') then call help(ip,iout) goto 27 endif C C Get band C if(ans.eq.'B'.or.ans.eq.'b') then write(*,193) 193 format('Approximate throughput for which of', 1 ' following bands?',/) do j=1,nband write(*,196) cband(j),int(wave(j)) 196 format(a2,' =',i6,' A for spectroscopy') enddo read(iunit,'(a2)') ans2 do j=1,nband if(ans2.eq.cband(j).or.ans2.eq.cbandlower(j)) iband=j enddo endif C C Get bandwidth C if(ans.eq.'W'.or.ans.eq.'w') then print *,'Current effective bandwidth for ',cband(iband) 1 ,' = ',band(iband) print *,'(= area under filter transmission curve)' print *,'New effective bandwidth for this band?' read(iunit,*) band(iband) endif C C Get apparent magnitude C if(ans.eq.'M'.or.ans.eq.'m') then write(6,104) cband(iband) 104 format(' Apparent magnitude in band ',a2,'?') print *, 1 '(number > 50 will be interpreted as 100 + Oke AB mag)' print *, 1 '(negative number will be interpreted as Jy', 1 ' (1E-26 W/Hz/m**2))' read(iunit,*) objm endif C C Get telescope/instrument combination C if(ans.eq.'I'.or.ans.eq.'i') then print *,'Telescope/instrument? (N/A = no longer offered)' do k=1,40 cinstrb(k)=' ' enddo do k=1,ninstr i=itrans(k) write(cinstrb(k),101) k,cinstr(i) 101 format(i3,2x,a25) enddo do k=1,20 write(6,'(a30,5x,a30)') cinstrb(k),cinstrb(k+20) enddo print *,' ' write(6,102) 102 format(' Choose a number ', 1 '(NB detector reverts to a default', 1 ' type on changing instrument)') read(iunit,*) k iinstr=itrans(k) idet=idefdet(iinstr) slit=defslit(iinstr) endif C C Get slit width C if(ans.eq.'S'.or.ans.eq.'s') then print *,'Slit width (or sqrt fibre/lenslet area) (arcsec)?' read(iunit,*) slit endif C C Get IPX (point or extended) C if(ans.eq.'P'.or.ans.eq.'p') then print *,'Mag is for point source (0) or', 1 ' per arcsec**2 for an extended source (1)?' read(iunit,*) ipx endif C C Get extinction C if(ans.eq.'E'.or.ans.eq.'e') then print *, 1 'New value of extinction per airmass for this band?' read(iunit,*) extin(iband) endif C C Get new value for telescope geometric area C if(ans.eq.'U'.or.ans.eq.'u') then print *,'Telescope geometric area = ', 1 telarea(itel(iinstr)) print *,'New value?' read(*,*) telarea(itel(iinstr)) endif C C Get airmass C if(ans.eq.'A'.or.ans.eq.'a') then print *,'Airmass?' read(iunit,*) airmass endif C C Get sky brightness C if(ans.eq.'K'.or.ans.eq.'k') then c729 continue print *,'Current value of sky brightness in ', 1 cband(iband),' = ',skymag(iband), 1 ' mags/arcsec**2' print *,'Enter new value for sky brightness' print *,'or enter 0 (zero) to run the sky-brightness', 1 ' calculator' print *,' (NB the calculator does not', 1 ' reset the current value)' print *,'or enter D, G or B to select', 1 ' typical sky brightness in', 1 ' dark, grey or bright of moon' print *,' (NB typical dark is not same as minimum)' print *,'or enter H for help' read(iunit,'(a80)') cline if(cline.eq.'H'.or.cline.eq.'h') then call helpsky(ip,iout) c goto 729 else if(cline.eq.'D'.or.cline.eq.'d') then if(iband.le.6) skymag(iband)=skymag0(iband)-0.4 else if(cline.eq.'G'.or.cline.eq.'g') then if(iband.le.3) skymag(iband)=skymag0(iband)-2.15 if(iband.ge.4.and.iband.le.6) 1 skymag(iband)=skymag0(iband)-1.3 else if(cline.eq.'B'.or.cline.eq.'b') then if(iband.le.3) skymag(iband)=skymag0(iband)-3.4 if(iband.ge.4.and.iband.le.6) 1 skymag(iband)=skymag0(iband)-2.7 else read(cline,*) dummy if(dummy.ne.0) then skymag(iband)=dummy else if(dummy.eq.0) then call skybrcalc(ip,iout,skymag0, 1 sol,eclat,glat,air,im,alpha,r,zm,ak) c goto 729 endif endif endif C C If NGRATING(IINSTR) > 0, get grating C if(ans.eq.'G'.or.ans.eq.'g') then print *, 1 'Grating? (NB for WYFFOS, the list below is', 1 ' 9 ISIS then 10 IDS gratings)' if(ngrating(iinstr).gt.0) then do i=1,ngrating(iinstr) write(6,105) i,cgrating(i,iinstr), 1 effgrat(iband,i,iinstr) 105 format(i2,2x,a20,' throughput in this band = ',f4.2) enddo write(6,172) 172 format('Choose a number (leftmost column above)') read(iunit,*) igrat endif endif C C Get detector C if(ans.eq.'D'.or.ans.eq.'d') then print *,'Detector?' print *, 1'(See ING web page for detailed detector info)' do i=1,ndet write(6,106) i,cdet(i),effdet(iband,i) 106 format(i3,2x,a40,' QE in this band = ',f4.2) enddo write(6,172) read(iunit,*) idet endif C C Get exposure time C if(ans.eq.'T'.or.ans.eq.'t') then print *,'Exposure time (sec)?' read(iunit,*) time endif C C Get FWHM C if(ans.eq.'F'.or.ans.eq.'f') then print *,'FWHM object * seeing?' read(iunit,*) fwhm endif 27 continue endif 500 continue end C****************************************************************************** C***** H E L P ************************************************************ C****************************************************************************** C subroutine help(ip,iout) C C C implicit none integer ip,iout write(iout,'('' '')') write(iout,101) 101 format( 1'Full details of the program''s operation are given in the',/, 1'comments preceding the code.',/, 1'More up-to-date help is given on the SIGNAL web page:',/, 1' http://www.ing.iac.es/~crb/signal.html',/, 1' ',/, 1'For imaging, the bandwidths used correspond to the',/, 1'50-mm or liq-CuSO4 U filters, the Harris B,V,R,I set,',/, 1'the RGO glass Z filter and the WHIRCAM J, H and Kshort',/, 1'filters. Note that the 125-mm glass U filter',/, 1'has a factor 3 lower throughput than the 50-mm U.',/, 1'Note also that the throughputs of the Harris B filters are',/, 1'up to 25% less than that of the KPNO B filters.',/, 1'To specify narrow-band filters or adjust the effective',/, 1'bandwidth, use command B.',/, 1' ',/, 1'Object brightnesses can also be given in Jy or Oke AB mag',/, 1'(use the M command).') if(ip.eq.6) call pause write(iout,102) 102 format( 1' ',/, 1'The program does not take into account:',/, 1'- Vignetting at large field radius',/, 1'- Vignetting by the spectrograph slit (the specified width',/, 1' of the slit is used only in calculating the contribution',/, 1' from the sky).',/, 1'- Loss of light to colour and ND filters or to polarisation',/, 1' optics in spectrograph',/, 1'- Loss of grating efficiency at large angles',/, 1' and loss of light to neighbouring orders',/, 1'- Brightening of the sky due to twilight, moonlight,',/, 1' solar activity, low ecliptic or galactic latitude,',/, 1' or light pollution (enhanced in dusty conditions)',/, 1' - see option K in the menu',/, 1'- For FOS, U,B are for second order, V,R,I are first order.') if(ip.eq.6) call pause write(iout,103) 103 format( 1' ',/, 1'The signal-to-noise calculation:',/, 1' Nobj = photons/A (per arcsec**2 if extended) from object',/, 1' Nsky = photons/A/arcsec**2 from sky',/, 1' BAND = equivalent width of filter in A (integral T(l)dl',/, 1' where T(l) is transmission, l is wavelength)',/, 1' P = number of pixels over which integration carried out',/, 1' READ = CCD readout noise (e-)',/, 1' FWHM = object fwhm (intrinsic and due to seeing) in arcsec',/, 1' ',/, 1' Imaging, point or extended sources:',/, 1' N = Nobj * BAND for point sources',/, 1' N = Nobj * BAND*(arcsec/pixel)**2 for extended sources',/, 1' S = Nsky * BAND*(arcsec/pixel)**2',/, 1' P = pi*(FWHM/(arcsec/pixel))**2 for point source',/, 1' P = 1 for an extended source') if(ip.eq.6) call pause write(iout,104) 104 format( 1' ',/, 1' Spectroscopy, point source:',/, 1' N = Nobj * A/pixel',/, 1' S = Nsky * slit-width(A) * arcsec/pixel * A/pixel',/, 1' P = 2 * FWHM (arcsec) / arcsec/pixel',/, 1' ',/, 1' Spectroscopy, extended source:',/, 1' N = Nobj * slit-width(A) * arcsec/pixel * A/pixel',/, 1' S = Nsky * slit-width(A) * arcsec/pixel * A/pixel',/, 1' P = 1',/, 1' ',/, 1' Signal-to-noise = N/sqrt(N+P*(S+READ**2))',/, 1' For point sources, the sky counts are those from a',/, 1' 2-FWHM-diameter circular aperture for imaging, and from',/, 1' a 2-FWHM * slit-width rectangle for spectroscopy.') write(iout,105) 105 format( 1' For extended sources, the calculations are per pixel.') end C****************************************************************************** C***** H E L P S K Y ****************************************************** C****************************************************************************** C subroutine helpsky(ip,iout) implicit none integer ip,iout character*1 pans character*70 prompt pans=' ' write(iout,101) 101 format( 1'SIGNAL''s default sky-brightness settings are median for',/, 1'dark-of-moon, solar minimum, at high galactic and ecliptic',/, 1'latitude and in the absence of twilight and moonlight.',/, 1'The BVR values are accurate to +- 0.1 mag and are similar to',/, 1'those measured at other dark sites (Chile, Hawaii etc.).',/, 1'U and I are accurate to 0.5 mag. R and I sky brightnesses',/, 1'vary randomly by several tenths of a mag with variations in',/, 1'the OH airglow. The V and R sky brightness include a',/, 1'contribution of about 0.1 mag due to NaD light pollution.',/, 1'Light pollution is negligible in other bands, and it is',/, 1'decreasing.',/, 1'The sky is markedly brighter (several tenths of a mag) under',/, 1'very dusty conditions (> 0.3 mag extinction).',/, 1'The sky is 0.4 mag brighter at solar maximum. Recent minima') prompt= 1'[Return to continue, Q to end help]' if(ip.eq.6) call pause2(pans,prompt) if(pans.eq.'Q'.or.pans.eq.'q') goto 22 write(iout,102) 102 format( 1'were in 1986.5, 1996.5. Variation is approximately',/, 1'sinusoidal with time.',/, 1'The sky is 0.4 mag brighter on the ecliptic than at the',/, 1'poles, varying as sine(b) approx.',/, 1'The airglow contribution (typically about 70% of the total',/, 1'in V) brightens approx as airmass. The sky is 0.3 mag',/, 1'brighter at airmass 1.5.',/, 1'Stars fainter than mag=20 contribute', 1' negligibly to the total.',/, 1'Starlight scattered by interstellar dust contributes about',/, 1'5% of the total, rising to about 30% on the galactic plane.',/, 1'The extragalactic contribution is negligible (< 1%).',/, 1'The brightness of the sky does not vary with time after',/, 1'astronomical twilight.',/, 1'For further details of the calculations of moonless sky',/, 1'brightness, see La Palma technical note 115',/, 1'(Benn & Ellison 1998).') if(ip.eq.6) call pause2(pans,prompt) if(pans.eq.'Q'.or.pans.eq.'q') goto 22 write(iout,103) 103 format( 1' ',/, 1'SKY BRIGHTNESS WITH MOON UP:',/, 1'When the moon is 60 deg from zenith,', 1' with extinction 0.15 mag,',/, 1'the zenith sky will brighten roughly by M as tabulated below:',/, 1' ',/, 1' New Crescent Quarter Gibbous Full',/, 1' Phase angle (deg) 180 135 90 45 0 ',/, 1' Approx day: 1 4 8 12 15',/, 1' D, G or B: D G G B B',/, 1' Illum. frac. % 0 25 50 75 100',/, 1' M (U, B, V) 0 0.5 2.0 3.1 4.3',/, 1' M (R) 0 0.3 1.3 2.4 3.5',/, 1' M (I) 0 0.2 1.1 2.2 3.3',/, 1' ') if(ip.eq.6) call pause2(pans,prompt) if(pans.eq.'Q'.or.pans.eq.'q') goto 22 write(iout,104) 104 format( 1'Sky brightness for other values of lunar phase, lunar zenith',/, 1'angle, sky position and extinction, can be estimated with',/, 1'SIGNAL''s sky-brightness calculator.',/, 1'Note that the quarter moon (i.e. half disc illuminated)',/, 1'is a factor of 10 (not 2) fainter than full, due to the',/, 1'opposition effect (also responsible for gegenschein on the',/, 1'ecliptic and dry heiligenschein on earth).',/, 1'The contribution of moonlight in V has been calculated',/, 1'according to the scattering formula of Krisciunas & Schaefer',/, 1'(1991, PASP, 103, 1033), normalised (multiplied by a factor',/, 1'of 2.4) to agree with measurements of sky brightness made at',/, 1'the JKT on a dust-free night in 7/98.',/, 1'The moonlight contribution in the other bands is calculated',/, 1'according to the U-B, B-V, V-R, R-I colours of moonlight',/, 1'measured on the same night in 7/98.') if(ip.eq.6) call pause2(pans,prompt) if(pans.eq.'Q'.or.pans.eq.'q') goto 22 write(iout,105) 105 format( 1'These values agree +-40% with measurements made by DHPJ in',/, 1'9/89, but the contribution of moonlight probably depends',/, 1'strongly on local conditions (e.g. dust, telescope baffling),',/, 1'and with current data, the contribution by moonlight on',/, 1'La Palma can probably only be predicted within a factor ~2.',/, 1' ',/, 1'Ian Steele (LJMU) has found that the background brightness',/, 1'at the JKT rises dramatically (factor >~5 brighter than',/, 1'the above numbers) if moonlight falls on the telescope',/, 1'structure (scattering within the telescope).') 22 continue end C****************************************************************************** C***** P A U S E ********************************************************** C****************************************************************************** C subroutine pause implicit none character*1 ans print *,'[Press return to continue]' read(*,'(a1)') ans end C****************************************************************************** C***** P A U S E 2 ******************************************************** C****************************************************************************** C subroutine pause2(pans,prompt) implicit none character*1 pans character*70 prompt print *,prompt read(*,'(a1)') pans end C****************************************************************************** C***** S K Y B R C A L C ************************************************** C****************************************************************************** C subroutine skybrcalc(ip,iout,skymag, 1sol,eclat,glat,air,im,alpha,r,zm,ak) C C SKYBRCALC calculates the brightness of the sky in each band, C in S10 units and in mag/arcsec^2 as a function of C ecl, gal latitudes, zd, phase of solar cycle, and, if moon up: C phase of moon, zenith distance of moon, separation between C sky position and moon, extinction C If IP=6, the parameters are requested from user, otherwise C those supplied in call are used. C C Variables C --------- C AIR Airmass (derived from ZD) C AK V-band extinction in mag C ALPHA Phase of moon (0 = full, 90 = 7-day, 180 = new) C ECLAT Ec lat (deg) C GLAT Gal lat (deg) C IM = 1 if moonlight to be considered C IP Unit number to which to send output (6 = screen) C NBAND Number of wavelength bands C QAIR Airglow contribution in S10 C QSTAR Diffuse starlight contribution in S10 C QZOD Zodiacal light contibution in S10 C Q3 QAIR+QSTAT+QZOD in S10 C QSKY(10) Moonless sky brightness in S10 derived from SKY C QMOON(10) Contrib of moon in S10 C QALL(10) QSKY+QMOON C R Ang sep of patch of sky and moon (deg) C SKY Moonless sky br in mag/arcsec^2 derived from Q3, C correcting for colour difference (band - V) C SKYALL Total sky br in mag/arcsec^2 C SOL Phase of solar sycle (0 = min, 1 = max) converted t C 0.8 + 1.2*SOL C ZM Zenith distance of moon (deg) C implicit none integer i,im,nband,j,ip,iout,loop real skymag(10),sol,eclat,glat,air,alpha,r,zm,ak real qair,qzod,qstar,q3 real qsky(10),sky(10),qmoon(10),qall(10),skyall(10) nband=5 loop=1 if(ip.eq.6) loop=10000 do i=1,loop if(ip.eq.6) then print *,'--------------------------------------------------' print *,'Running sky-brightness calculator' print *, 1 'Solar activity (0 = minimum, 1 = maximum, 9 to finish', 1 ' with calculator)?' print *, 1 '(Dates of last three solar minima: 1976.8, 1986.8, 1996.5)' read(*,*) sol if(sol.eq.9) goto 22 print *, 1 'Ecliptic latitude in deg (longitude', 1 ' assumed to be > 90 deg from sun)?' read(*,*) eclat print *,'Galactic latitude in deg?' read(*,*) glat print *,'Zenith distance (deg)?' read(*,*) air air=1./sqrt(1-0.96*(sin(air/57.3))**2.) print *,'Phase angle of moon (-1 for moonless sky)?' print *,'(e.g. 0 = full, 90 = 7-day old moon, 180 = new moon)' read(*,*) alpha im=0 if(alpha.ne.-1) then im=1 print *,'Zenith distance of moon (deg)?' read(*,*) zm print *, 1 'Angular separation of moon and target patch of sky (deg)?' read(*,*) r print *,'V-band extinction (mag)?' read(*,*) ak endif endif ! if IP=6 C C Convert SOL to MJy C sol=0.8+1.2*sol call skybrcalc2(im,nband,skymag,sol,eclat,glat,air, 1 alpha,r,zm,ak, 1 qair,qzod,qstar,q3, 1 qsky,sky,qmoon,qall,skyall) C C If IP=7|8, write html header C if(ip.ne.6) write(iout,494) 494 format('',/,'',/, 1' SIGNAL ',/,'',/, 1'',/, 1'

Predicted brightness of La Palma sky', 1'

',/,'

')

      write(iout,116) (sol-0.8)/1.2,eclat,glat,acos(1/air)*57.3,air
116   format(
     1 'Solar activity =',f4.1,', ecliptic and galactic',
     1 ' latitudes =',f5.1,',',f5.1,' deg,',/,
     1 ' ZD = ',f4.1,' deg, airmass = ',f5.3)
      if(im.eq.1) write(iout,117) alpha,zm,r,ak
117   format('Lunar phase = ',f5.1,' deg, lunar zenith distance =',
     1 f5.1, ' deg,',/,
     1 '  distance from moon = ',f5.1,' deg, extinction =',
     1 f5.2,' mag')
      write(iout,101)
     1 q3,qair,qzod,qstar,
     1 skymag(1)-skymag(2),skymag(2)-skymag(3),skymag(3)-skymag(4),
     1 skymag(4)-skymag(5),
     1 (qsky(j),j=1,nband),(sky(j),j=1,nband)
101   format(
     1 'Moonless V sky brightness, in S10 units ',
     1 '(equivalent 10th-mag stars/deg^2):',/,
     1 f6.1,' ( =',f6.1,' airglow +',f6.1,' zodiacal light +',f6.1,
     1 ' starlight)',/,
     1 'The (moonless) sky brightness in other bands is obtained by',
     1 ' assuming ',/,
     1 'colours U-B, B-V, V-R, R-I = ',4f6.2,':',//,
     1 '  Band                 ',
     1 'U        B        V        R        I',/,
     1 '  Total (no moon)',5f9.1,'  S10 units',/,
     1 '  Total (no moon) ',5f9.2,' mag/arcsec^2')
      if(im.eq.1) 
     1 write(iout,102) 
     1 (qmoon(j),j=1,nband),(qall(j),j=1,nband),
     1 (skyall(j),j=1,nband),(sky(j)-skyall(j),j=1,nband)
102   format(/,
     1 '  Moonlight      ',5f9.1,'  S10 units',/,
     1 '  Total          ',5f9.1,'  S10 units',/,
     1 '  Total           ',5f9.2,' mag/arcsec^2',/,
     1 ' (Change with moon',5f9.2,' mag)',/)
      write(iout,109)
 109  format('Accuracy:',
     1 ' +- 0.1 mag without moon (in B,V),',
     1 ' +- 1 mag with moon',/)
C
C   If IP=7|8, write closing html tags
C
      if(ip.ne.6) write(iout,495)
495   format('
',/,'',/,'') enddo 22 continue end C****************************************************************************** C***** S K Y B R C A L C 2 ************************************************ C****************************************************************************** C subroutine skybrcalc2(im,nband,skymag,sol,eclat,glat,air, 1 alpha,r,zm,ak, 1 qair,qzod,qstar,q3, 1 qsky,sky,qmoon,qall,skyall) C C Calculates sky brightness for given conditions, using LP tech note C 115 (Benn & Ellison) for dark-sky brightness, and Krisciunas & Schaefer C 1991 PASP 103 1033 p.1038 for moonlight contribution. C All working in S10 units until final conversion to mag/arcsec**2 C C Variables not used in SKYBRCALC: C S, FR, XZ, XZM, BNL - terms at start Sect 5 of Krisciunas C & Schaefer for calculating moonlight contr. C BNL is moonlight sky brightness in V in nanoLamberts C FUDGE(I) normalise the Krisciunas prediction to values measured C in 6 bands at the JKT 5/7/98 C BS10 = BNL * 3.8, the V moonlight sky br in S10 units C implicit none integer i,im,nband real sol,eclat,glat,air,qair,qzod,qstar,q3 real qsky(10),sky(10),qmoon(10),qall(10),skyall(10) real skymag(10),fudge(10) real alpha,s,fr,r,xz,zm,xzm,bnl,ak,bs10 data fudge/2.0,1.6,2.4,3.0,5.9,5*0.0/ C C Calculate V sky brightness without moon C qair=(145+130*(sol-0.8)/1.2)*air if(eclat.lt.60) then qzod=140-90.*sin(eclat/57.3) else qzod=60 endif qstar=100.*exp(-abs(glat)/10.) q3=qair+qzod+qstar C C Assume sky colour implied by values in SKYMAG, and calculate C dark-sky brightness in the other bands C do i=1,nband sky(i)=27.78-2.5*log10(q3)+skymag(i)-skymag(3) qsky(i)=10.**((27.78-sky(i))/2.5) enddo C C Calculate contribution of moon in each band C if(im.eq.1) then s=10.**(-0.4*(3.84+0.026*alpha+4.e-9*alpha**4.)) fr=10.**5.36*(1.06+(cos(r/57.3))**2)+10.**(6.15-r/40.) xz =air xzm=1./sqrt(1-0.96*(sin(zm/57.3))**2.) bnl=s*fr*10.**(-0.4*ak*xzm)*(1-10.**(-0.4*ak*xz)) bs10=bnl*3.8 C print *,'s,fr,xz,xzm,bnl,bs10 ',s,fr,xz,xzm,bnl,bs10 do i=1,nband qmoon(i)=fudge(i)*bs10 qall(i)=qsky(i)+qmoon(i) skyall(i)=27.78-2.5*log10(qall(i)) enddo endif end C****************************************************************************** C***** S L I T V I G N **************************************************** C****************************************************************************** C subroutine slitvign(iinstr,itype,slit,fwhm,vign) C C Calculates VIGN, the fraction of light entering slit for given C pointlike object FWHM. C implicit none integer iinstr,itype(50) real slit,fwhm,vign,sigma,r,rat vign=1 C C Slit spectroscopy C Formulae below are approx (better than 5%) C fit to dependence of light fraction in slit vs slit/fwhm C given by light_in_slit.f numerical simulations C if(itype(iinstr).eq.2) then if(fwhm.eq.0) then vign=1.0 else rat=slit/fwhm if(rat.lt.0.76) vign=0.868*rat if(rat.ge.0.76.and.rat.lt.1.40) vign=0.37+0.393*rat if(rat.ge.1.40.and.rat.lt.2.30) vign=1.00-0.089*(2.3-rat) if(rat.ge.2.3) vign=1.0 endif endif C C Fibre/lenslet spectroscopy C r=slit/2. sigma=fwhm/2.35 if(itype(iinstr).eq.3) then vign=1.0-exp(-r*r/2./sigma/sigma) endif end C****************************************************************************** C***** W E L C O M E ****************************************************** C****************************************************************************** C subroutine welcome(version) implicit none character*28 version print *,' ' print *,'SIGNAL '//version 1 //' crb@ing.iac.es' print *,' ' write(*,101) 101 format( 1'SIGNAL calculates photons/A detected from a source of given',/, 1'apparent magnitude (per arcsec**2 if extended), as:',/, 1' 10.**(-mag/2.5)',/, 1' * photons/sec/A/m**2 giving mag = 0 at top of atm.',/, 1' * transmission of atm. at given airmass',/, 1' * exposure time in sec',/, 1' * unobstructed area of main mirror in m**2',/, 1' * measured throughput of telescope/instrument',/, 1' * quantum efficiency of detector in the given band',/, 1'Accuracy +-20%. Signal-to-noise is per pixel step in',/, 1'wavelength for point sources,', 1' and per pixel for extended ones.',/, 1'For imaging, the counts obtained above are multiplied by the',/, 1'effective bandwidth of the filters in A, and signal-to-noise',/, 1'is calculated within a 2-FWHM-diameter aperture for a point',/, 1'source and per pixel for extended sources.') write(*,102) 102 format( 1'Signal-to-noise = Nobj/sqrt(Nobj+pixels*(Nsky+readnoise**2))',/, 1' ',/, 1'First time users, type H for more information.') end C****************************************************************************** C***** W R I T E O U T **************************************************** C****************************************************************************** C C subroutine writeout( 1 iw,iwb,iwd,iwi, 1 ninstr,ndet, 1 itel,imirror,itype,ngrating,iaxis,itrans, 1 idefdet,idefdet2,idefgrat,idefgrat2, 1 npixels,irot, 1 telarea,arcmm,defslit,spread,dispmm, 1 read,adu,pixel, 1 effmirror,wave,band,fluxjy,photons,skymag,extin, 1 effinstr,effgrat,effdet,truethru,zp,zmeas, 1 zmeasc,zmeasd, 1 version,cband,cinstr,comment,cgrating,cdet 1 ) implicit none integer i,j,k,kk,ninstr,ndet integer itrans(50),itel(50),itype(50),imirror(50),ngrating(50) integer iaxis(50),npixels(30,2),idefdet(50),idefdet2(50) integer iw,iwd(30),iwi(50),irot(30) integer iwb,idefgrat(50),idefgrat2(50),iid,iig real defslit(50),spread(50),telarea(50),effmirror real effinstr(10,50),effgrat(10,20,50),effdet(10,30) real zmeas(10,50,5),zp(10,50,5) real truethru(10,50),dispmm(20,50),arcmm(50) real pixel(30),read(30),adu(30) real band(10),fluxjy(10),photons(10),skymag(10), 1 extin(10),wave(10) character*2 cband(10) character*20 cgrating(20,50) character*25 cinstr(50) character*28 version character*40 cdet(30) character*80 comment(50),zmeasc(50),zmeasd(50) character*5 czmeas(10,50,5) do k=1,5 do j=1,50 do i=1,10 czmeas(i,j,k)=' ----' if(zmeas(i,j,k).gt.0.0) 1 write(czmeas(i,j,k),'(f5.1)') zmeas(i,j,k) enddo enddo enddo open(unit=11,name='outsum') write(11,1104) version,ninstr,ndet, 1 telarea(3),telarea(2),effmirror 1104 format( 1 'Instrument and detector parameters used by SIGNAL',/, 1 '=================================================',//, 1 a28,//, 1 'The counts and S:N predicted by SIGNAL are based on the ', 1 'numbers given below.',//, 1 'Total number of instruments in SIGNAL = ',i2,/, 1 'Total number of detectors in SIGNAL = ',i2,/, 1 '(Not all are shown below.)',//, 1 'Assumed collecting areas of WHT and INT (in m^2) = ', 1 2f5.1,/, 1 'Assumed reflectivity of aluminium telescope mirrors = ', 1 f5.2, 1 //,'EFFINSTR is the a priori predicted throughput ', 1 'of each', 1 ' instrument in each band',/,'(excluding atm, tel, gratings,', 1 ' filters).',/, 1 /,'TRUETHRU is an empirical correction factor,', 1 ' to make SIGNAL''s predictions',/,'match the observed', 1 ' instrumental zeropoints in each band.',/, 1 'For most instruments, the TRUETHRU values are close to 1.', 1 //,'''Assumed zeropoint'' is calculated on the basis of'/, 1 'EFFINSTR * grating throughput * detector QE * TRUETHRU.', 1 //'The measured zeropoints are empirical on-sky values.',/) write(11,1115) 1 (cband(k),k=1,iwb), 1 (int(wave(k)),k=1,iwb), 1 (int(band(k)),k=1,iwb), 1 (int(fluxjy(k)),k=1,iwb), 1 (int(photons(k)),k=1,iwb), 1 (skymag(k),k=1,iwb), 1 (extin(k),k=1,iwb) 1115 format(/,'Bandpasses used, sky brightness, extinction',/, 1 '-------------------------------------------',/, 1 'Band ',9(a2,4x),/, 1 'Wavel. (A) ',9(i5,1x),/, 1 'Bandwidth (A)',9(i5,1x),/, 1 'Mag 0 in Jy ',9(i5,1x),/, 1 'Mag 0 in ',9(i5,1x),/, 1 ' ph/s/A/cm^2',/, 1 'Dark sky mag ',9(f6.2),/, 1 ' /arcsec^2',/, 1 'Extinction ',9(f6.2),/, 1 ' (airm. 1.0) ',/) write(11,1132) 1132 format('Instruments included in SIGNAL',/, 1 '------------------------------',/) do j=1,ninstr write(11,1131) j,cinstr(j) 1131 format(i3,2x,a25) enddo write(11,1133) 1133 format(//,'Detectors included in SIGNAL',/, 1 '----------------------------',/) do j=1,ndet write(11,1134) j,cdet(j) 1134 format(i3,2x,a30) enddo do j=1,50 if(iwi(j).ne.0) then iw=iwi(j) write(11,1101) cinstr(iw),comment(iw)(1:70), 1 itel(iw),imirror(iw),itype(iw),ngrating(iw), 1 arcmm(iw),arcmm(iw)*pixel(idefdet(iw)), 1 defslit(iw),spread(iw),idefdet(iw), 1 cdet(idefdet(iw))(1:35),iaxis(iw),itrans(iw),iw 1101 format(/,15('====='),/,a25,/,15('====='),/,a70,/, 1 'ITEL ',i3,' (3 = WHT, 2 = INT, 1 = JKT)',/, 1 'IMIRROR ',i3,' (no. alum tel mirrors en route)',/, 1 'ITYPE ',i3, 1 ' (1 = imag, 2 = spec, 3 = ifu/fibre)',/, 1 'NGRATING ',i3,' (no. gratings available)',/, 1 'ARCMM ',f7.2,' "/mm (spatial scale)', 1 ' = ',f7.2,' "/pixel',/, 1 'DEFSLIT ',f7.2,' " (default slit width, for spec)',/, 1 'SPREAD ',f7.2,' (no. pix over which spec spread,', 1 ' for IFU)',/, 1 'IDEFDET ',i3,' (def. CCD) = ',a35,/, 1 'IAXIS ',i3,' (1 = spec disp is usu in x,', 1 ' 2 is usu in y)',/, 1 'ITRANS ',i3,' (order in menu)',/, 1 'IW ',i3,' (instrument number)',/) write(11,1102) (cband(k),k=1,iwb),(effinstr(k,iw),k=1,iwb) 1102 format('Throughput',/,'----------',/,13x,9(a2,3x), 1 ' A/mm A/pix',//, 1 'EFFINSTR ',9f5.2,/) do i=1,ngrating(iw) write(11,1103) cgrating(i,iw)(1:10), 1 (effgrat(k,i,iw),k=1,iwb), 1 dispmm(i,iw),dispmm(i,iw)*pixel(idefdet(iw)) 1103 format(a10,9f5.2,f8.1,f7.3) enddo write(11,1105) cdet(idefdet(iw))(1:6), 1 (effdet(k,idefdet(iw)),k=1,iwb) 1105 format(/,a6,4x,9f5.2) if(idefdet2(iw).ne.0) then write(11,1105) cdet(idefdet2(iw))(1:6), 1 (effdet(k,idefdet2(iw)),k=1,iwb) endif write(11,1109) 1 (truethru(k,iw),k=1,iwb) 1109 format(/,'Empirical / predicted thruput:',/, 1 'TRUETHRU ',9f5.2) write(11,1110) 1 (effinstr(k,iw)*truethru(k,iw),k=1,iwb) 1110 format(/,'EFFINSTR * empirical/predicted:',/, 1 10x,9f5.2) write(11,1111) cdet(idefdet(iw))(1:6), 1 (effdet(k,idefdet(iw))*truethru(k,iw)* 1 effinstr(k,iw),k=1,iwb) 1111 format(/,'EFFINSTR * empirical/prediced * QE of ',a6,':',/, 1 10x,9f5.2) C C Write out assumed zeropoints C do kk=1,2 C C Loop over up to 2 default detector (and/or grating) combinations C if(kk.eq.1) then iid=idefdet(iw) iig=idefgrat(iw) else if(kk.eq.2) then iid=idefdet2(iw) iig=idefgrat2(iw) endif if(iid.eq.0) goto 1357 if(itype(iw).eq.1) then write(11,1112) cdet(iid)(1:6),(zp(k,iw,kk),k=1,iwb) 1112 format(/,'Assumed zeropoint mags for 1 ph/sec ', 1 'with ',a6,':',/, 1 9x,9f5.1) else write(11,1113) cgrating(iig,iw)(1:10), 1 cdet(iid)(1:6),(zp(k,iw,kk),k=1,iwb) 1113 format(/,'Assumed zeropoint mags for 1 ph/sec/A ', 1 'with ',a10,' / ',a6,':',/, 1 9x,9f5.1) endif 1357 continue enddo C C Write out measured zeropints C if(itype(iw).eq.1) then write(11,1118) zmeasc(iw)(1:50), 1 (czmeas(k,iw,1),k=1,iwb),zmeasd(iw)(1:50) 1118 format(/,'Measured zeropoint mags for 1 ph/sec ', 1 'with ',/,10x,a50,/, 1 9x,9a5,/,10x,a50) else write(11,1119) zmeasc(iw)(1:50), 1 (czmeas(k,iw,1),k=1,iwb),zmeasd(iw)(1:50) 1119 format(/,'Measured zeropoint mags for 1 ph/sec/A ', 1 'with ',/,10x,a50,/, 1 9x,9a5,/,10x,a50) endif endif enddo C C Write summary of throughput parameters C write(11,1120) 1120 format(/,15('====='),/, 1 'Summary of instrument-throughput parameters',/, 1 '-------------------------------------------') do i=1,4 if(i.eq.1) write(11,1121) if(i.eq.2) write(11,1123) if(i.eq.3) write(11,1124) if(i.eq.4) write(11,1125) 1121 format(/, 1 'Predicted throughput EFFINSTR of instrument (excl atm, tel,', 1 ' gratings, filters):',/) 1123 format(/,'Observed / predicted throughput TRUETHRU:',/) 1124 format(/,'Actual instrument throughput EFFINSTR * TRUETHRU:',/) 1125 format(/,'Actual instrument throughput * detector QE:',/) do j=1,50 if(iwi(j).ne.0) then iw=iwi(j) if(i.eq.1) write(11,1122) cinstr(iw)(1:20), 1 (effinstr(k,iw),k=1,iwb) 1122 format(a20,9f5.2) if(i.eq.2) 1 write(11,1122) cinstr(iw)(1:20),(truethru(k,iw),k=1,iwb) if(i.eq.3) write(11,1122) cinstr(iw)(1:20), 1 (effinstr(k,iw)*truethru(k,iw),k=1,iwb) if(i.eq.4) write(11,1126) cinstr(iw)(1:20), 1 (effinstr(k,iw)*truethru(k,iw)* 1 effdet(k,idefdet(iw)),k=1,iwb),cdet(idefdet(iw))(1:6) 1126 format(a20,9f5.2,2x,a6) endif enddo enddo C C Write out CCD QE C write(11,1106) (cband(k),k=1,9) 1106 format(/,15('====='),//,'CCD parameters',/, 1 '--------------',/, 1 13x,' CCD QE in each band',/, 1 13x,8(a2,3x),a2) do j=1,30 if(iwd(j).ne.0) then iw=iwd(j) write(11,1107) cdet(iw)(1:6), 1 (effdet(k,iw),k=1,9) 1107 format(a6,4x,9f5.2) endif enddo C C Write out other CCD parameters C write(11,1108) 1108 format(/, 1 'CCD Read noise Gain Pix size Pixels IROT',/, 1 ' e- e-/ADU microns x * y ',/) do j=1,30 if(iwd(j).ne.0) then iw=iwd(j) write(11,1116) cdet(iw)(1:6), 1 read(iw),adu(iw),pixel(iw), 1 npixels(iw,1),npixels(iw,2),irot(iw) 1116 format(a6,1x,f5.1,f12.1,f10.4,i10,i5,i3) endif enddo end