Last weekend I wrote 4 Python programs which concern raw FITS image sharpness:
All programs, and examples, are available in the Sharpness project. This is an XCode project, however XCode is only being used as an editor and file manager: there are no targets or executables in the project. The programs are run under Python 2.5.1 (either with 'python program.py' or from the '#!/usr/bin/env python' line in the file), and you can use BBEdit (etc...) as the program editor.
Usage: [python] sharpness.py directory [> output-file], or
[python] sharpness.py fits-file [min] [max]
To run on a folder of files, use:
> sharpness.py contk1 contk1/data1050.a.fits 2004-07-13 16:17:23.972260 0.04539629 contk1/data1047.a.fits 2004-07-13 16:17:16.547470 0.04383086 contk1/data1048.a.fits 2004-07-13 16:17:19.018101 0.04270498 contk1/data1049.a.fits 2004-07-13 16:17:21.525703 0.04204120 contk1/data1045.a.fits 2004-07-13 16:17:11.548302 0.03994618 contk1/data1042.a.fits 2004-07-13 16:17:03.969390 0.03979242 contk1/data1043.a.fits 2004-07-13 16:17:06.517683 0.03971004 contk1/data1044.a.fits 2004-07-13 16:17:08.984247 0.03949092 contk1/data1046.a.fits 2004-07-13 16:17:13.994263 0.03925742 contk1/data1041.a.fits 2004-07-13 16:17:01.523765 0.03888214
The program reads all .fits files in the folder, estimates their sharpness, and prints the sorted results. You can also redirect this output to a text file. The program takes less than 0.5 seconds to process each image, on a 2x2GHz MacBook.
To estimate the sharpness of a single image, use:
> sharpness.py contk1/data1041.a.fits Min radius Max radius Max wavelength Min wavelength Count Score ---------- ---------- -------------- -------------- ----- ----- 0.00000000 8.00000000 512.00000000 64.00000000 197 0.73594810 8.00000000 16.00000000 64.00000000 32.00000000 600 0.13945106 16.00000000 24.00000000 32.00000000 21.33333333 996 0.05120139 24.00000000 32.00000000 21.33333333 16.00000000 1416 0.02514732 32.00000000 40.00000000 16.00000000 12.80000000 1816 0.01377158 40.00000000 48.00000000 12.80000000 10.66666667 2188 0.00819796 48.00000000 56.00000000 10.66666667 9.14285714 2632 0.00515955 56.00000000 64.00000000 9.14285714 8.00000000 3008 0.00346567 64.00000000 72.00000000 8.00000000 7.11111111 3388 0.00242458 72.00000000 80.00000000 7.11111111 6.40000000 3840 0.00167302 80.00000000 88.00000000 6.40000000 5.81818182 4232 0.00115471 88.00000000 96.00000000 5.81818182 5.33333333 4604 0.00090070 96.00000000 104.00000000 5.33333333 4.92307692 5032 0.00072927 104.00000000 112.00000000 4.92307692 4.57142857 5432 0.00057161 112.00000000 120.00000000 4.57142857 4.26666667 5844 0.00044305 120.00000000 128.00000000 4.26666667 4.00000000 6208 0.00039044 128.00000000 136.00000000 4.00000000 3.76470588 6656 0.00033251 136.00000000 144.00000000 3.76470588 3.55555556 7012 0.00029439 144.00000000 152.00000000 3.55555556 3.36842105 7432 0.00025345 152.00000000 160.00000000 3.36842105 3.20000000 7848 0.00021591 160.00000000 168.00000000 3.20000000 3.04761905 8252 0.00019948 168.00000000 176.00000000 3.04761905 2.90909091 8640 0.00019078 176.00000000 184.00000000 2.90909091 2.78260870 9080 0.00017908 184.00000000 192.00000000 2.78260870 2.66666667 9428 0.00018725 192.00000000 200.00000000 2.66666667 2.56000000 9848 0.00018150 200.00000000 208.00000000 2.56000000 2.46153846 10248 0.00019118 208.00000000 216.00000000 2.46153846 2.37037037 10628 0.00019263 216.00000000 224.00000000 2.37037037 2.28571429 11080 0.00022758 224.00000000 232.00000000 2.28571429 2.20689655 11464 0.00029652 232.00000000 240.00000000 2.20689655 2.13333333 11868 0.00041801 240.00000000 248.00000000 2.13333333 2.06451613 12288 0.00082729 248.00000000 256.00000000 2.06451613 2.00000000 12654 0.00518241 512 - 64: 0.73594810 64 - 32: 0.13945106 32 - 16: 0.07634871 16 - 8: 0.03059476 8 - 4: 0.00828738 4 - 2: 0.00936998 sharpness = 0.03888214 contk1/data1041.a.fits 2004-07-13 16:17:01.523765 0.03888214
The program also prints a table of frequency/wavelength values, and the percent energy in the image at each wavelength. The sharpness is the energy in the 16 - 4 pixel band.
To get the sharpness of a single image, plus plots, use:
> sharpness.py contk1/data1041.a.fits 0 3000
where the last two numbers on the line are the min and max scale values of the image. Note that the plots shown are of the log spectrum and average log spectrum, which are not the actual data used to estimate the sharpness. See the sharpness.py program text for more details. Also note that you must close all the windows to exit the program.
Usage: [python] plotsharp.py sharpness-file
To run on a .txt file of sharpness results, use:
> plotsharp.py sharpness.txt 2400 entries Min time: contk/data1041.a.fits 2004-07-13 16:17:01.523765 0.03888214 Max time: contk/data3540.a.fits 2004-07-13 18:30:32.704249 0.01595912 Min sharp: contk/data3021.a.fits 2004-07-13 17:56:32.344494 0.00218952 Max sharp: contk/data1161.a.fits 2004-07-13 16:22:02.456419 0.04652096
Note that sharpness is estimated on the raw FITS image, full-frame. Therefore, any negative values, noise (readout, Poisson, or bad pixels), and the presence of the slit will affect this result. Also, bad images (usually meaning planet partially off-screen or multiple planet images) have not been removed from this set. Ideally, sharpness should be estimated after the image has been aligned and preprocessed, and should only be performed on a subset of the clouds which has been passed through a 2d Blackman filter. Once alignment and preprocessing have been performed, the sharpness can be recomputed in this manner.
Nevertheless, the results presented here correspond well with visual evaluation. Here are 6 images from best to worst, equally spaced in the sharpness.txt file output for 2400 7/13/07 images:
Usage: [python] plotfits.py fits-file [min] [max] [header-flag (anything)]
To run on a .fits file, use:
> plotfits.py contk/data1041.a.fits 0 3000 1
SIMPLE = T / DATA IS IN FITS FORMAT
BITPIX = 16 / 16 BITS TWOS COMPLEMENT INTEGERS
NAXIS = 2 / NUMBER OF AXIS
NAXIS1 = 512 / PIXELS ON 1st MOST VARYING AXIS
NAXIS2 = 512 / PIXELS ON 2nd MOST VARYING AXIS
DATAMIN = -4960 / MIN DATA VALUE IN FILE
DATAMAX = 20779 / MAX DATA VALUE IN FILE
DATAMEAN= 0.00 / MEAN DATA VALUE IN FILE
DIVISOR = 4 / Normalization value
ORIGIN = 'Institute for Astronomy'
TELESCOP= 'NASA IRTF'
INSTRUME= 'SpeX, IRTF Guider/Camera'
OBSERVER= 'Eliot Young and Mark Bullock'
OBJECT = 'Venus'
COMMENT = 'Morning Venus'
IRAFNAME= 'data1041.a.fits'
DSPTMFLE= 'none' / DSPTimingInfo File
BEAM = 'A' / Object(A) or sky(B)
TIME_OBS= '16:17:01.523765' / UT TIME OF ACQISTION ('hh:mm:ss.ss')
DATE_OBS= '2004-07-13' / UT DATE OF ACQUISITION ('yyyy/mm/dd')
TIME_GPS= '17:01.594634' / UT TIMESTAMP from GPS ('mm:ss.usec')
CAMMODE = 0 / CameraMode is Basic
ITIME = 0.2500 / INTEGRATION TIME IN SECONDS
CO_ADDS = 4 / NUMBER OF INTEGRATIONS
CYCLES = 400 / Number of cycles
NUMARRAY= 1 / Number of data sub-arrays
ARRAY0 = 0,0,512,512 / x,y,wid,hgt of data sub-arrays
BEAMPAT = 0 / Beam Pattern
BBMODE = 2 / BB Mode is BBMODE_D16
CBMODE = 1 / CB Mode is ARC_D
NDR = 1 / Number of Non-Destructive Reads
SLOWCNT = 12 / Number of Slow Counts
GRSTCNT = 800 / Global Reset Count. 1 cnt = 25 nsec
GPSTMFLG= 1 / GPS_Time_Flag
BGRESET = 1,2000,1600 / Backgroud Reset's flag, msec, cnt
RMEX206 = 0 / RemoveExtra206MuxData. MUX206=0
TABLE_MS= 235.508 / msec to clock out array
CALMIR = Out / menu=0,pos=252000
DIT = open / menu=2,pos=160000
OSF = Open / menu=0,pos=21333
ROTATOR = 90.00 / Rotator_Mechanical_Angle; pos=180000
POSANGLE= 0.00 / Sky_Position_Angle
SLIT = 0.3x15 / menu=2,pos=89000
GRAT = LowRes15 / menu=5,pos=179350
GFLT = contK / menu=10,pos=448000
AFOC = 50000 / pos=50000,volts= 0.24
PLATE_SC= 0.12 / PlateScale arcsec/pixel
FLTZP = 0.0 / GFlt Filter Zero Point
EXT_COF = 0.000 / GFlt Extinction coefficient
TC330A = '30.00,74.13,30.00,Med30.00%' / Spectragraph ChA,ChB,SetPt,Heater
TC330B = '30.00,11.89,30.00,Med10.00%' / Guider ChA,ChB,SetPt,Heater
TC208 = ' 49.43 15.67 74.30 73.81 74.56 77.93 94.25 277.70'
QTH_LAMP= Off / Quart Lamp
INC_LAMP= Off / Incandesce Lamp
IR_SRC = Off / IR Source
ARG_SRC = Off / Argon Lamp
SHUTTER = Close / Cal Box shutter
TPD = ' 04:48:10.56 17:39:54.3 -03:25:29.47 1.511 2000.0 -O'
RA = 04:48:10.56 / Right Ascension
DEC = 17:39:54.3 / Declination
HA = -03:25:29.47 / Hour Angle
AIRMASS = 1.511 / Air Mass
EPOCH = 2000.0 / Epoch
To plot an image without showing the header, use:
> plotfits.py contk/data1041.a.fits 0 3000
To plot an image using auto-scale, omit the min and max values, or pass 0 for both:
> plotfits.py contk/data1041.a.fits
Note that you can display several images by using multiple plotfits.py commands, each terminated with the '&' character:
20: ~/Projects/Python/Sharpness > plotfits.py contk/data1041.a.fits & [1] 283 21: ~/Projects/Python/Sharpness > plotfits.py contk/data1042.a.fits & [2] 284 22: ~/Projects/Python/Sharpness > plotfits.py contk/data1043.a.fits & [3] 285
Usage: [python] mergesharp.py sharpness-file1 sharpness-file2 [> output-file]
You can combine 2 sharpness .txt files with:
27: ~/Projects/Python/Sharpness > cat sharpness1.txt contk/data1050.a.fits 2004-07-13 16:17:23.972260 0.04539629 contk/data1047.a.fits 2004-07-13 16:17:16.547470 0.04383086 contk/data1048.a.fits 2004-07-13 16:17:19.018101 0.04270498 contk/data1049.a.fits 2004-07-13 16:17:21.525703 0.04204120 contk/data1045.a.fits 2004-07-13 16:17:11.548302 0.03994618 contk/data1042.a.fits 2004-07-13 16:17:03.969390 0.03979242 contk/data1043.a.fits 2004-07-13 16:17:06.517683 0.03971004 contk/data1044.a.fits 2004-07-13 16:17:08.984247 0.03949092 contk/data1046.a.fits 2004-07-13 16:17:13.994263 0.03925742 contk/data1041.a.fits 2004-07-13 16:17:01.523765 0.03888214 28: ~/Projects/Python/Sharpness > cat sharpness2.txt contk/data1051.a.fits 2004-07-13 16:17:26.529575 0.04340547 contk/data1053.a.fits 2004-07-13 16:17:31.528383 0.04169112 contk/data1052.a.fits 2004-07-13 16:17:28.970800 0.04164434 contk/data1058.a.fits 2004-07-13 16:17:43.974738 0.04061172 contk/data1060.a.fits 2004-07-13 16:17:49.016979 0.03716876 contk/data1059.a.fits 2004-07-13 16:17:46.548851 0.03468087 contk/data1057.a.fits 2004-07-13 16:17:41.512497 0.03426781 contk/data1055.a.fits 2004-07-13 16:17:36.507754 0.03303268 contk/data1056.a.fits 2004-07-13 16:17:38.970658 0.03005856 contk/data1054.a.fits 2004-07-13 16:17:33.970351 0.02946241 29: ~/Projects/Python/Sharpness > mergesharp.py sharpness1.txt sharpness2.txt contk/data1050.a.fits 2004-07-13 16:17:23.972260 0.04539629 contk/data1047.a.fits 2004-07-13 16:17:16.547470 0.04383086 contk/data1051.a.fits 2004-07-13 16:17:26.529575 0.04340547 contk/data1048.a.fits 2004-07-13 16:17:19.018101 0.04270498 contk/data1049.a.fits 2004-07-13 16:17:21.525703 0.04204120 contk/data1053.a.fits 2004-07-13 16:17:31.528383 0.04169112 contk/data1052.a.fits 2004-07-13 16:17:28.970800 0.04164434 contk/data1058.a.fits 2004-07-13 16:17:43.974738 0.04061172 contk/data1045.a.fits 2004-07-13 16:17:11.548302 0.03994618 contk/data1042.a.fits 2004-07-13 16:17:03.969390 0.03979242 contk/data1043.a.fits 2004-07-13 16:17:06.517683 0.03971004 contk/data1044.a.fits 2004-07-13 16:17:08.984247 0.03949092 contk/data1046.a.fits 2004-07-13 16:17:13.994263 0.03925742 contk/data1041.a.fits 2004-07-13 16:17:01.523765 0.03888214 contk/data1060.a.fits 2004-07-13 16:17:49.016979 0.03716876 contk/data1059.a.fits 2004-07-13 16:17:46.548851 0.03468087 contk/data1057.a.fits 2004-07-13 16:17:41.512497 0.03426781 contk/data1055.a.fits 2004-07-13 16:17:36.507754 0.03303268 contk/data1056.a.fits 2004-07-13 16:17:38.970658 0.03005856 contk/data1054.a.fits 2004-07-13 16:17:33.970351 0.02946241
This program reads both files and resorts the combined results. You can also redirect this output to another .txt file. Note that although the sharpness.py program is fairly fast to run, you can get almost double the throughput on a dual-core machine by running two copies of sharpness.py, one in each of two terminal windows. You can then use mergesharp.py to combine the results.
İSky Coyote 2007