diff options
Diffstat (limited to 'noao/imred/crutil/src/crnebula.cl')
| -rw-r--r-- | noao/imred/crutil/src/crnebula.cl | 116 |
1 files changed, 116 insertions, 0 deletions
diff --git a/noao/imred/crutil/src/crnebula.cl b/noao/imred/crutil/src/crnebula.cl new file mode 100644 index 00000000..c5f214f9 --- /dev/null +++ b/noao/imred/crutil/src/crnebula.cl @@ -0,0 +1,116 @@ +# CRNEBULA -- Cosmic ray cleaning for images with fine nebular structure. + +procedure crnebula (input, output) + +file input {prompt="Input image"} +file output {prompt="Output image"} +file crmask {prompt="Cosmic ray mask"} +file residual {prompt="Residual median image"} +file rmedresid {prompt="Residual ring median image"} +real var0 = 1. {prompt="Variance coefficient for DN^0 term", min=0.} +real var1 = 0. {prompt="Variance coefficient for DN^1 term", min=0.} +real var2 = 0. {prompt="Variance coefficient for DN^2 term", min=0.} +real sigmed = 3 {prompt="Sigma clip factor for residual median"} +real sigrmed = 3 {prompt="Sigma clip factor for residual ring median"} + +int mbox = 5 {prompt="Median box size"} +real rin = 1.5 {prompt="Inner radius for ring median"} +real rout = 6 {prompt="Outer radius for ring median"} +bool verbose = no {prompt="Verbose"} + +begin + file in, out, med, sig, resid, rmed + struct expr + + # Query once for query parameters. + in = input + out = output + + # Check on output images. + if (out != "" && imaccess (out)) + error (1, "Output image already exists ("//out//")") + if (crmask != "" && imaccess (crmask)) + error (1, "Output mask already exists ("//crmask//")") + if (residual != "" && imaccess (residual)) + error (1, "Output residual image already exists ("//residual//")") + if (rmedresid != "" && imaccess (rmedresid)) + error (1, + "Output ring median difference already exists ("//rmedresid//")") + + # Create median results. + med = mktemp ("cr") + sig = mktemp ("cr") + resid = residual + if (resid == "") + resid = mktemp ("cr") + if (verbose) + printf ("Creating CRMEDIAN results\n") + crmedian (in, "", crmask="", median=med, sigma=sig, residual=resid, + var0=var0, var1=var1, var2=var2, lsigma=100., hsigma=sigmed, + ncmed=mbox, nlmed=mbox, ncsig=25, nlsig=25) + + # Create ring median filtered image. + rmed = mktemp ("cr") + rmedian (in, rmed, rin, rout, ratio=1., theta=0., zloreject=INDEF, + zhireject=INDEF, boundary="wrap", constant=0., verbose=verbose) + + # Create output images. + if (rmedresid != "") { + printf ("(a-b)/c\n") | scan (expr) + imexpr (expr, rmedresid, med, rmed, sig, dims="auto", + intype="auto", outtype="real", refim="auto", bwidth=0, + btype="nearest", bpixval=0., rangecheck=yes, verbose=no, + exprdb="none") + imdelete (rmed, verify-) + imdelete (sig, verify-) + if (out != "") { + if (verbose) + printf ("Create output image %s\n", out) + printf ("((a<%.3g)||(abs(b)>%.3g)) ? c : d\n", + sigmed, sigrmed) | scan (expr) + imexpr (expr, out, resid, rmedresid, in, med, dims="auto", + intype="auto", outtype="auto", refim="auto", bwidth=0, + btype="nearest", bpixval=0., rangecheck=yes, verbose=no, + exprdb="none") + } + if (crmask != "") { + if (verbose) + printf ("Create cosmic ray mask %s\n", crmask) + printf ("((a<%.3g)||(abs(b)>%.3g)) ? 0 : 1\n", + sigmed, sigrmed) | scan (expr) + set imtype = "pl" + imexpr (expr, crmask, resid, rmedresid, dims="auto", + intype="auto", outtype="short", refim="auto", bwidth=0, + btype="nearest", bpixval=0., rangecheck=yes, verbose=no, + exprdb="none") + } + imdelete (med, verify-) + } else { + if (out != "") { + if (verbose) + printf ("Create output image %s\n", out) + printf ("((a<%.3g)||(abs((b-c)/d)>%.3g)) ? e : b\n", + sigmed, sigrmed) | scan (expr) + imexpr (expr, out, resid, med, rmed, sig, in, dims="auto", + intype="auto", outtype="auto", refim="auto", bwidth=0, + btype="nearest", bpixval=0., rangecheck=yes, verbose=no, + exprdb="none") + } + if (crmask != "") { + if (verbose) + printf ("Create cosmic ray mask %s\n", crmask) + printf ("((a<%.3g)||(abs((b-c)/d)>%.3g)) ? 0 : 1\n", + sigmed, sigrmed) | scan (expr) + set imtype = "pl" + imexpr (expr, crmask, resid, med, rmed, sig, dims="auto", + intype="auto", outtype="short", refim="auto", bwidth=0, + btype="nearest", bpixval=0., rangecheck=yes, verbose=no, + exprdb="none") + } + imdelete (med, verify-) + imdelete (sig, verify-) + imdelete (rmed, verify-) + } + if (residual == "") + imdelete (resid, verify-) +end |