EnPM
EnPM

Reputation: 1

aperture photometry with astropy photutils using a set of predefined apertures

I am using photutils to perform isophotal fitting on the image of a galaxy.


#First a dummy iteration to get the right centre
geometry = EllipseGeometry(x0=301, y0=307, sma=20, eps=0.8, pa=np.deg2rad(110))
ellipse = Ellipse(masked_img, geometry) 
isolist0 = ellipse.fit_image(sma0 = x00, minsma=4, maxsma=210, step=8, maxit=100, minit=10, linear=True, fix_center=False) 

#second iteration with center fixed
geometry = EllipseGeometry(x0=np.median(isolist0.x0), y0=np.median(isolist0.y0), sma=20, eps=0.8, pa=np.deg2rad(110))
ellipse = Ellipse(masked_img, geometry) 
isolist = ellipse.fit_image(sma0 = x00, minsma=4, maxsma=210, step=8, maxit=200, minit=20, linear=True, nclip=3, integrmode = 'median', fix_center=True)

So far, so good. The result is an isophote list that contains information on the aperture parameters (e.g., radius, center, position angle, ellipticity), as well as information on the flux density profile (specifically the mean intensity and its error, see https://photutils.readthedocs.io/en/stable/api/photutils.isophote.IsophoteList.html#photutils.isophote.IsophoteList).

Now, I have a second image (the same galaxy in a different band), say masked_img_2. I want to obtain the flux profile (isolist.intens and isolist.int_err) for masked_img_2, but without doing the isophotal fitting again; instead, I want to use the same apertures as in my existing isolist object. In other words, I want to do aperture photometry matching the apertures obtained before with the isophotal fitting.

I have tried a few different things but nothing quite works. What came the closest was this:

n_intensities = []
mean_intensities_err = []


# Loop over each isophote in the isolist
for isophote in isolist:
    # Extract parameters
    sma = isophote.sma  # Semi-major axis
    eps = isophote.eps  # Ellipticity
    pa = isophote.pa    # Position angle in radians
    x0 = isophote.x0    # X center
    y0 = isophote.y0    # Y center

    # Calculate semi-minor axis
    sma_minor = sma * (1 - eps)

    # Create the EllipticalAperture for the current isophote
    ellipse = EllipticalAperture((x0, y0), sma, sma_minor)

    # Calculate points along the ellipse
    theta = np.linspace(0, 2 * np.pi, 100)  # 100 points around the ellipse
    x_coords = x0 + sma * np.cos(theta) * np.cos(pa) - sma_minor * np.sin(theta) * np.sin(pa)
    y_coords = y0 + sma * np.cos(theta) * np.sin(pa) + sma_minor * np.sin(theta) * np.cos(pa)

    # Clip coordinates to be within image boundaries
    x_coords = np.clip(x_coords, 0, masked_img_2.shape[1] - 1).astype(int)
    y_coords = np.clip(y_coords, 0, masked_img_2.shape[0] - 1).astype(int)

    # Extract pixel values along the elliptical path
    ellipse_data = masked_img_2[y_coords, x_coords]

    # Calculate mean intensity along the elliptical path
    if ellipse_data.size > 0:
        mean_intensity = np.nanmean(ellipse_data)
        int_err = np.nanstd(ellipse_data) / np.sqrt(ellipse_data.size)  # RMS error   
        
    else:
        mean_intensity = 0  # Or some value to indicate no data
        int_err = 0

    mean_intensities.append(mean_intensity)
    mean_intensities_err.append(int_err)

To test it, I input masked_img instead of masked_img_2, just to see that I can recover the same as when using the isophotal fitting. The results are similar for the intensities but they are not perfect, and can lead to difference of ~0.2 in the surface brightness profiles.

Plot showing that the recovery is not perfect

So, my method does not fully work. Moreover, I would think that there is a much simpler way of using direct photutils to perform matched aperture photometry. Any help is greatly appreciated!

Upvotes: 0

Views: 190

Answers (0)

Related Questions