diff options
Diffstat (limited to 'noao/obsutil/src/ccdtime')
-rw-r--r-- | noao/obsutil/src/ccdtime/Revisions | 79 | ||||
-rw-r--r-- | noao/obsutil/src/ccdtime/ccddb.x | 222 | ||||
-rw-r--r-- | noao/obsutil/src/ccdtime/ccdtime.par | 15 | ||||
-rw-r--r-- | noao/obsutil/src/ccdtime/mkpkg | 17 | ||||
-rw-r--r-- | noao/obsutil/src/ccdtime/t_ccdtime.x | 307 | ||||
-rw-r--r-- | noao/obsutil/src/ccdtime/x_ccdtime.x | 1 |
6 files changed, 641 insertions, 0 deletions
diff --git a/noao/obsutil/src/ccdtime/Revisions b/noao/obsutil/src/ccdtime/Revisions new file mode 100644 index 00000000..44405b23 --- /dev/null +++ b/noao/obsutil/src/ccdtime/Revisions @@ -0,0 +1,79 @@ +.help revisions Jun88 noao.obsutil +.nf + +======= +V2.12.2 +======= + +ccdtime.par +t_ccdtime.x + 1. The minimum seeing is now 0.001. + 2. The formating of the seeing was changed in case the seeing is set + less than 0.1. + (4/24/03, Valdes) + +======= +V2.12.1 +======= + +===== +V2.12 +===== + +t_ccdtime.x +../doc/ccdtime.hlp + It is now an error if time<0, time>10-000, abs(mag)>40, snr<0 or + snr>100000. (8/24/00, Valdes) + +../doc/ccdtime.hlp + In the formula for r(sky) was pixel area term was in the wrong place. + (3/9/99, Valdes) + +t_ccdtime.x + For the case where SNR is very large and a time is specified the + iteration on the magnitude might not complete. The iteration is now + capped at 100 and the test for convergence is now normalized. + (11/6/98, Valdes) + +t_ccdtime.x +../doc/ccdtime.hlp + 1. The calculation of exposure time given a SNR was changed from an + interative solution to an analytic solution. + 2. The times are printed to 0.01s. + 3. The photometry aperture is now the rounded-up integer with a minimum + of 9 pixels. + (9/8/98, Valdes) + +t_ccdtime.x +ccddb.x +../doc/ccdtime.hlp + 1. The database keywords can now be index by reference to the telescope, + filter, and/or telescope. + 2. A new filter keyword, "extinction", was added to specify the + extinction. + 3. The extinction is now used to fixe the previous incorrect behavior + that used 1 mag/airmass extinction. The old results are preserved + by making the default extinction be 1 if missing. However the + database files should be updated to have correct extinctions. + (8/19/98, Valdes) + +ccddb.x +../doc/ccdtime.hlp + 1. The code would not work with database entries containing whitespace. + 2. The help was not correct in describing how the number of pixels used + in the photometry is calculated from the seeing FWHM. + (4/5/94, Valdes) + +t_ccdtime.x + Modified CCDTIME to use a plate scale instead of the f/ratio and to + include an airmass term. (10/23/93, Valdes) + +t_ccdtime.x + +ccddb.x + +ccdtime.x - +ccdtime.par +doc/ccdtime.hlp + Revised CCDTIME to use a telescope/filter/detector database and to + compute and print additional information. (8/16/93, Valdes) + +.endhelp diff --git a/noao/obsutil/src/ccdtime/ccddb.x b/noao/obsutil/src/ccdtime/ccddb.x new file mode 100644 index 00000000..e0f4dd1d --- /dev/null +++ b/noao/obsutil/src/ccdtime/ccddb.x @@ -0,0 +1,222 @@ +# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc. + +# General text database routines. + +# Symbol table definitions. +define LEN_INDEX 10 # Length of symtab index +define LEN_STAB 512 # Length of symtab +define SZ_SBUF 512 # Size of symtab string buffer +define SYMLEN 40 # Length of symbol structure +define SZ_DBVAL 79 # Size of database value string + +# Symbol table structure +define DBVAL Memc[P2C($1)] # Database value string + + +# DBOPEN -- Open database and store the requested information in symbol table. + +pointer procedure dbopen (dname, fname, kname, ename) + +char dname[ARB] #I Directory name +char fname[ARB] #I File name +char kname[ARB] #I Key name +char ename[ARB] #I Entry name +pointer db #O Database symbol table pointer + +int fd, found, open(), fscan(), nscan() +pointer sp, pname, name, key, str, sym +pointer stopen(), stenter() +bool streq(), strne() +errchk open, stopen, stenter, fscan, dberror + +begin + call smark (sp) + call salloc (pname, SZ_FNAME, TY_CHAR) + call salloc (name, SZ_LINE, TY_CHAR) + call salloc (key, SZ_FNAME, TY_CHAR) + call salloc (str, SZ_LINE, TY_CHAR) + + # Open database. + call sprintf (Memc[pname], SZ_FNAME, "%s%s") + call pargstr (dname) + call pargstr (fname) + fd = open (Memc[pname], READ_ONLY, TEXT_FILE) + + # Strip entry name whitespace and convert to lower case. + call strcpy (ename, Memc[name], SZ_LINE) + call xt_stripwhite (Memc[name]) + call strlwr (Memc[name]) + + # List entries in database. + if (Memc[name] == '?') { + Call printf ("Entries for %s in database %s:\n") + call pargstr (kname) + call pargstr (Memc[pname]) + while (fscan (fd) != EOF) { + call gargwrd (Memc[key], SZ_FNAME) + call gargwrd (Memc[str], SZ_LINE) + call gargwrd (Memc[str], SZ_LINE) + if (nscan()<3 || Memc[key]=='#' || strne (Memc[key], kname)) + next + call printf ("\t%s\n") + call pargstr (Memc[str]) + } + call close (fd) + call sfree (sp) + return (NULL) + } + + # Find entry. + found = 0 + while (fscan (fd) != EOF) { + call gargwrd (Memc[key], SZ_FNAME) + call gargwrd (Memc[str], SZ_LINE) + call gargwrd (Memc[str], SZ_LINE) + if (nscan()<3 || Memc[key]=='#' || strne (Memc[key], kname)) + next + found = 1 + call strlwr (Memc[str]) + if (streq (Memc[str], Memc[name])) { + found = 2 + break + } + } + + # Check if entry was found. + if (found != 2) { + call close (fd) + if (found != 1) + call dberror ("DBOPEN: Database entry not found", kname) + else + call dberror ("DBOPEN: Database entry not found", ename) + } + + # Create symbol table. + db = stopen (ename, LEN_INDEX, LEN_STAB, SZ_SBUF) + + # Read the file and enter the parameters in the symbol table. + sym = stenter (db, Memc[key], SYMLEN) + call strcpy (ename, DBVAL(sym), SZ_DBVAL) + while (fscan(fd) != EOF) { + call gargwrd (Memc[key], SZ_FNAME) + call gargwrd (Memc[str], SZ_LINE) + call gargwrd (Memc[str], SZ_LINE) + if (nscan()>0 && (streq(Memc[key],"end") || streq(Memc[key],kname))) + break + if (nscan() < 3 || Memc[key] == '#') + next + sym = stenter (db, Memc[key], SYMLEN) + call strcpy (Memc[str], DBVAL(sym), SZ_DBVAL) + } + + call close (fd) + call sfree (sp) + + return (db) +end + + +# DBCLOSE -- Close the database symbol table pointer. + +procedure dbclose (db) + +pointer db # Database symbol table pointer + +begin + if (db != NULL) + call stclose (db) +end + + +# DBGETD -- Get double database parameter. + +double procedure dbgetd (db, param, arg1, arg2) + +pointer db # Database symbol table pointer +char param[ARB] # Database parameter +char arg1[ARB], arg2[ARB] # Optional arguments + +char str[SZ_LINE] +int ip, ctod() +double dval +errchk dbgstr + +begin + call dbgstr (db, param, arg1, arg2, str, SZ_LINE) + + ip = 1 + if (ctod (str, ip, dval) <= 0) + call dberror ("DBGETD: Database parameter not double", param) + return (dval) +end + + +# DBGSTR -- Get string valued parameter. + +procedure dbgstr (db, param, arg1, arg2, str, maxchar) + +pointer db # Database symbol table pointer +char param[ARB] # Database parameter +char arg1[ARB], arg2[ARB] # Optional arguments +char str[maxchar] # Database parameter value +int maxchar # Maximum characters for string + +pointer sp, param1, sym, stfind() +errchk dberror + +begin + call smark (sp) + call salloc (param1, SZ_LINE, TY_CHAR) + + sym = NULL + if (arg1[1] != EOS && arg2[1] != EOS) { + call sprintf (Memc[param1], SZ_LINE, "%s(%s,%s)") + call pargstr (param) + call pargstr (arg1) + call pargstr (arg2) + sym = stfind (db, Memc[param1]) + if (sym == NULL) { + call sprintf (Memc[param1], SZ_LINE, "%s(%s,%s)") + call pargstr (param) + call pargstr (arg2) + call pargstr (arg1) + sym = stfind (db, Memc[param1]) + } + } + if (sym == NULL && arg1[1] != EOS) { + call sprintf (Memc[param1], SZ_LINE, "%s(%s)") + call pargstr (param) + call pargstr (arg1) + sym = stfind (db, Memc[param1]) + } + if (sym == NULL && arg2[1] != EOS) { + call sprintf (Memc[param1], SZ_LINE, "%s(%s)") + call pargstr (param) + call pargstr (arg2) + sym = stfind (db, Memc[param1]) + } + if (sym == NULL) + sym = stfind (db, param) + + call sfree (sp) + + if (sym == NULL) + call dberror ("DBGSTR: Database parameter not found", param) + call strcpy (DBVAL(sym), str, maxchar) +end + + +# DBERROR -- Print database error. + +procedure dberror (errstr, param) + +char errstr[ARB] # Error string +char param[ARB] # Parameter +char errmsg[SZ_LINE] # Error message + +begin + call sprintf (errmsg, SZ_LINE, "%s (%s)") + call pargstr (errstr) + call pargstr (param) + call error (1, errmsg) +end diff --git a/noao/obsutil/src/ccdtime/ccdtime.par b/noao/obsutil/src/ccdtime/ccdtime.par new file mode 100644 index 00000000..d7a36c24 --- /dev/null +++ b/noao/obsutil/src/ccdtime/ccdtime.par @@ -0,0 +1,15 @@ +time,r,h,INDEF,,,Time (seconds) +magnitude,r,h,20.,,,Magnitude +snr,r,h,10.,0.1,,Signal-to-noise ratio +database,s,h,"ccdtime$kpno.dat",,,Database file +telescope,s,h,,,,Telescope name (? for list) +detector,s,h,,,,Detector name (? for list) +sum,i,h,1,,,CCD summing factor +seeing,r,h,1.5,0.001,10.0,Seeing (arcsec) +airmass,r,h,1.2,1.,,Airmass +phase,r,h,0.,0.,28.,Moon phase (0-28) +f1,s,h,U,,,Filter 1 +f2,s,h,B,,,Filter 2 +f3,s,h,V,,,Filter 3 +f4,s,h,R,,,Filter 4 +f5,s,h,I,,,Filter 5 diff --git a/noao/obsutil/src/ccdtime/mkpkg b/noao/obsutil/src/ccdtime/mkpkg new file mode 100644 index 00000000..d468ef93 --- /dev/null +++ b/noao/obsutil/src/ccdtime/mkpkg @@ -0,0 +1,17 @@ +# Make the CCDTIME task. + +$checkout libpkg.a ../ +$update libpkg.a +$checkin libpkg.a ../ +$exit + +standalone: + $update libpkg.a + $omake x_ccdtime.x + $link x_ccdtime.o libpkg.a -lxtools -o xx_ccdtime.e + ; + +libpkg.a: + ccddb.x + t_ccdtime.x + ; diff --git a/noao/obsutil/src/ccdtime/t_ccdtime.x b/noao/obsutil/src/ccdtime/t_ccdtime.x new file mode 100644 index 00000000..0c0e8822 --- /dev/null +++ b/noao/obsutil/src/ccdtime/t_ccdtime.x @@ -0,0 +1,307 @@ +define NF 25 # Maximum number of filters +define SZ_FILTER 7 # Maximum length of filter name +define TOL 0.0001 # Convergence tolerance + + +# T_CCDTIME -- Compute the time, magnitude, and signal-to-noise for CCD +# exposures for a given telescope, detector, and filters. The telescope +# detector, and filter data come from a text database. The computed +# quantities fix two parameters and compute the third. + +procedure t_ccdtime() + +pointer database # Telescope, detector, filter database +pointer telescope # Telescope +pointer detector # Detector +pointer fltr[5] # Filters +real time # Target time (sec) +real mag # Target magnitude +real snr # Target SNR +int sum # CCD summing +real seeing # Seeing (arcsec) +real phase # Moon phase (0-28) +real airmass # Airmass + +int i, j, k, nf, n +real a, b, c +real aper, scale, trans, nread, dark, pixsize, lum, p +real star[NF], sky[NF], t, m, s, nstar, nsky, ndark, noise, npix +real dqe, ext, starmag, counts, sky0, sky1, sky2 +pointer sp, tdb, ddb, filter[NF], fdb[NF] + +int clgeti() +real clgetr() +double dbgetd() +pointer dbopen() +errchk dbopen + +define done_ 10 + +begin + call smark (sp) + call salloc (database, SZ_FNAME, TY_CHAR) + call salloc (telescope, SZ_FNAME, TY_CHAR) + call salloc (detector, SZ_FNAME, TY_CHAR) + do i = 1, 5 + call salloc (fltr[i], SZ_LINE, TY_CHAR) + do i = 1, NF + call salloc (filter[i], SZ_FILTER, TY_CHAR) + + # Get task parameters. + call clgstr ("database", Memc[database], SZ_FNAME) + call clgstr ("telescope", Memc[telescope], SZ_FNAME) + call clgstr ("detector", Memc[detector], SZ_FNAME) + call clgstr ("f1", Memc[fltr[1]], SZ_LINE) + call clgstr ("f2", Memc[fltr[2]], SZ_LINE) + call clgstr ("f3", Memc[fltr[3]], SZ_LINE) + call clgstr ("f4", Memc[fltr[4]], SZ_LINE) + call clgstr ("f5", Memc[fltr[5]], SZ_LINE) + + time = clgetr ("time") + mag = clgetr ("magnitude") + snr = clgetr ("snr") + sum = clgeti ("sum") + seeing = clgetr ("seeing") + phase = clgetr ("phase") + airmass = clgetr ("airmass") + + # The input filter strings may be lists of filters and here + # we expand them into all filters. + + nf = 1; k = 0 + do i = 1, 5 { + for (j = fltr[i]; Memc[j] != EOS; j = j + 1) { + Memc[filter[nf]+k] = Memc[j] + if (Memc[j] == ',') { + Memc[filter[nf]+k] = EOS + if (k > 0) { + nf = nf + 1; k = 0 + } + } else + k = k + 1 + } + Memc[filter[nf]+k] = EOS + if (k > 0) { + nf = nf + 1; k = 0 + } + } + nf = nf - 1 + + i = 0 + if (IS_INDEFR(time)) + i = i + 1 + else if (time <= 0.) + call error (1, "Requested time must be greater than zero") + else if (time > 100000.) + call error (1, "Requested time must be less than 100,000") + + if (IS_INDEFR(mag)) + i = i + 1 + else if (mag > 40.) + call error (1, "Requested magnitude must be less than 40") + else if (mag < -40.) + call error (1, "Requested magnitude must be greater than -40") + + if (IS_INDEFR(snr)) + i = i + 1 + else if (snr <= 0.) + call error (1, "Requested SNR must be greater than zero") + else if (snr > 100000.) + call error (1, "Requested SNR must be less than 100,000") + + if (i > 1) { + call sfree (sp) + call error (1, + "At least two of time, magnitude, and snr must be specified") + } + + if (phase < 14) + p = phase + else + p = 28 - phase + + # Open database entries. + # If '?' this will print list and then the task will exit. + # If an error occurs in the telescope or detector abort. + # If an error occurs in the filter ignore the filter. + + tdb = dbopen ("", Memc[database], "telescope", Memc[telescope]) + ddb = dbopen ("", Memc[database], "detector", Memc[detector]) + + n = 0 + do i = 1, nf { + iferr (fdb[n+1]=dbopen("",Memc[database],"filter",Memc[filter[i]])) + next + if (fdb[n+1] == NULL) + goto done_ + n = n + 1 + filter[n] = filter[i] + } + if (tdb == NULL || ddb == NULL || n == 0) + goto done_ + + # Get star and sky rates for telescope/detector/filter combination. + # Convert to a standard rate at 20th magnitude at the given airmass. + + do i = 1, n { + # Get telescope parameters. + aper = dbgetd (tdb, "aperture", Memc[filter[i]], Memc[detector]) + scale = dbgetd (tdb, "scale", Memc[filter[i]], Memc[detector]) + trans = dbgetd (tdb, "transmission", Memc[filter[i]], + Memc[detector]) + + # Get detector parameters. + dqe = dbgetd (ddb, Memc[filter[i]], Memc[filter[i]], + Memc[telescope]) + nread = dbgetd (ddb, "rdnoise", Memc[filter[i]], Memc[telescope]) + dark = dbgetd (ddb, "dark", Memc[filter[i]], Memc[telescope]) + pixsize = dbgetd (ddb,"pixsize",Memc[filter[i]],Memc[telescope]) * + (scale / 1000) * sum + npix = max (9, nint (1.4 * (seeing / pixsize)**2+0.5)) + + # Get filter parameters. + iferr (ext = dbgetd (fdb[i], "extinction", Memc[telescope], + Memc[detector])) + ext = 1 + starmag = dbgetd (fdb[i], "mag", Memc[telescope], Memc[detector]) + counts = dbgetd (fdb[i], "star", Memc[telescope], Memc[detector]) + sky0 = dbgetd (fdb[i], "sky0", Memc[telescope], Memc[detector]) + sky1 = dbgetd (fdb[i], "sky1", Memc[telescope], Memc[detector]) + sky2 = dbgetd (fdb[i], "sky2", Memc[telescope], Memc[detector]) + + lum = 10. ** (0.4 * (starmag - 20.)) + star[i] = counts * lum * aper ** 2 * trans * + 10 ** (0.4 * (1 - airmass) * ext) + star[i] = star[i] * dqe + sky[i] = sky0 + sky1 * p + sky2 * p * p + sky[i] = star[i] * pixsize ** 2 * 10. ** ((20. - sky[i]) / 2.5) + } + if (!IS_INDEFR(mag)) + lum = 10. ** (0.4 * (20. - mag)) + + # Print header and column labels. + call printf ("Database: %-20s Telescope: %-10s Detector: %-10s\n") + call pargstr (Memc[database]) + call pargstr (Memc[telescope]) + call pargstr (Memc[detector]) + call printf (" Sum: %-2d Arcsec/pixel: %-4.2f Pixels/star: %-4.1f\n") + call pargi (sum) + call pargr (pixsize) + call pargr (npix) + call printf (" Seeing: %-.3g Airmass: %-4.2f Phase: %-4.1f\n") + call pargr (seeing) + call pargr (airmass) + call pargr (phase) + call printf ("\n%7s %7s %7s %7s %7s %7s %s\n") + call pargstr ("Filter") + call pargstr ("Time") + call pargstr ("Mag") + call pargstr ("SNR") + call pargstr ("Star") + call pargstr ("Sky/pix") + call pargstr (" Noise contributions") + call printf ("%47w %7s %7s %7s\n") + call pargstr ("Star") + call pargstr ("Sky") + call pargstr ("CCD") + + # Compute exposure time through each filter. + + if (!IS_INDEFR(mag) && !IS_INDEFR(snr)) { + call printf ("\n") + do i = 1, n { + a = ((star[i] * lum) / snr) ** 2 + b = -(star[i] * lum + npix * (sky[i] + dark)) + c = -npix * nread ** 2 + t = (-b + sqrt (b**2 - 4 * a * c)) / (2 * a) + + nstar = star[i] * lum * t + nsky = sky[i] * t + ndark = dark * t + noise = sqrt(nstar + npix * (nsky + ndark + nread**2)) + s = nstar / noise + m = mag + + call printf ( + "%7s %7.2f %7.1f %7.1f %7.1f %7.1f %7.2f %7.2f %7.2f\n") + call pargstr (Memc[filter[i]]) + call pargr (t) + call pargr (m) + call pargr (s) + call pargr (nstar) + call pargr (nsky) + call pargr (sqrt (nstar)) + call pargr (sqrt (npix * nsky)) + call pargr (sqrt (npix * (ndark + nread**2))) + } + } + + # Compute magnitude through each filter. + # Use resubstitution to iterate for SNR. + + if (!IS_INDEFR(time) && !IS_INDEFR(snr)) { + call printf ("\n") + do i = 1, n { + m = 20 + s = 0 + do j = 1, 100 { + t = time + nstar = star[i] * 10**(0.4*(20.0-m)) * t + nsky = sky[i] * t + ndark = dark * t + noise = sqrt(nstar + npix * (nsky + ndark + nread**2)) + m = 20 - 2.5 * log10 (snr * noise / (t * star[i])) + s = nstar / noise + if (abs(1-s/snr) <= TOL) + break + } + + call printf ( + "%7s %7.2f %7.1f %7.1f %7.1f %7.1f %7.2f %7.2f %7.2f\n") + call pargstr (Memc[filter[i]]) + call pargr (t) + call pargr (m) + call pargr (s) + call pargr (nstar) + call pargr (nsky) + call pargr (sqrt (nstar)) + call pargr (sqrt (npix * nsky)) + call pargr (sqrt (npix * (ndark + nread**2))) + } + } + + # Compute SNR through each filter. + + if (!IS_INDEFR(time) && !IS_INDEFR(mag)) { + call printf ("\n") + do i = 1, n { + t = time + m = mag + nstar = star[i] * lum * t + nsky = sky[i] * t + ndark = dark * t + noise = sqrt(nstar + npix * (nsky + ndark + nread**2)) + s = nstar / noise + + call printf ( + "%7s %7.2f %7.1f %7.1f %7.1f %7.1f %7.2f %7.2f %7.2f\n") + call pargstr (Memc[filter[i]]) + call pargr (t) + call pargr (m) + call pargr (s) + call pargr (nstar) + call pargr (nsky) + call pargr (sqrt (nstar)) + call pargr (sqrt (npix * nsky)) + call pargr (sqrt (npix * (ndark + nread**2))) + } + } + call printf ("\n") + +done_ + call dbclose (tdb) + call dbclose (ddb) + do i = 1, n + call dbclose (fdb[i]) + call sfree (sp) +end diff --git a/noao/obsutil/src/ccdtime/x_ccdtime.x b/noao/obsutil/src/ccdtime/x_ccdtime.x new file mode 100644 index 00000000..cccfa908 --- /dev/null +++ b/noao/obsutil/src/ccdtime/x_ccdtime.x @@ -0,0 +1 @@ +task ccdtime = t_ccdtime |