title: Radial Profile from a .fits image
body: <p>I have been trying to plot a radial profile of a fits image using a modified script I found on-line. I always get y axis units which are completely different to what's expected. I'm not even sure what the y axis units are. I have attached the fits file and a profile I keep getting and the correct radial profile I plotted using another program. </p>
  
  <p>I am very new to python so I have no idea why this keeps happening. Any help to fix this will be so greatly appreciated.  </p>
  
  <p>This is the code I've been using: </p>
  
  <pre><code>import numpy as np
  import pyfits
  import matplotlib.pyplot as plt
  import matplotlib.image as mpimg
  
  def azimuthalAverage(image, center=None):
      """
      Calculate the azimuthally averaged radial profile.
  
      image - The 2D image
      center - The [x,y] pixel coordinates used as the center. The default is 
               None, which then uses the center of the image (including 
               fracitonal pixels).
  
      """
      # Calculate the indices from the image
      y, x = np.indices(image.shape)
  
  
      if not center:
          center = np.array([(x.max()-x.min())/2.0, (y.max()-y.min())/2.0])
  
      r = np.hypot(x - center[0], y - center[1])
  
      # Get sorted radii
      ind = np.argsort(r.flat)
      r_sorted = r.flat[ind]
      i_sorted = image.flat[ind]
  
      # Get the integer part of the radii (bin size = 1)
      r_int = r_sorted.astype(int)
  
      # Find all pixels that fall within each radial bin.
      deltar = r_int[1:] - r_int[:-1]  # Assumes all radii represented
      rind = np.where(deltar)[1]       # location of changed radius
      nr = rind[1:] - rind[:-1]        # number of radius bin
  
      # Cumulative sum to figure out sums for each radius bin
      csim = np.cumsum(i_sorted, dtype=float)
      tbin = csim[rind[1:]] - csim[rind[:-1]]
  
      radial_prof = tbin / nr
      print center
      print i_sorted
      print radial_prof
      return radial_prof
  
  #read in image
  hdulist = pyfits.open('cit6ndf2fitsexample.fits')
  scidata = np.array(hdulist[0].data)[0,:,:]
  center = None
  radi = 10
  rad = azimuthalAverage(scidata, center)
  
  plt.xlabel('radius(pixels?)', fontsize=12)
  plt.ylabel('image intensity', fontsize=12)
  plt.xlim(0,10)
  plt.ylim(0, 3.2)
  plt.plot(rad[radi:])
  plt.savefig('testfig1.png')
  plt.show()
  </code></pre>
  
  <p>Profile with wrong y axis units</p>
  
  <p><a href="https://i.stack.imgur.com/kgJk9.png" rel="noreferrer"><img src="https://i.stack.imgur.com/kgJk9.png" alt="enter image description here"></a></p>
  
  <p>Profile with expected correct units created using Celtech Aperture Photometry Tool. </p>
  
  <p><a href="https://i.stack.imgur.com/RFyZi.png" rel="noreferrer"><img src="https://i.stack.imgur.com/RFyZi.png" alt="enter image description here"></a></p>