I think it's finally done.

My contrast curve plotter/point source detector finally works as intended!

To detect a point source companion or plot a contrast curve requires the program to first determine the background noise in the image. The program does so by first determining the aperture of the star (radius of one FWHM) and then determining the pixels that make up annuli that are 2 FWHM wide around the star. The radius of each annulus is defined as the midpoint between the inner and outer radii that define the annulus. Each annulus increases in radius by FWHM/2. This overlap of annuli allow us to find point sources that might have been on the edge of an annulus and otherwise missed.

For the point source detector, we further break up each annulus into spheres arranged in a ring in the annulus. Each sphere has a radius of 1 FWHM and the spheres overlap each other by half a FWHM. We determine the background of the spherical regions and then find point sources according to the darkpix, thresh, and adj values I've mentioned in previous posts. By breaking up the annuli into regions, we can more accurately determine the local background noise in each region and see if any region contains a point source.


(To refresh: darkpix is a decimal that determines the percent of darkest pixels used to find the mean noise of the background. thresh is the threshhold, or the multiplier we use with sigma, the standard deviation of the flux. Only pixels with intensity values greater than thresh*sigma+background mean will be labeled as potential point sources. adj determines the number of adjacent pixels greater than thresh*sigma+background mean a potential point source has to have in order for a pixel to be labeled as the position of the point source. )


For the contrast curve graph, we measure the difference in magnitude in each annulus with the peak flux of the central star. We call this difference contrast, or delta M. The graph plots contrast with respect to the radius of each annulus (the midpoint between the inner and outer radii). The magnitude value in each annulus is calculated in this way:

After getting the mean background of the regions in each annulus, we average the values of each region to get the average intensity of the annulus itself. We then convert this to a delta mag number by using the equation:

Contrast = -2.5*log10(flux1/flux2)

Flux 1 is the average flux of the annulus; flux 2 is the peak flux of the star. Finally, we plot the contrast with respect to magnitude on our graph.

If there is no companion point source, the program will plot a contrast curve; otherwise, the program will label the location of the main point source and its companion(s). The user can opt to save either of these figures by pressing "w". A text file with some important information will be saved along with the figure.
-----
This serves as a general description of what my program does. There are some specific details that I also need to mention:

- The darkpix value in my program, when used to detect extra point sources, has been set as 0.8 to increase robustness of the point source detector. The darkpix value is generally set to 0.95 when determining which pixels constitute the "background" of the FITS file.
- The ratio of arcseconds to pixels changes for each fits file. However, due to the variety of FITS files I am working with (some of them without headers as well), I cannot create subclasses for each FITS file (as I originally intended) and automatically make my program look for the arcseconds to pixels ratio. Thus the user needs to manually type it in into the phot4.py module. Hopefully this process will become more user-friendly as I clean up the program a bit.
-When the user enters multiple fits files into the command line when working with IPhot, the user needs to make sure that all files are working correctly. If not, the program will exit as it reaches that file. I have decided as of now to not let my program simply "pass" over these files as it would let a silent error through.
-I will post more as I find more errors.

Matching Justin's Results

So it turns out that I actually wasn't able to match Justin's results because we forgot that Justin had sent me images that had intensities on a square-root scale (square-root because it makes fainter companions stick out). Thus, all I needed to do to correct for this would be to square the original FITS array before making some contrast curves.

So: here are the new results for two different contrast curves:

(Properties: threshold=3 sigma)

KIC109:

KIC118:

Compare to Justin's results: 
At radial distances of 0.25", 0.5", 1", 2", 4" respectively:
KIC109   3.7, 5.6, 7.8, 8.6, 8.6  
KIC118   4.2, 5.7, 7.2, 7.6, 7.6


Unfortunately, I can't usually go past 1" on my graphs yet due to the way my point source finder works (it'll crash my computer if I try). However you can see that this is generally sufficient distance to make the contrast curve flatten out by the end. 


You may also have noticed that my delta mags tend to be slightly lower than Justin's values. This could be due to the fact that Justin took a thin slice of the picture and calculated the delta mags based on distance in one direction, while I averaged the noise at a certain distance away in all directions. But overall, the graphs seem to agree with Justin's results pretty well (It seems accurate up to 1").


To-Do:
-I will attempt to work out some kinks in my contrast curve calculator to make it run a bit more efficiently.
-Will continue working on a way to find point sources, as fitting 2D-gaussians still doesn't seem to work.

Accepting Directories

Alright, now that I've finally adjusted to moving from Caltech to home (and learnt how to deal with my dust allergies), IPhot now accepts directories as well as individual .FITS files in one directory. (I guess I have to change the usage() description of IPhot as well.) Here are the three different ways IPhot can process and find .FITS files in command line:

1. By default, if no additional argument is written after 'IPhot', IPhot will assume the directory is the current directory. Any .FITS files in this directory will then be processed. Let's say you have a file called fitsimage.fits in this directory. Typing in IPhot will give you this:
   
   C:\[Your path here]> IPhot
   fitsimage.fits
   
   The image will open up and you'll be able to click on it and such just as before.
   
2. You can specify a directory after IPhot. IPhot will go through this directory and process all the .FITS images in that directory.
   
   C:\[Your path here]> IPhot "C:\\Example\\Path\\To\\Your\\Images\\" -d 0.95 -t 4
   ImageFromThisPath.fits
   
   ImageFromThisPath.fits and any other .fits images will then show up on the screen like before, and will create contrast curves in the way specified by the options following the path. You can still press 'q' to quite or any other key for the next image. Note the specific way you have to type the paths - it has to have quotes around it (since it's a string, and also because some path names have spaces inside them, like "Documents and Settings" that PC's use), and you should also write double \\'s to be safe.
   
3. This is the way I originally implemented this program to handle .FITS files.
    
   C:\[Your path here]> IPhot CK00012_snap_20110412.fits -d 0.97 -t 5
   CK00012_snap_20110412.fits
   
   The program will go through the .FITS files one by one and create contrast curves for them (if you make it do so) according to the options following the .FITS files. Note that you can still type in multiple individual .FITS file names. 

I gotta say, this was pretty fun to work on. But I'm still not done with my point sources finder -- it's rather reluctant to actually find some point sources (simulated or 'real' ones) -- so I still can't test this code very extensively to make sure that it works brilliantly. For now though, I'll say that it seems to work fine.

Contrast Curve Function 2.0 III and Point Source Detector

So it turns out that actually, because of various factors like dithering and such, my contrast curves are probably correct.
...
Yay! That's actually pretty exciting stuff.

However, now my point source detector is far too insensitive. I am trying to search for regions that have point sources, and then use the fitgaussian(...) function to detect those point sources. However, the fitgaussian(...) runs very slowly over large arrays and doesn't seem to be all that effective at pinpointing the point sources. Even after I've put down an artificial point source, it still doesn't detect anything in any of the regions. Perhaps it's the way I've implemented it, and I just have to keep working towards a solution.

Learned:
-What dithering means.
To-Do:
-Find a way to make my point source detector work correctly, whether by using the fit gaussian function, or some other function.

Contrast Curve Function v2.0 II

I finished the new contrastcurve(...) function today.


def contrastcurve(self,ann_space=0,darkpix=0.95,thresh=5): Returns list of radii and corresponding noise (thresh*sigma+mean background noise) and point sources (if any). It detects point sources by fitting gaussians over regions that have, on average, greater intensity than the noise and searches for peaks in those regions.


And it........... still doesn't work right. 


Here's an example though, of what the contrast curve looks like for KIC109:

Interestingly, however, if I calculate magnitude using log instead of log base 10, I get this result:

This agrees much more closely with Justin's data. Weird! I still haven't figured out why. 


By the way, Justin's data is this:
KIC109   3.7, 5.6, 7.8, 8.6, 8.6  at 0.25'', 0.5'', 1'', 2'', 4'' respectively.
(This is for thresh=3; mine are produced with thresh=5, but I ran the code again and this barely makes any difference). 


To-Do:
-Start testing how my code behaves around images with bright point sources (yes I said I would do this today, but then the contrast curves looked funny...) while I try to figure out why my contrast curves are behaving this way.
Learned:
-I actually learned something interesting today: If you have a function that returns multiple things in an array, it's much more efficient to do thing1, thing2, thing3=function() 
as opposed to thing1, thing2 = function()[0], function[1] even if you only want to save thing1 and thing2. This is probably because in the first case, Python will run the function once and save the corresponding variables, whereas in the second case, Python will run the function once, save the 0th element, and then run the function again and save the first. So I guess that ultimately what I learned is that if you have a function that takes a while to run, it's better to use the first case to save time, but if you have a function that's relatively fast for large amounts of data, but that returns a bunch of stuff, it's better to use the second case to save memory. If you have a slow function that returns lots of data.... well then good luck to you.

Contrast Curve Function v2.0 I

Today I finished most of the contrast curve function I talked about in my previous post. I haven't gotten it to return a list of detected point sources yet (it only returns radii and respective 5*sigma as of now), but it should integrate relatively well into all the other code I've written so far.

I have also been given a rough timeline of what I'm supposed to do for the rest of the summer (yes, although officially ten weeks, my project will continue till the end of summer/beginning of the fall term). So here it is:

Goals:
1) Rewrite my code and try and get the same results that Justin did.
2) Start using my code for all the fits files I have.
3) Detect problems and check them.
4) Repeat as necessary (This is not a fun step).
5) Make my code more user-friendly.
6) At some point start putting contrast curves up on a website.
.
.
.
n) Somehow help with the writing of a paper.

To-Do:
-Tomorrow I'll have to figure out how this new contrast curve function will actually detect point sources. Since we're talking about regions now, I might have to also change the way the program points out detected sources to the user.
-I also have to start writing code so that I can upload directories of .FITS files into the program. That way, the user won't have to spend a long time typing in the names of all the .FITS images.
-At some point I'll have to see what I can do about adapting my program to various arcsecond/pixel ratio. That is, I want my program to be able to recognize a type of image by its file name and automatically assign the correct arcsecond/pixel ratio.
-Finally, if possible/not too aesthetically ugly, I should try and make all the pictures axes display units in terms of arcseconds instead of pixels.

New Implementation

I talked with Tim and Justin about my code. The problem seems to be that in my current program, I am overestimating the amount of background noise since I'm calculating noise by averaging the flux from individual pixels. What I should do instead was to average noise from little regions and calculate the noise from the average of those regions. It seems like I'm just working with slightly bigger pixels, but by first averaging the flux from each region, I am able to smooth out my background and ultimately decrease the noise in each image.

So now, I have to code a different approach to finding point sources. It will still use much of what I have already written, but it's definitely the biggest change yet to my code.

Here's a rough outline of what I'm about to do:
I will partition each annulus into little circular regions of data. I will then average the flux in each region, and then find the standard deviation from the regions. This will be used to create sigma. The background mean will also be calculated from the regions. Detection of a point source apart from the main source will determined by if the a region has an average intensity greater than thresh*sigma. My program might pinpoint the specific location of the point source less precisely, but at least it won't detect fake point sources that turned out to just be speckles. The contrast curve will still be created in much the same way. It's just the way I determine background noise that will change the most significantly. 


I'll start working on this next week.