|
Oh the joys of having to rewrite other people's code so it works
|
View this Thread in Original format
| DJ Kenosis |
Written in IDL as a part of the data reduction package for the infrared instrument I used last month in Arizona. It wasn't doing what I was wanting it to do (nor what it should do), so I had to rewire it (just finished today) But still, grr....
---------
pro imageproc6, input_file_list, bad_pixel_mask, dark_image, flat_image, output_directory, NOGUI=nogui^M
;^M
; does dark and flat corrections of science images (already linearity corrected)^M
;^M
; respects NUMCOADDS for both science and dark images^M
;^M
; DPC 20051009 original code^M
;^M
;------------------------------------------------------------- read input lists -------------------------------------------^M
np = N_PARAMS()^M
if(N_ELEMENTS(NOGUI) eq 0) then NOGUI = 0^M
;^M
if(np eq 0) then begin^M
if(NOGUI eq 0) then begin^M
science_files = dialog_pickfile(Title="Select Science Images (*_LC.fits)",/MULTIPLE_FILES)^M
science_list = ' '^M
bad_file = dialog_pickfile(Title="Select Bad Pixel Mask")^M
dark_file = dialog_pickfile(Title="Select Dark Count Image")^M
flat_file = dialog_pickfile(Title="Select Flat Field Image")^M
output_dir = dialog_pickfile(Title="Select the Directory for the Output Science Images",/DIRECTORY)^M
endif else begin^M
science_files = strarr(4000)^M
science_list = ''^M
bad_file = ''^M
dark_file = ''^M
flat_file = ''^M
output_dir = ''^M
print, "Enter pathname to list of input Science images (*_LC.fits):"^M
read, science_list^M
print, "Enter pathname to the Bad Pixel Mask:"^M
read, bad_file^M
print, "Enter pathname to the Dark Count Image:"^M
read, dark_file^M
print, "Enter pathname to the Flat Field Image:"^M
read, flat_file^M
print, "Enter a pathname for the Directory to hold the output Science images:"^M
read, output_dir^M
endelse^M
endif else begin^M
science_list = input_file_list^M
science_files = strarr(4000)^M
bad_file = bad_pixel_mask^M
dark_file = dark_image^M
flat_file = flat_image^M
output_dir = output_directory^M
endelse^M
;^M
; if a science list was indicated, read the list^M
;^M
if(science_list ne ' ') then begin^M
openr, lu_in, science_list, /GET_LUN^M
onefile = ''^M
nfiles = 0L^M
while (~EOF(lu_in)) do begin^M
readf, lu_in, format='(a)', onefile^M
science_files[nfiles] = STRTRIM(onefile,2)^M
nfiles++^M
endwhile^M
FREE_LUN, lu_in^M
science_files = science_files[0:nfiles-1]^M
endif else begin^M
nfiles = size(science_files)^M
nfiles = nfiles[1]^M
endelse^M
;^M
;---------------------------------------------------------------------------------------------------------------^M
;^M
; "discovery" part of code -^M
;^M
; read dark, BPM, flat files^M
; capture dark Texp and NUMCOADDS^M
;^M
image_bad = readfits(bad_file, header_bad)^M
image_dark = readfits(dark_file, header_dark)^M
image_flat = readfits(flat_file, header_flat)^M
;^M
dark_Texp = FXPAR(header_dark, "EXPTIME")^M
dark_COADDING = FXPAR(header_dark, "COADDING")^M
dark_NUMCOADD = FXPAR(header_dark, "NUMCOADD")^M
;^M
; if dark COADDING true, generate a scaled dark image by NUMCOADD, using safe_math^M
;^M
if(dark_COADDING) then begin^M
image_dark_scaled = safe_math(image_dark, dark_NUMCOADD, image_bad, "/")^M
endif else begin^M
image_dark_scaled = image_dark^M
endelse^M
;^M
; temp, disable dark^M
;^M
;image_dark_scaled = image_dark_scaled * 0.0d0^M
;^M
;----------------------------------------------------------------------------------------------------------------^M
;^M
; Science Correction Loop^M
;^M
for ifile = 0, nfiles-1 do begin^M
;^M
; read science frame^M
;^M
image_science = readfits(science_files[ifile],header_science)^M
;^M
; read exposure time, COADDING, NUMCOADD from header^M
;^M
science_Texp = FXPAR(header_science, "EXPTIME")^M
science_COADDING = FXPAR(header_science, "COADDING")^M
science_NUMCOADD = FXPAR(header_science, "NUMCOADD")^M
;^M
; tests^M
;^M
; first, verify that both dark and science have identical exposure times^M
;^M
; if(science_Texp ne dark_Texp) then begin^M
; print, "ERROR - mismatch of Texp for dark and science images for frame ",FILE_BASENAME(science_files[ifile])^M
; print, " Correction of this file is aborted..."^M
; goto, next_science_image^M
; endif^M
;^M
; if science COADDING true, scale science image by NUMCOADD, using safe_math^M
;^M
if(science_COADDING) then begin^M
image_science_scaled = safe_math(image_science, science_NUMCOADD, image_bad, "/")^M
endif else begin^M
image_science_scaled = image_science^M
endelse^M
;^M
; perform the dark correction^M
image_science_dc = safe_math(image_science_scaled, image_dark_scaled, image_bad, '-')^M
;^M
; and the flat correction^M
;^M
image_science_dfc = safe_math(image_science_dc, image_flat, image_bad, '/')^M
;^M
; and if the science image was scaled, return its scaling^M
;^M
if(science_COADDING) then begin^M
image_science_final = safe_math(image_science_dfc, science_NUMCOADD, image_bad, '*')^M
endif else begin^M
image_science_final = image_science_dfc^M
endelse^M
; hot/cold pixel removal & smoothing TC 12/8/05
badmask2=fltarr(1024,1026)
badmask2(*,*)=0.
; lose_bad,image_science_final,image_science_final,1,4,bad
image_science_final=fixpix_rs(image_science_final,image_bad,iter=8)
lose_bad,image_science_final,image_science_final,1,1,bad
badcorr=where(image_science_final eq 1)
badnum=floor(median(image_science_final))
image_science_final[badcorr]=floor(median(image_science_final))
lose_bad,image_science_final,image_science_final,badnum,4
; index=where(image_science_final le max(image_science_final))
; dims=size(image_science_final,/dimensions)
; ncol=dims[0]
; col_index=index mod ncol
; row_index=index/ncol
badmask2(*,512:514)=1
image_science_final=fixpix_rs(image_science_final,badmask2,iter=5)
; lose_bad,image_science_final,image_science_final,1,0,bad
; image_science_final=fixpix_rs(image_science_final,bad,iter=2)
; badcorr=where(image_science_final eq 1)
; image_science_final=fixpix_rs(image_science_final,badcorr,iter=2)
;^M
; create output science filename^M
;^M
science_base_name = FILE_BASENAME(science_files[ifile])^M
; parts = STRSPLIT(science_base_name, "LC", /EXTRACT)^M
; new_name = parts[0] + "LDFC" + parts[1]^M
new_pathname = output_dir + science_base_name^M
;^M
; update science header^M
;^M
dcomment = STRING(format='(" Dark Corrected using image ",a)',FILE_BASENAME(dark_file))^M
fcomment = STRING(format='(" Flat Corrected using image ",a)',FILE_BASENAME(flat_file))^M
FXADDPAR, header_science, "COMMENT", dcomment^M
FXADDPAR, header_science, "COMMENT", fcomment^M
FXADDPAR, header_science, "FILENAME", science_base_name, "Fully Corrected Filename"^M
;^M
; and write output file^M
;^M
writefits, new_pathname, image_science_final, header_science^M
;^M
print, "Finished correcting and saving image ",science_base_name^M
next_science_image:^M
;^M
; do next science image^M
;^M
endfor^M
;^M
; all done, quit^M
;^M
print, "Science Dark and Flat Corrector is Done!"^M
;^M
end^M |
|
|
| philippe |
Hehe, sounds like your thesis is starting (with IDL) like mine did (with Fortran) :wtf:
Oh, and: sed "s/^M//g" ;)
BTW, do you remember this olive green T-shirt?
$ cd /pub
$ more beer
Paavo--such a geek--loved it!! |
|
|
| DJ Kenosis |
| quote: | Originally posted by philippe
Hehe, sounds like your thesis is starting (with IDL) like mine did (with Fortran) :wtf:
Oh, and: sed "s/^M//g" ;)
BTW, do you remember this olive green T-shirt?
$ cd /pub
$ more beer
Paavo--such a geek--loved it!! |
HAHA, yeah. I remember that.
You did your thesis work in cosmology, right? What were you working on?
So Philippe, the problem is that the data reduction package for MIMIR on the Perkins telescope at Lowell automatically sets the 'bad pixel values' to -1.e6. Of course that's not good if you're wanting to do photometry on a star that happens to be there and it also changes (slightly) the median value of the pixels. I had to feed the program a 'fake' bad pixel mask instead of the real one in order to partially work around this.
Complicating matters are the fact that MIMIR has a crack running right through the middle of it at around the 512th column. The count values here are not only bad but *variably bad*. So when you construct a sky frame out of all your frames and then subtract the sky frame from each of the images you have a row of pixels whose values are pretty much a crapshoot. I added a hot/cold pixel interpolation scheme to correct this, but the scheme was built for interpolating over individual bad pixels, not rows of them: it's not as good in this regard as the value for a given 'bad pixel' is always going to include another 'bad pixel' value since the subroutine only rejects the hottest and coldest neighboring pixels.
The routine here is the original image processing routine, altered by me to interpolate over bad pixels (image = object_frame - dark/(flat - dark)/median(flat-dark)). I've had to write another one for building a sky frame as well as sky subtraction. |
|
|
| philippe |
| Sorry, can't help you there :( I did numerical simulations and had nothing to do with data reduction :D |
|
|
| DJ Kenosis |
| quote: | Originally posted by philippe
Sorry, can't help you there :( I did numerical simulations and had nothing to do with observations :D |
haha, ok. Numerical simulations of...? |
|
|
| philippe |
| quote: | Originally posted by DJ Kenosis
haha, ok. Numerical simulations of...? |
The whole universe :p So, yes, cosmology :) |
|
|
| DJ Kenosis |
| quote: | Originally posted by philippe
The whole universe :p So, yes, cosomology :) |
Oh yah? I simulate a universe when I get bored with everything else :D :haha: |
|
|
| DJ Kenosis |
| quote: | Originally posted by Xtracktor
*takes notes* |
there's a test on friday. If you fail then 'no trance for you.'
- the trance nazi |
|
|
| Xtracktor |
| quote: | Originally posted by DJ Kenosis
there's a test on friday. If you fail then 'no trance for you.'
- the trance nazi |
Well, my last final is on wednesday...so just gonna have to add this to the schedule, at least this is WAY more interesting than poly sci :whip: |
|
|
| Electrophile |
Oh the joys of the von Richter reaction. Is the mechanism right? Is the mechanism wrong? Does it exist? I guess so. Let's just "accept" it.
 |
|
|
| dj_bas |
| Who wants to see my SQL code? |
|
|
|
|
