Changeset 20527

Timestamp:

Nov 4, 2008, 3:42:28 AM (18 years ago)

Author:

Sebastian Jester

Message:

add logic for putting plots on same page using sm windows; catch change from FWHM_X to FWHM_MAJ etc.; propagate blankrad properly in gauss-fitting

File:

• branches/sj_ippTests_branch_20080929/ippTests/compIPPphoto.py (modified) (13 diffs)

branches/sj_ippTests_branch_20080929/ippTests/compIPPphoto.py

@@ -62 +62 @@
 # XXX Todo: labels! group by runs? different colours for different
 # bands; include seeing? plot things against seeing?
+
+# Also: the plotcol lists need to be grouped into things that are to
+# appear on the same page. For the beginning, that's the entire list
+# (and that's probably quite enough anyway).
+
 
 plotcol_1frame_tlist = [
@@ -107 +112 @@
 def compIPPphoto(summaryTable,mode,plotcol_1frame_tlist=plotcol_1frame_tlist,\
                      plotcol_summary_tlist=plotcol_summary_tlist,\
-                 cmfDir = '/IPP/data/SDSS/stripe82/coadd/prod/run_ipp_20080815/',\
-                 skip=True):
+                     cmfDir = '/IPP/data/SDSS/stripe82/coadd/prod/run_ipp_20080815/',\
+                     copyfields_list = ['RUN','RERUN','CAMCOL','FIELD','FILTER','FWHM_MAJ','FWHM_MIN'],\
+                     skip=True):
     """summaryTable: .fits table for output
     mode: new or append for creating summaryTable new or appending current run's output to it.
@@ -144 +150 @@
             filters.append(filter_name)
             bandindex_hash[filter_name] = bandindex
-        stats_hash, column_hash, goodrow_hash = computeStatistics(matchtable)
+        stats_hash, column_hash, goodrow_hash = computeStatistics(matchtable,copyfields_list = copyfields_list)
         # Sort res_hash by its column names to make sure order is
         # identical in all rows of the table, and output is more
@@ -284 +290 @@
     os.remove(tmp_outname)
     return all_outname
-        
+
+def setWindow(i,Nplots):
+    """Set sm window for plotting the i-th out of Nplots. The logic is
+    that N_x * N_y >= Nplots and N_y-N_x is 0 or 1. So Ny is the sqrt
+    of N or the next-greatest square; Nx is either the same or one
+    fewer if that's still enough."""
+    
+    from scipy import sqrt
+    if int(sqrt(Nplots))**2 == Nplots:
+        Ny = int(sqrt(Nplots))
+        Nx = Ny
+    else:
+        Ny = int(sqrt(Nplots)+1)
+        Nx = Ny-1
+        if Nx * Ny < Nplots:
+            Nx+=1
+
+    # That's nx, ny; now count up to total Nx ny. This is like
+    # representing i in a number system where the y-window number is
+    # the 10s and the x-window number the 1s, but the 10s are base Ny
+    # and the 1s are base Nx, and counting starts at 1 not 0. I.e. x =
+    # (i % Nx) or Nx if that's 0; y = i / Nx + 1
+    winx = i% Nx
+    if winx == 0:
+        winx = Nx
+    winy = int(i/Nx)+1
+    if int(i/Nx) == i/float(Nx):
+        winy -= 1
+    return Nx,Ny,winx,winy
 
 def valuesKeysSortedByKeys(hash):
@@ -379 +413 @@
     # "reference columns"
     outhash = {}
+    if 'FWHM_MAJ' in copyfields_list and 'FWHM_MAJ' not in h.ascardlist().keys():
+        copyfields_list.remove('FWHM_MAJ')
+        copyfields_list.remove('FWHM_MIN')
+        copyfields_list.append('FWHM_X')
+        copyfields_list.append('FWHM_Y')
     for f in copyfields_list:
         outhash[f] = h[f]
@@ -557 +596 @@
     filters = ['u','g','r','i','z']
     # Read primary  header of PS1 file to work out band
-    ps1copyfields_hash = headerfieldHash(['FWHM_MAJ','FWHM_MIN','FILTER'],PS1cmf,HDU=0)
+    try:
+        ps1copyfields_hash = headerfieldHash(['FWHM_MAJ','FWHM_MIN','FILTER'],PS1cmf,HDU=0)
+    except KeyError:
+        ps1copyfields_hash = headerfieldHash(['FWHM_X','FWHM_Y','FILTER'],PS1cmf,HDU=0)
     sdssbandstr = ps1copyfields_hash['FILTER']
     bandindex = filters.index(sdssbandstr)
@@ -832 +874 @@
     """Fit Gaussians to image in infile's HDU near the
     positions and with fluxes and sigmas given in initialGuessList"""
-    # XXX Todo: Make frame with residuals. Make blank array, then insert residuals.    
     import pyfits
     from numpy import nan
@@ -839 +880 @@
     f_handle=pyfits.open(infile)
     img = f_handle[HDU].data
-    resid_img = 1e-7*ones(img.shape,Float32)
+    resid_img = zeros(img.shape,Float32)
     fitlist = []
     for initialGuessTuple in initialGuessList:
@@ -860 +901 @@
             res = subim
         (scale,x,xsig,y,ysig,theta)=par
-        win.showArray(res)
+        # win.showArray(res)
         resid_img[int(y0)-halfwinsize:int(y0)+halfwinsize+1,int(x0)-halfwinsize:int(x0)+halfwinsize+1] = res
         # print initialGuessTuple, scale,x,xsig,y,ysig,theta
         fitlist.append(par)
     f_handle.close()
-    outname = 'gaussfit_resid_'+infile
+    outname = 'gaussfit_resid_r%ipix_%s' % (blankrad,infile)
     primHDU = pyfits.PrimaryHDU(resid_img)
     hdulist = pyfits.HDUList([primHDU])
@@ -871 +912 @@
     return fitlist
         
-def doEfbComparison(simfile='simtest.1.00a.fits',efbfile='efbsim.1.00a.fits',simtab='simtest.1.00a.cmf'):
+def doEfbComparison(simfile='simtest.1.00a.fits',efbfile='efbsim.1.00a.fits',simtab='simtest.1.00a.cmf',blankrad=9):
     import pyfits
     from numpy import array, rec
@@ -886 +927 @@
         initlist.append((counts[i],x[i],2.,y[i],2.,0.))
         inparlist.append((counts[i],x[i],y[i]))
-    efbparlist = doGaussPsfPhotometry(efbfile,initlist)
-    simparlist = doGaussPsfPhotometry(simfile,initlist)
+    efbparlist = doGaussPsfPhotometry(efbfile,initlist,blankrad=blankrad)
+    simparlist = doGaussPsfPhotometry(simfile,initlist,blankrad=blankrad)
 
     savelist = [ inparlist[i] + efbparlist[i] + simparlist[i] for i in range(nrows)]
@@ -900 +941 @@
     saveHDU = tabHDUfromRecArray(saverows)
     fitsend = re.compile('(\.[sc]mf|\.fits?)$')
-    saveHDU.writeto('gaussfit_'+fitsend.sub('',efbfile)+'___'+simfile,clobber=True)
+    saveHDU.writeto('gaussfit_r%ipix_%s___%s' %(blankrad,fitsend.sub('',efbfile),simfile),clobber=True)
     
 
@@ -931 +972 @@
         fpObjcDir += '/'
     # For testing:
-    # return ['1056-0192.421/1056-0192.421.ch.421.CHIP1.cmf'],['1056-0192.421/fpObjc-001056-1-0192.fit']
+    return ['1056-0192.421/1056-0192.421.ch.421.CHIP1.cmf'],['1056-0192.421/fpObjc-001056-1-0192.fit']
 
     from glob import glob