diff options
Diffstat (limited to 'noao/imred/src/fibers/arcrefs.cl')
| -rw-r--r-- | noao/imred/src/fibers/arcrefs.cl | 326 |
1 files changed, 326 insertions, 0 deletions
diff --git a/noao/imred/src/fibers/arcrefs.cl b/noao/imred/src/fibers/arcrefs.cl new file mode 100644 index 00000000..d0d6fcf1 --- /dev/null +++ b/noao/imred/src/fibers/arcrefs.cl @@ -0,0 +1,326 @@ +# ARCREFS -- Determine dispersion relation for reference arcs. + +procedure arcrefs (arcref1, arcref2, extn, arcreplace, apidtable, response, + crval, cdelt, done, log1, log2) + +file arcref1 +file arcref2 +string extn +file arcreplace +file apidtable +file response +string crval = "INDEF" +string cdelt = "INDEF" +file done +file log1 +file log2 + +struct *fd + +begin + string imtype + string arcref, arcrefms, arc, arcms + string temp, temp1, temp2, str1, str2 + int i, n, nspec, dc + real w + bool log + struct str3 + + imtype = "." // envget ("imtype") + i = stridx (",", imtype) + if (i > 0) + imtype = substr (imtype, 1, i-1) + n = strlen (imtype) + + temp = mktemp ("tmp$iraf") + temp1 = mktemp ("tmp$iraf") + temp2 = mktemp ("tmp$iraf") + + # Extract the primary arc reference spectrum. Extract and replace + # any replacement arcs defined in the arcreplace file. Determine the + # dispersion function with IDENTIFY/REIDENTIFY. Set the wavelength + # parameters with MSDISPCOR. + + arcref = arcref1 + arcrefms = arcref1 // extn + if (!access (arcrefms//imtype)) { + print ("Extract arc reference image ", arcref) | tee (log1) + apscript (arcref, output=arcrefms, ansrecenter="NO", + ansresize="NO", ansedit="NO", anstrace="NO", + nsubaps=params.nsubaps, background="none", clean=no, + weights="none") + sapertures (arcrefms, apertures="", apidtable=apidtable, + wcsreset=no, verbose=no, beam=INDEF, dtype=INDEF, w1=INDEF, + dw=INDEF, z=INDEF, aplow=INDEF, aphigh=INDEF, title=INDEF) + if (response != "") { + if (params.nsubaps == 1) + sarith (arcrefms, "/", response, arcrefms, w1=INDEF, + w2=INDEF, apertures="", bands="", beams="", apmodulus=0, + reverse=no, ignoreaps=no, format="multispec", + renumber=no, offset=0, clobber=yes, merge=no, + errval=0, verbose=no) + else { + blkrep (response, temp, 1, params.nsubaps) + sarith (arcrefms, "/", temp, arcrefms, w1=INDEF, w2=INDEF, + apertures="", bands="", beams="", apmodulus=0, + reverse=no, ignoreaps=yes, format="multispec", + renumber=no, offset=0, clobber=yes, merge=no, + errval=0, verbose=no) + imdelete (temp, verify=no) + } + } + + if (arcreplace != "") { + if (!access (arcreplace)) + error (1, "Can't access file "//arcreplace) + fd = arcreplace + while (fscan (fd, arc, str1, str2) != EOF) { + i = strlen (arc) + if (i > n && substr (arc, i-n+1, i) == imtype) + arc = substr (arc, 1, i-n) + if (arc != arcref) + next + arc = str1 + if (i > n && substr (arc, i-n+1, i) == imtype) + arc = substr (arc, 1, i-n) + arcms = arc // extn + + if (access (arcms//imtype)) + imdelete (arcms, verify=no) + + print ("Extract arc reference image ", arc) | tee (log1) + apscript (arc, output=arcms, ansrecenter="NO", + ansresize="NO", ansedit="NO", anstrace="NO", + nsubaps=params.nsubaps, background="none", clean=no, + weights="none") + sapertures (arcms, apertures="", apidtable=apidtable, + wcsreset=no, verbose=no, beam=INDEF, dtype=INDEF, + w1=INDEF, dw=INDEF, z=INDEF, aplow=INDEF, aphigh=INDEF, + title=INDEF) + if (response != "") { + if (params.nsubaps == 1) + sarith (arcms, "/", response, arcms, w1=INDEF, + w2=INDEF, apertures="", bands="", beams="", + apmodulus=0, reverse=no, ignoreaps=no, + format="multispec", renumber=no, offset=0, + clobber=yes, merge=no, errval=0, verbose=no) + else { + blkrep (response, temp, 1, params.nsubaps) + sarith (arcms, "/", temp, arcms, w1=INDEF, + w2=INDEF, apertures="", bands="", beams="", + apmodulus=0, reverse=no, ignoreaps=yes, + format="multispec", renumber=no, offset=0, + clobber=yes, merge=no, errval=0, verbose=no) + imdelete (temp, verify=no) + } + } + scopy (arcms, arcrefms, w1=INDEF, w2=INDEF, apertures=str2, + bands="", beams="", apmodulus=1000, offset=0, + format="multispec", clobber=yes, merge=yes, renumber=no, + verbose=yes) | tee (log1, > log2) + imdelete (arcms, verify=no) + } + fd = "" + } + } + + # Get the dispersion parameters from the header. These are + # used for all further spectra and also flag whether this + # spectrum has been processed. If the parameters are missing + # the spectrum needs to have the dispersion function and + # wavelength scale determined. The HEDIT is needed because + # in some cases the user may exit IDENTIFY without updating + # the database (if the image was deleted but the database + # entry was not). + + hselect (arcrefms, "dclog1", yes) | scan (str1) + if (nscan () != 1) { + print ("Determine dispersion solution for ", arcref) | tee (log1) + #delete (database//"/id"//arcrefms//"*", verify=no) + printf ("%s %s\n", crval, cdelt) | scan (str3) + if (str3 == "INDEF INDEF") + identify (arcrefms, section="middle line", database=database, + coordlist=params.coordlist, nsum=1, match=params.match, + maxfeatures=50, zwidth=100., ftype="emission", + fwidth=params.fwidth, cradius=params.cradius, + threshold=params.threshold, minsep=2., + function=params.i_function, order=params.i_order, + sample="*", niterate=params.i_niterate, + low_reject=params.i_low, high_reject=params.i_high, + grow=0., autowrite=yes) + else + autoidentify (arcrefms, crval, cdelt, + coordlist=params.coordlist, + interactive="YES", section="middle line", nsum="1", + ftype="emission", fwidth=params.fwidth, + cradius=params.cradius, threshold=params.threshold, + minsep=2., match=params.match, function=params.i_function, + order=params.i_order, sample="*", + niterate=params.i_niterate, low_reject=params.i_low, + high_reject=params.i_high, grow=0., dbwrite="YES", + overwrite=yes, database="database", verbose=yes, + logfile=logfile, plotfile=plotfile, + reflist="", refspec="", crpix="INDEF", cddir="unknown", + crsearch="-0.5", cdsearch="INDEF", aidpars="") + + hedit (arcrefms, "refspec1", arcrefms, add=yes, + show=no, verify=no) + + nspec = 1 + hselect (arcrefms, "naxis2", yes) | scan (nspec) + if (nspec > 1) + reidentify (arcrefms, "", interactive=yes, + section="middle line", shift=0., step=1, nsum=1, + cradius=params.cradius, threshold=params.threshold, + nlost=100, refit=params.refit, trace=no, override=yes, + addfeatures=params.addfeatures, newaps=no, + coordlist=params.coordlist, match=params.match, + maxfeatures=50, minsep=2., database=database, + plotfile=plotfile, logfiles=logfile, verbose=yes) + + # Dispersion correct the reference arc. This step is required to + # use the confirm option of MSDISPCOR to set the wavelength scale + # for all further spectra. Set the newdisp flag. + + print ("Dispersion correct ", arcref) | tee (log1) + dispcor (arcrefms, "", linearize=params.linearize, + database=database, table="", w1=INDEF, w2=INDEF, dw=INDEF, + nw=INDEF, log=params.log, flux=params.flux, samedisp=yes, + global=no, ignoreaps=no, confirm=yes, listonly=no, verbose=no, + logfile=logfile) + if (params.nsubaps > 1) { + imrename (arcrefms, temp, verbose=no) + scopy (temp, arcrefms, w1=INDEF, w2=INDEF, apertures="1-999", + bands="", beams="", apmodulus=0, offset=0, + format="multispec", clobber=no, merge=no, renumber=no, + verbose=no) + blkavg (temp, temp, 1, params.nsubaps, option="sum") + imcopy (temp, arcrefms//"[*,*]", verbose=no) + imdelete (temp, verify=no) + } + proc.newdisp = yes + } + if (extn == ".ms") + print (arcref, >> done) + + # Extract the alternate shift arc reference. Transfer the dispersion + # function from the primary arc reference and then identify shift + # lines. + + if (arcref2 != "") { + arcref = arcref2 + arcrefms = arcref2 // extn + if (proc.newdisp && access (arcrefms//imtype)) + imdelete (arcrefms, verify=no) + if (!access (arcrefms)) { + print ("Extract arc reference image ", arcref) | tee (log1) + apscript (arcref, output=arcrefms, ansrecenter="NO", + ansresize="NO", ansedit="NO", anstrace="NO", + nsubaps=params.nsubaps, background="none", clean=no, + weights="none") + sapertures (arcrefms, apertures="", apidtable=apidtable, + wcsreset=no, verbose=no, beam=INDEF, dtype=INDEF, w1=INDEF, + dw=INDEF, z=INDEF, aplow=INDEF, aphigh=INDEF, title=INDEF) + if (response != "") { + if (params.nsubaps == 1) + sarith (arcrefms, "/", response, arcrefms, w1=INDEF, + w2=INDEF, apertures="", bands="", beams="", + apmodulus=0, reverse=no, ignoreaps=no, + format="multispec", renumber=no, offset=0, + clobber=yes, merge=no, errval=0, verbose=no) + else { + blkrep (response, temp, 1, params.nsubaps) + sarith (arcrefms, "/", temp, arcrefms, w1=INDEF, + w2=INDEF, apertures="", bands="", beams="", + apmodulus=0, reverse=no, ignoreaps=yes, + format="multispec", renumber=no, offset=0, + clobber=yes, merge=no, errval=0, verbose=no) + imdelete (temp, verify=no) + } + } + } + + hselect (arcrefms, "dclog1", yes) | scan (str1) + if (nscan () != 1) { + print ("Determine dispersion solution for ", arcref) | + tee (log1) + #delete (database//"/id"//arcrefms//"*", verify=no) + + print (":r ", arcref1//extn, "\na\nd") | + identify (arcrefms, section="middle line", database=database, + coordlist="", nsum=1, match=params.match, maxfeatures=50, + zwidth=100., ftype="emission", fwidth=params.fwidth, + cradius=params.cradius, threshold=params.threshold, + minsep=2., function=params.i_function, + order=params.i_order, sample="*", + niterate=params.i_niterate, low_reject=params.i_low, + high_reject=params.i_high, grow=0., autowrite=yes, + cursor="STDIN", >G "dev$null", >& "dev$null") + identify (arcrefms, section="middle line", database=database, + coordlist="", nsum=1, match=params.match, maxfeatures=50, + zwidth=100., ftype="emission", fwidth=params.fwidth, + cradius=params.cradius, threshold=params.threshold, + minsep=2., function=params.i_function, + order=params.i_order, sample="*", + niterate=params.i_niterate, low_reject=params.i_low, + high_reject=params.i_high, grow=0., autowrite=yes) + print (":feat ", temp) | + identify (arcrefms, section="middle line", database=database, + coordlist="", nsum=1, match=params.match, maxfeatures=50, + zwidth=100., ftype="emission", fwidth=params.fwidth, + cradius=params.cradius, threshold=params.threshold, + minsep=2., function=params.i_function, + order=params.i_order, sample="*", + niterate=params.i_niterate, low_reject=params.i_low, + high_reject=params.i_high, grow=0., autowrite=yes, + cursor="STDIN", >G "dev$null", >& "dev$null") + print (":r ", arcref1//extn, "\na\nd", > temp1) + fd = temp + while (fscan (fd, i, w, w, w) != EOF) { + if (nscan() == 4) { + print (w, 1, 1, "m", >> temp1) + print (w, >> temp2) + } + } + print ("g", >> temp1) + fd = ""; delete (temp, verify=no) + + nspec = 1 + hselect (arcrefms, "naxis2", yes) | scan (nspec) + for (i = 1; i <= nspec; i+=1) + identify (arcrefms, section="line "//i, + database=database, coordlist="", nsum=1, + match=params.match, maxfeatures=50, zwidth=100., + ftype="emission", fwidth=params.fwidth, + cradius=params.cradius, threshold=params.threshold, + minsep=2., function=params.i_function, + order=params.i_order, sample="*", + niterate=params.i_niterate, + low_reject=params.i_low, high_reject=params.i_high, + grow=0., autowrite=yes, cursor=temp1, < temp2, + >G "dev$null", >>& temp) + delete (temp1, verify=no); delete (temp2, verify=no) + system.match ("Coordinate shift", temp, stop=no, print_file_n=yes, + metacharacte=yes) | tee (log1, > log2) + delete (temp, verify=no) + + dispcor (arcrefms, "", linearize=params.linearize, + database=database, table="", w1=INDEF, w2=INDEF, + dw=INDEF, nw=INDEF, log=params.log, flux=params.flux, + samedisp=yes, global=no, ignoreaps=no, confirm=no, + listonly=no, verbose=yes, logfile=logfile, > log2) + if (params.nsubaps > 1) { + imrename (arcrefms, temp, verbose=no) + scopy (temp, arcrefms, w1=INDEF, w2=INDEF, apertures="1-999", + bands="", beams="", apmodulus=0, offset=0, + format="multispec", clobber=no, merge=no, renumber=no, + verbose=no) + blkavg (temp, temp, 1, params.nsubaps, option="sum") + imcopy (temp, arcrefms//"[*,*]", verbose=no) + imdelete (temp, verify=no) + } + } + if (extn == ".ms") + print (arcref, >> done) + } +end |