# commands extracted from tutorial writeup.
# Not all may be used during presentation

#-----------------------
# starting a interpreter
#-----------------------

python  # starts standard Python interpreter
ipython # starts ipython (enhanced interactive features)
pyraf   # starts PyRAF 

#--------------------------------------------
# some quick examples using pyfits and arrays
#--------------------------------------------

>>> import pyfits
>>> from numarray import *

>>> pyfits.info('pix.fits')
Filename: pix.fits
No.    Name         Type      Cards   Dimensions   Format
0    PRIMARY     PrimaryHDU      71  (512, 512)    Int16

>>> im = pyfits.getdata('pix.fits')
>>> import numdisplay
>>> numdisplay.display(im)
>>> numdisplay.display(im,z1=0,z2=300)
>>> fim = 1.*im
>>> bigvals = where(fim > 10)
>>> fim[bigvals] = 10*log(fim[bigvals]-10) + 10
>>> numdisplay.display(fim)
>>> hdr = pyfits.getheader('pix.fits')
>>> print hdr
SIMPLE  =                    T / Fits standard
BITPIX  =                   16 / Bits per pixel
NAXIS   =                    2 / Number of axes
NAXIS1  =                  512 / Axis length
NAXIS2  =                  512 / Axis length
EXTEND  =                    F / File may contain extensions
[...]
CCDPROC = 'Apr 22 14:11 CCD processing done'
AIRMASS =    1.08015632629395   / AIRMASS
HISTORY 'KPNO-IRAF'
HISTORY '24-04-87'
HISTORY 'KPNO-IRAF'           /
HISTORY '08-04-92'            /

>>> date = hdr['date']
>>> date
'2004-06-05T15:33:51'
>>> hdr['date'] = '4th of July'
>>> hdr.update('flatfile','flat17.fits')
>>> pyfits.writeto('newfile.fits',fim)    # Creates barebone header
>>> pyfits.writeto('newfile.fits',fim,hdr) # User supplied header
>>> pyfits.append('existingfile.fits',fim, hdr)
>>> pyfits.update('existingfile.fits',fim, hdr, ext=3)

#-----------------------------------------------------------
# Some comments about python variables, copies vs views, etc
#-----------------------------------------------------------
>>> value = 3
>>> print value
3
>>> value = "hello"

>>> x = im # the image we read in
>>> im[5,5] = -999
>>> x = im.copy()

#-----------------------------------
# And now a word about using objects
#-----------------------------------
>>> f = open('myfile.txt')
>>> f.name
>>> f.mode
>>> f.readline()
'a line from the text file'
f.seek(1024) # move 1024 bytes from the file beginning.

#---------------------------------
# An example of an error traceback
#--------------------------------- 
>>> f = pyfits.open(3)
Traceback (most recent call last):
  File "<stdin>", line 1, in ?
  File "/usr/stsci/pyssg/py/pyfits.py", line 3483, in open
    ffo = _File(name, mode=mode, memmap=memmap)
  File "/usr/stsci/pyssg/py/pyfits.py", line 2962, in __init__
    self.__file = __builtin__.open(name, python_mode[mode])
TypeError: coercing to Unicode: need string or buffer, int found

#-------------
# Array basics
#-------------

#----------------
# Creating arrays
#----------------
>>> x = zeros((20,30))
>>> arange(10)
array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9])
>>> arange(10.)
array([ 0., 1., 2., 3., 4., 5., 6., 7., 8., 9])
>>> arange(3,10)
array([3, 4, 5, 6, 7, 8, 9])
>>> arange(1., 10., 1.1) # note trickiness
array([1. , 2.1, 3.2, 4.3, 5.4, 6.5, 7.6, 8.7, 9.8])
>>> print array([3,1,7])
[3 1 7]
>>> print array([[2,3],[4,4]])
[[2 3]
 [4 4]]

#------------ 
# Array types
#------------
>>> zeros(3, Int8) # Signed byte
>>> zeros(3, type=UInt8) # Unsigned byte
>>> array([2,3], type=Float32)
>>> arange(4, type=Complex64)
>>> arr.type()
Float32
>>> a = arr.astype(Float64)

#----------------------
# Printing large arrays
#----------------------
>>> x = arange(1000000)
print x
[    0     1     2 ..., 999997 999998 999999]

#----------------
# Indexing arrays
#----------------
>>> x = arange(10)
>>> x
array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9])
>>> x[2] # 3rd element
2
>>> x[-2] # -1 represents the last element, -2 next to last...
8
>>> x[2:5]
array([2, 3, 4])
>>> y = x[2:5]
>>> y[0] = 99
>>> y
array([99, 3, 4])
>>> x
array([0, 1, 99, 3, 4, 5, 6, 7, 8, 9])

>>> x[:5] # presumes start from beginning
array([ 0, 1, 99, 3, 4])
>>> x[2:] # presumes goes until end
array([99, 3, 4, 5, 6, 7, 8, 9])
>>> x[:] # selects whole dimension
array([0, 1, 99, 3, 4, 5, 6, 7, 8, 9])
>>> x[2:8:3] # "Stride" every third element
array([99, 5])
>>> x[[4,2,4,1]]
array([4, 99, 4, 1])
>>> x[where(x>5)]
array([99, 6, 7, 8, 9])
>>> m = x > 5
>>> m
array([0,0,1,0,0,0,1,1,1,1], type=Bool)
>>> x[m]
array([99, 6, 7, 8, 9])

#------------------------
# Multidimensional arrays
#------------------------
>>> im = arange(24)
>>> im.shape=(4,6)
>>> im
array([[ 0,  1,  2,  3,  4,  5],
       [ 6,  7,  8,  9, 10, 11],
       [12, 13, 14, 15, 16, 17],
       [18, 19, 20, 21, 22, 23]])
>>> im[1, 2]
8
>>> im[1]
array([6, 7, 8, 9, 10, 11])

#-------------------
# Broadcasting rules
#-------------------
>>> x[:5] = 0
>>> x[:5] = array([0,1,2,3,4]) # OK
>>> x[:5] = array([0,1,2,3]) 
!Error
>>> x = zeros((5,4))
>>> x[:,:] = [2,3,2,3]
>>> x
array([[2, 3, 2, 3],
       [2, 3, 2, 3],
       [2, 3, 2, 3],
       [2, 3, 2, 3],
       [2, 3, 2, 3]])


#-------
# Ufuncs
#-------
>>> x = arange(9)
>>> x.shape = (3,3)
>>> x
array([0, 1, 2],
      [3, 4, 5],
      [6, 7, 8]])
>>> add.reduce(x) # sums along the first index
array([9, 12, 15])
>>> add.reduce(x, axis=1) # sums along the 2nd index
array([3, 12, 21])
>>> add.accumulate(x) # cumulative sum along the first index
array([[0,  1,  2],
       [3,  5,  7],
       [9, 12, 15]])
>>> multiply.outer(arange(3),arange(3))
array([[0, 0, 0],
       [0, 1, 2],
       [0, 2, 4]])
>>> x.sum() # the L indicates a Python Long integer
36L 

#-----------------------------
# Practical example
# Summing flux within a radius
#   of a position
#-----------------------------
# first find location of maximum in image
y, x = indices(im.shape, type=Float32) # note order of x, y!
# after finding peak at 257,258 using ds9
x = x-257 # make peak have 0 value
y = y-258
radius = sqrt(x**2+y**2)
mask = radius < 50
display mask*im
(mask*im).sum() # sum total of masked image
# or
im[where(mask)].sum() # sum those points within the radius
# look Ma, no loops!