<jupyter_start><jupyter_text>the stars module : interactive star marking, aperture photometry, WCS, plus photometric transformations Start with the usual setup<jupyter_code>%matplotlib qt
import matplotlib.pyplot as plt
import numpy as np
import glob
import pdb
from pyvista import imred, tv, stars<jupyter_output>pyautogui does not seem to be available, disabling arrow key cursor moves<jupyter_text>Instantiate a reducer object, and a TV object<jupyter_code>red=imred.Reducer('TM61',dir='/home/tmo/20211020')
t=tv.TV()
t.tvclear()
a=red.reduce(851)
t.tv(a)
red.scale<jupyter_output><empty_output><jupyter_text>Interactive marking of stars with stars.mark(). Note that photutils has routines for automatic detection of stars.<jupyter_code>seeing=3
tab=stars.mark(t,rad=seeing/red.scale,new=True)
print(tab)<jupyter_output>id    x       y        MJD      EXPTIME FILTER  AIRMASS 
--- ------- ------- ------------ ------- ------ ---------
  1 3344.90 1448.30 59507.536177    10.0     SR 1.0339353
  2 1996.00 1514.28 59507.536177    10.0     SR 1.0339353<jupyter_text>Note that if you provide an input table to stars.mark(), along with auto=True, it will re-centroid around position in table Aperture photometry with stars.photom (front end for photutils). This will do aperture photometry on image specified in first argument, for objects specified in the second argument ( a Table). Will do multiple aperture radii as specified with rad= kekyword, using sky annulus as specified by skyrad= keyword.  If display= is given, apertures and annuli will be displayed on the image.Returns the table with columns added for aperture photometry, uncertainties and sky values. Magnitudes are instrumental magnitudes, not corrected for exposure time, i.e. -2.5 log (counts)<jupyter_code>stars.photom(a,tab,rad=np.arange(3,10),skyrad=[20,30],display=t)
print(tab)<jupyter_output>id    x       y        MJD      EXPTIME ... aper9.0 aper9.0err   sky    skysig
--- ------- ------- ------------ ------- ... ------- ---------- -------- ------
  1 3344.90 1448.30 59507.536177    10.0 ... -15.123      0.002 12024.45 125.17
  2 1996.00 1514.28 59507.536177    10.0 ... -13.997      0.006 12960.36 129.39<jupyter_text>World coordinate system (WCS) With solve=True, Reducer will try to solve for a WCS using a local instllation of astrometry.net Example for autofinding stars on a TMO image. Here's an example of using a master file with a set of standard stars, using WCS to get pixel positions and selecting only those stars on the current frame.<jupyter_code>#reduce an image and add WCS using astrometry.net routines
# note that I haven't even included flat fielding! so results will be poor
red.dir=('/home/tmo/20211116')
a=red.reduce('SA112-S001-SR-R001-SR.fts',solve=True)

# ugriz_smith file has positions and standard mags
tab=Table.read('ugriz_smith.fits')

# populate x and y columns with pixel locations in this image
tab['x'],tab['y']=a.wcs.wcs_world2pix(tab['RA'],tab['DEC'],0)
#table formatting
for col in ['x','y','u','g','r','i','z'] : tab[col].info.format='.2f'
for col in ['RA','DEC'] : tab[col].info.format='.6f'

# select out stars that fall on the image
nrow,ncol=a.data.shape
gd=np.where((tab['x']>0)&(tab['y']>0)&(tab['x']<ncol)&(tab['y']<nrow))[0]
tab=tab[gd]

# display them
t.tvclear()
t.tv(a)
for xx,yy in zip(tab['x'],tab['y']) : 
    t.tvcirc(xx,yy,20)
    
# show the table
tab<jupyter_output>Reading file: /home/tmo/20211116/SA112-S001-SR-R001-SR.fts
  subtracting overscan:  512.8631917962498
  plate solving with local astrometry.net....<jupyter_text>This has the nice advantage that you can just add the photometry onto this table and facilitate getting the photometric transformations!Now we use the x,y positions from the WCS, recentroid using mark() with the auto=True option (exit=True means it doesn't then wait for you to mark additional stars), and then do the aperture photometry.<jupyter_code>t.tv(a)
phot=stars.mark(t,tab,rad=3/red.scale,auto=True,exit=True)
stars.photom(a,phot,rad=[10],skyrad=[20,30],display=t)<jupyter_output>/home/local/Anaconda3-2020.02/lib/python3.7/site-packages/astropy/table/table.py:3064: FutureWarning: elementwise == comparison failed and returning scalar instead; this will raise an error or perform elementwise comparison in the future.
  result = self.as_array() == other<jupyter_text>Note that the photom() routine does not include the exposure time in the instrumental magnitudes. Photometric transformations The general transformation equations follow:$$g = SG + k_g X + t_g(g-r) + z_g$$$$r = SR + k_r X + t_r(g-r) + z_r$$$$...$$where there would be an equation for every filter. For each filter, we need to associate the standard magnitude to which we want to transform (e.g., SG to g, SR to r, etc.), and also associate a standard color to use for the color term, which is typically a color index between the magnitude that we are transforming to, and an adjacent one (e.g., g-r for both g and r above).Using standard stars, we solve for the extinction coefficient, the transformation coefficient, and the zero point for each filter. Generally, this is set up as a linear least squares problem for each filter:\begin{equation}\begin{pmatrix}X_1 & color_1 & 1 \\X_2 & color_2 & 1 \\... & ... & ... \\\end{pmatrix}\begin{pmatrix}k_g\\t_g\\z_g\\\end{pmatrix} =\begin{pmatrix}g-SG_1\\g-SG_2\\...\\\end{pmatrix}\end{equation}where there is a row for every observation in the filter.<jupyter_code>%matplotlib inline
from tools import plots
filts=['SR']        #observed filter names
mags=['r']          # standard mag to transform to
colors=[['g','r']]  # pair for each filter
aper='aper10'
for filt,mag,color in zip(filts,mags,colors) :
    j=np.where(tab['FILTER'] == filt)[0]
    mjd=tab['MJD'][j]
    airmass=tab['AIRMASS'][j]
    exptime=tab['EXPTIME'][j]
    col=tab[color[0]][j]-tab[color[1]][j]
    zero=np.ones(len(j))
    # include the EXPTIME correction for instrumental mag
    inst=tab[aper]+2.5*np.log10(exptime)
    stan=tab[mag]
    # right hand side is standard mag - instrumental mag
    rhs=stan-inst
    # left hand side matrix is the airmass,color,1 for each observation
    # for my example, I don't have multiple airmass, so just including 2 terms
    #mat=np.array([airmass,color,zero]).T
    mat=np.array([col,zero]).T
    
    #plot difference between standard and observe a function of things
    #plots (from holtz-tools) contains convenience routines for making and labelling plots
    fig,ax=plots.multi(2,2,figsize=(8,5))
    # plotc plots points at x,y positions, colored by z
    colname=color[0]+'-'+color[1]
    plots.plotc(ax[0,0],airmass,rhs,col,
                xt='airmass',yt=mag+'-'+filt,zt=colname,colorbar=True,size=20)
    plots.plotc(ax[0,1],col,rhs,airmass,
                xt='color',yt=mag+'-'+filt,zt='airmass',colorbar=True,size=20)
    plots.plotc(ax[1,0],mjd,rhs,airmass,
                xt='MJD',yt=mag+'-'+filt,zt='airmass',colorbar=True,size=20)
    fig.tight_layout()
    
    # do the fit for the coefficients
    fit,resid,rank,s=np.linalg.lstsq(mat,rhs)
    print('fit: ',fit)
    print('resid ',resid)
    
    # now use the solution to get transformed magnitudes from the obs
    trans=np.matmul(mat,fit)+obs
    plots.plotc(ax[1,1],mjd,stan-trans,airmass,
                xt='MJD',yt='stan-trans',zt='airmass',colorbar=True,size=20)<jupyter_output>/home/local/Anaconda3-2020.02/lib/python3.7/site-packages/ipykernel_launcher.py:38: FutureWarning: `rcond` parameter will change to the default of machine precision times ``max(M, N)`` where M and N are the input matrix dimensions.
To use the future default and silence this warning we advise to pass `rcond=None`, to keep using the old, explicitly pass `rcond=-1`.
{'fraction_non_alphanumeric': 0.1350407450523865, 'fraction_numerical': 0.04203854611305135, 'mean_word_length': 4.05359477124183, 'pattern_counts': {'":': 0, '<': 25, '<?xml version=': 0, '>': 25, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '7', 'dir': 'jupyter-structured-clean-dedup', 'id': '172800', 'license': 'MIT', 'n_tokens_mistral': 2782, 'n_tokens_neox': 2510, 'n_words': 783, 'path': 'doc/stars.ipynb', 'repo_name': 'holtzmanjon/pyvista', 'size': '7731'}