Sunday, July 10, 2011

Automated Astrophotography with Python - Part 7

FOCUSING WITH FOCUSMAX
Achieving great focus is undoubtedly one of the most critical operations for any astronomical imaging application. Using the free application, FocusMax, this step can be reduced to a simple, repeatable, and highly effective operation. Fortunately, FocusMax provides all of the necessary interfaces to quickly automate this important step. The following listing shows a new class named 'cFocuser'. This class is responsible for controlling all aspects of the focus operation. The listing for the 'cFocuser' class is shown below:
import time
import pythoncom
import win32com.client

FOCUSMAXTIME = 150
LASTFOCHFD = -1

ERROR   = True
NOERROR = False

##------------------------------------------------------------------------------
## Class: cFocuser
##------------------------------------------------------------------------------
class cFocuser:
    def __init__(self,printlog):
        printlog.log(0,"Connecting to FocusMax...")
        self.__FM = win32com.client.Dispatch("FocusMax.FocusControl")

    def checkFocusStar(self,printlog):
        printlog.log(0,"Checking focus by measuring HFD...")
        HFD = []
        count = 0
        # take 7 HFD measurements, sort the list, report the median value
        while count < 7:
            busy = self.__FM.SingleExposeAsyncStatus        
            self.__FM.SingleExposeAsync()
            # wait for measurement to end
            startTime = time.time()
            busy = -1
            while busy == -1 and (time.time() - startTime) < FOCUSMAXTIME:
                time.sleep(1)
                busy = self.__FM.SingleExposeAsyncStatus
            if self.__FM.SingleExposeAsyncStatus != 1:
                printlog.log(0,"HFD Measurement failed or timed out")
                return -1        
            result = self.__FM.HalfFluxDiameter
            flux = self.__FM.TotalFlux
            time.sleep(0.5)
            printlog.log(1,"HFD Measurement: %0.2f  Total Flux: %d" %
                         (result,flux))
            if flux > 20000:
                HFD.append(result)
                count = count + 1
        HFD.sort()
        j = len(HFD)
        if not j%2:
            median = (HFD[(j/2)-1]+HFD[j/2])/2.0
        else:
            median = HFD[j/2]
        focuserPos = self.__FM.Position
        printlog.log(0,"Median HFD = %0.2f at Pos = %d" % (median,focuserPos))
        return median

    def focusStar(self,printlog):
        global LASTFOCHFD
        printlog.log(0,"Starting FocusMax Focus method...")
        try:
            busy = self.__FM.FocusAsyncStatus        
            start = self.__FM.FocusAsync()
        except pythoncom.com_error, (hr, msg, exc, arg):
            printlog.log(0,"ERROR: %s" % exc[2])
            return ERROR
        else:
            if start:
                # wait for focus method to end
                startTime = time.time()
                busy = -1
                while busy == -1 and (time.time() - startTime) < FOCUSMAXTIME:
                    time.sleep(1)
                    busy = self.__FM.FocusAsyncStatus
                if self.__FM.FocusAsyncStatus != 1:
                    printlog.log(0,"Focus method failed or timed out.")
                    return ERROR        
                printlog.log(0,"FocusMax Focus method complete...")
                tgtHFD = self.__FM.HalfFluxDiameter
                if tgtHFD > 0.5:
                    median = self.checkFocusStar(printlog)
                    if median == -1:
                        return ERROR
                    elif median < 0.5 or median > 4.0:
                        printlog.log(0,"ERROR: Focus results out of range.")
                        return ERROR
                    else:
                        LASTFOCHFD = median
                        return NOERROR
                else:
                    printlog.log(0,"ERROR: Focus method failed.")
                    return ERROR
            else:
                printlog.log(0,"ERROR: FocusMax Focus method failed to start.")
                return ERROR

##
## END OF 'cFocuser' Class
##
The constructor for 'cFocuser' is called whenever a new object of the class is instantiated and is used to create a new object (__FM) that is bound to the FocusMax object (FocusControl) required by this class. Following the constructor, the 'checkFocusStar()' method performs a check of focus quality by invoking the 'SingleExposeAsync()' method. This method is called multiple times to take a series of measurement of a suitable focus star that has previously been centered in the field of view. After waiting for the exposure to complete, the method checks the properties that holds the measurement's HFD and total flux (brightness). The measurements are written to the screen and to the log file and the HFD value is appended to a list data structure. After seven measurements are collected, the HFD list is sorted and the median value is extracted and returned to the method's calling routine.

The other method in this class is named 'focusStar()' and performs the actual focus operation (nearly equivalent to clicking on the 'Focus' button on FocusMax's front panel). This method invokes the FocusMax 'FocusAsync()' method to kick off the focus operation. The method then polls the 'FocusAsyncStatus' property to determine when the focus operation is complete. If the focus star's HFD is a reasonable value, the 'checkFocusStar()' method is called to determine the median of seven consecutive HFD measurements. The median value is then assigned to a global variable, 'LASTFOCHFD'.

UNIT TESTING
Unlike unit tests for modules presented up to this point, the unit test for this class can not be simulated. Therefore, actual stars must be available to adequately verify the operation of this module. Due to poor weather and equipment failures I have been unable to generate and verify a unit test for the 'cFocuser' class. As soon as I am able, I will update this post with a unit test.

WHAT'S NEXT?
In the next post, I will demonstrate the method that I use to automatically select focus stars. This method will work with previous modules and result in a means to select a suitable focus star near the object to image, slew to the star, focus and determine median HFD, then slew to the main object to begin imaging. This procedure, along with periodic slews back to the focus star to check and ,if necessary, re-focus, is the main sequence of events that I use for my astronomical imaging script.

Tuesday, June 28, 2011

M21 - Open Cluster

(Click image to go to info page)

The open cluster M21 is a loose collection of about 50 stars located less than one degree northeast of the Trifid nebula (M20) in Sagittarius. The cluster lies approximately 4000 light years away and has a diameter of about 15 arc-minutes. Together, M20 and M21 make up an asterism known as "Webb's Cross" with M21 at the northern end and M20 at the southern end of the cross. M21 was first observed by Charles Messier in 1764. [Source: O'Meara, Stephen James. "The Messier Objects". Cambridge: University of Cambridge and Sky Publishing, 1998. Print.]

The data for this LRGB photo was captured between 3 June and 5 June 2011. All images were taken at the Chiefland Astronomy Village in Chiefland, Florida. Complete details regarding this image can be found by clicking on the image above or by clicking here.

Monday, June 27, 2011

M20 - Trifid Nebula

(Click image to go to info page)

The Trifid Nebula is a large, bright nebula and open cluster located in the constellation Sagittarius. It consists of four bright lobes separated by dark lanes of obscuring matter. As seen in the photo above and in the full-frame photo found here, M20 is located in a densely populated region of the Milky Way. At a magnitude of 5.9 and a distance of approximately 5000 light-years, the Trifid was first observed by Charles Messier in 1764. [Source: O'Meara, Stephen James. "The Messier Objects". Cambridge: University of Cambridge and Sky Publishing, 1998. Print.]

The data for this LRGB photo was captured between 3 June and 5 June 2011. All images were taken at the Chiefland Astronomy Village in Chiefland, Florida. Complete details regarding this image can be found by clicking on the image above or by clicking here.

Sunday, June 26, 2011

Automated Astrophotography with Python - Part 6

PLATE SOLVING WITH CCDSOFT AND THE SKY6
One of the most important features of any automated imaging script includes the ability to determine exactly where the imaging system is pointing at any given time. This is accomplished by astrometric plate solves that correlates positions of stars in an image with predetermined, known star positions. The combination of CCDSoft and the Sky6 (or MaxIm DL Pro and Pinpoint) provide a software path to implement this feature. For my astro-imaging purposes, I've chosen a more convoluted path to achieving plate solves. Specifically, I acquire my plate solve images using MaxIm DL (which pulls approximate RA and Dec coordinates of the image's center from the Sky6 via an ASCOM driver), then load the image into CCDSoft to perform the actual plate solve operation.

In order to capture the initial image to be plate solved it is necessary to add some code to the 'cCamera' class. A new class constant that specifies the save directory of the plate solve images is added at the top of the module:
PLATESOLVE_PATH = r"c:\\imaging scripts\\plate_solve_images\\"
Also, two new methods are added to the 'cCamera' class. These new methods are shown below:
def exposePlateSolve(self,length,filterSlot,name,printlog):
    self.__CAMERA.Expose(length,1,filterSlot)
    while not self.__CAMERA.ImageReady:
        time.sleep(0.5)
    # save exposure into plate solve folder
    printlog.log(0,"Saving plate solve image...")
    self.__CAMERA.SaveImage(PLATESOLVE_PATH + name + "_platesolve.fit")

def setQuarterSubFrame(self):
    self.__CAMERA.StartX = self.CAMERA.CameraXSize / 4
    self.__CAMERA.StartY = self.CAMERA.CameraYSize / 4
    self.__CAMERA.NumX = self.CAMERA.CameraXSize / 2
    self.__CAMERA.NumY = self.CAMERA.CameraYSize / 2
The method 'exposePlateSolve()' is used to capture the image, generate a filename, and save the image to the directory specified by PLATESOLVE_PATH. (Typically, all my plate solve images are unbinned, 10-second exposures taken through the luminance filter.) Since plate solves using CCDSoft and the Sky6 do not execute reliably on large images, I always use sub-framing to reduce the image size to one-fourth the area of a full-size frame. The 'setQuarterSubFrame()' method is responsible for setting up the reduced sub-frame so that the center of the sub-frame and the center of the full-size frame coincide.

The next listing shows a new class named 'cPlateSolver'. This class is responsible for initiating the plate solve operation through CCDSoft. CCDSoft, in turn, calls on the Sky6 to assist with the matching of stars in the database with the stars in the plate solve image. The listing for 'cPlateSolver' is shown below:
import pythoncom
import win32com.client

IMAGESCALE = 1.808
PLATESOLVE_PATH = r"c:\\Imaging_Scripts\\plate_solve_images\\"

class cPlateSolver:
    def __init__(self,printlog):
        printlog.log(0,"Connecting to CCDSoft...")
        self.__IMAGE = win32com.client.Dispatch("CCDSoft.Image")

    def insertWCS(self,name,targetRA,targetDec,mount,printlog):
        plateSolveResults = []
        # open platesolve.fit in CCDSoft and perform plate solve
        self.__IMAGE.Path = PLATESOLVE_PATH + name + "_platesolve.fit"
        self.__IMAGE.Open()
        # set image scale in CCDSoft
        self.__IMAGE.ScaleInArcsecondsPerPixel = IMAGESCALE
        try:
            self.__IMAGE.InsertWCS()
        except pythoncom.com_error, (hr, msg, exc, arg):
            errorstr = exc[2]
            errorstr = errorstr.replace('\n',' ')
            printlog.log(0,"Plate Solve %s" % errorstr)
            self.__IMAGE.Close()
            plateSolveResults = [0.0,0.0,-1,-1]
        else:
            # if no errors, report results
            printlog.log(0,"Plate Solve Results:")
            printlog.log(1,"Image Scale: %0.3f asp" % IMAGESCALE)
            printlog.log(1,"North Angle: %0.2f deg" % self.__IMAGE.NorthAngle)
            # get 2000 Epoch RA and Dec of center of plate solved image
            RADec = self.__IMAGE.XYToRADec(self.__IMAGE.Width/2.0,
                                         self.__IMAGE.Height/2.0)
            plateSolveResults.append(RADec[0])
            plateSolveResults.append(RADec[1])
            J2000RA = mount.getAngleToDMS(RADec[0])
            printlog.log(1,"Solved RA (J2000) : %02dh %02dm %0.2fs" % \
                         (J2000RA[0],J2000RA[1],J2000RA[2]))
            J2000DEC = mount.getAngleToDMS(RADec[1])
            printlog.log(1,"Solved DEC (J2000): %02d deg %02d' %0.2f\"" % \
                         (J2000DEC[0],J2000DEC[1],J2000DEC[2]))
            # Compute angular separation and position angle from
            # plate solved location to target location
            angSeparation = mount.getAngularSeparation(targetRA,targetDec,
                                                       RADec[0],RADec[1])
            positionAngle = mount.getPositionAngle(targetRA,targetDec,
                                                   RADec[0],RADec[1])
            angSeparation = angSeparation * 3600
            plateSolveResults.append(angSeparation)
            plateSolveResults.append(positionAngle)
            printlog.log(1,"Pointing error = %0.2f arcsec at %0.2f deg" % \
                         (angSeparation,positionAngle))
            self.__IMAGE.Close()
        return plateSolveResults        

    def execPlateSolve(self,camera,mount,name,targetRA,targetDec,printlog,
                       mode="Image"):
        # set binning to unbinned
        printlog.log(0,"Setting imager bin mode to 1x1")
        camera.setBinning(1,printlog)
        # set up camera for quarter-size subframe in center
        camera.setQuarterSubFrame()
        # take 10 second exposure
        printlog.log(0,"Exposing 10 second image for Plate Solve...")
        if mode.upper() == "IMAGE":
            camera.exposePlateSolve(10.0,0,PLATESOLVE_PATH+name,printlog)
        else:
            camera.exposePlateSolve(10.0,0,PLATESOLVE_PATH+name+'_test',
                                    printlog)
        #insert WCS
        printlog.log(0,"Performing plate solve with CCDSoft and The Sky6...")
        return self.insertWCS(name,targetRA,targetDec,mount,printlog)
            
##
## END OF 'cPlateSolver' Class
##
Two constants are declared and initialized at the top of the listing. 'IMAGESCALE' contains the arc-seconds per pixel resolution of the image to be solved and 'PLATESOLVE_PATH' is the path to the directory where the image has been stored. The class constructor creates an instance of the 'Image' object of CCDSoft.

The plate solve operation is initiated by a call to the 'execPlateSolve()' method of 'cPlateSolver'. This method sets up the binning and sub-framing of the plate solve image then calls the 'exposePlateSolveImage()' method of the 'camera' instance of the 'cCamera' object. The image is then exposed and saved under a name that is generated using the 'name' argument of the method (the save filename is modified if the mode is determined to be something other than 'Image'). After the plate solve image is exposed and saved, the plate solve operation is performed by a call to the 'insertWCS()' method.

The 'insertWCS()' method is fairly simple. After opening the plate solve image in CCDSoft and setting the 'ScaleInArcsecondsPerPixel' attribute, the 'insertWCS()' method of CCDSoft's 'Image' object is called from inside a try...except...else construct. If the plate solve fails, an informative error message is printed to the screen and log file. Otherwise, the RA and Dec coordinates of the plate solve image's center is determined using the 'XYToRADec()' method of CCDSoft's 'Image' object. The RA and Declination are reformatted, written to screen and log file, and appended to the 'plateSolveResults' list data structure. Additionally, the separation (in arc-seconds) and position angle from the solved coordinates to the intended coordinates are found, written to screen and log file, and also appended to the 'plateSolveResults' list. This is is then returned to the calling routine ('execPlateSolve()') which, in turn, returns the list to its calling routine.

The unit test for this class demonstrates use of the methods in the 'cPlateSolver' class and is shown below:
if __name__ == '__main__':

    import cLogger_6
    import cCamera_6
    import cMount_6

    # set up path and name of the logfile
    cCameraLogFile = r"C:\Imaging_Scripts\logfiles\cplatesolver.log"

    # instantiate the class with the path and name
    printer = cLogger_6.cLogger(cCameraLogFile)
    printer.logHeader('Unit Test cPlateSolver Class')

    # Create an instance of the cCamera class
    testCamera = cCamera_6.cCamera(printer)

    # Create an instance of the cMount class
    testMount = cMount_6.cMount(printer)

    # Create an instance of the cPlateSolver class
    testPlateSolve = cPlateSolver(printer)

    objects = ['M92','M107']
    for obj in objects:
        printer.log(0,"")
        if not testMount.findObject(obj):
            printer.log(0,"%s located" % (obj,))
            coords = testMount.getCoordinates()
            targetRA = coords['RA J2000']
            targetDec = coords['DEC J2000']
            DMScoords = testMount.getAngleToDMS(targetRA)
            printer.log(0,"%s RA  (J2000) = %02.0fh %02.0fm %02.0fs" % \
                        (obj,DMScoords[0],DMScoords[1],DMScoords[2]))
            DMScoords = testMount.getAngleToDMS(targetDec)
            printer.log(0,"%s Dec (J2000) = %02.0fd %02.0f' %02.0f\"" % \
                        (obj,DMScoords[0],DMScoords[1],DMScoords[2]))
            psResults = testPlateSolve.execPlateSolve(testCamera,testMount,
                                                      obj,targetRA,targetDec,
                                                      printer,'test')
        else:
            printer.log(0,"%s could not be found." % obj)

    printer.logFooter()
    printer.blank()
In the unit test shown above, after importing the 'cLogger', 'cCamera', and 'cMount' modules, an instance of the class contained in each module is created. An instance of the 'cPlateSolver' class is then created and two object names are placed in a list. Next, a loop is entered that iterates through the object list. The first step of the loop is to find the object in the Sky6's database and, if found, extract the object's RA and Declination coordinates which are used as the target or intended coordinates. The 'execPlateSolve()' method of the 'cPlateSolver' object is then called with mode set to something other than 'Image'. This causes the plate solve image's name to be mangled by inserting '_test' in the middle of the file name. This is so as to not overwrite the actual plate solve image that is used to generate a successful plate solve after image acquisition. The process is then repeated for each object in the list.

TESTING THE CLASSES
Like in previous posts, before testing classes in this section, the Sky6 simulator must be manually connected to the telescope control system using the usual procedure. Also, the imaging and guide cameras in MaxIm DL should be first connected to the simulated cameras. The 'cPlateSolver' module can then be executed from the IDE in the normal manner. The source files described in this chapter can be downloaded from chapter_6.zip. Note: the .fit files included in the .zip must be placed in the same directory that the plate solve images are saved to.

As an additional exercise, edit the unit test for the 'cPlateSolver' module and completely remove the 'test' argument in the call to the 'execPlateSolve()' method then run the module. This will cause the plate solve image's name to NOT be mangled and will thus overwrite the real plate solve image included in the .zip file (so, you should make a copy of the actual .fit image prior to this test.). After the plate solve image is taken, the image generated from the simulator will be used for a plate solve. Obviously, this will fail but the manner of the failure is instructive. If the error message returned is 'Plate Solve Error: command failed. Error code = 206 (0xce).' this most likely means the RA and Declination coordinates of the center of the image were not placed into the FITS header of the plate solve image (look for keywords 'OBJCTRA' and 'OBJCTDEC'). Conversely, if the error message is 'Plate Solve Image link failed because there are not enough stars in the image. Try adjusting the image's background and range. Error code = 651 (0x28b).', this indicates that the image just can't be solved (and this is the expected behavior).

WHAT'S NEXT?
In the next post, I will present a new class that demonstrates the use of 'FocusMax' using MaxIm DL Pro. The module will include a unit test which can't be executed from in a simulator mode but must have real stars, real images, and a properly configured 'FocusMax' to execute successfully.

CLICK HERE FOR NEXT POST IN THIS SERIES

Tuesday, June 21, 2011

M5 - Globular Cluster Now In Color

(Click image to go to info page)

As described in a previous post found here, the high-resolution (unbinned) luminance data for this showcase globular cluster was captured in early May 2011 from my favorite dark site at the Chiefland Astronomy Village. On my return trip to CAV in June 2011, I collected the low-resolution (binned 2x2) red, green, and blue data. Normally, I adjust the exposure times for each filter per results of manual G2V calibration that I performed when I first acquired the SBIG ST-8300M camera and then use 1/1/1 weighting when color combining. This time I used the same exposure (5-minutes) for each sub-frame for each filter then used the free application eXcalibrator to determine the RGB weighting for color combine. In short, it worked great! eXcalibrator found 13 G2V stars in the image field and from these determined a RGB weighting of 1.000/1.128/1.301 (very close to the results of my manual G2V calibration).

The color balance of the image above has not been manually adjusted from the initial color combine. Complete details regarding this image and more color images of M5 at various resolutions can be found by clicking on the image above or by clicking here.