Trying to overplot TESS and 2MASS images using WCSaxes and Astroquery

Question

I'm trying to plot a 2MASS image on the same projection as a TESS image. I used to do this with pywcsgrid2, but I can't seem to install it anymore. So I'm trying with Astropy WCSAxes.

Both the TESS image and 2MASS image are retrieved using Astroquery functions (TESScut and SkyView, respectively). I can create a plot of each image individually on its appropriate WCS axes. But when I try to plot the 2MASS image onto an axis with the TESS WCS (either as a contour or an image itself), it shrinks the image down into the lower left corner. Can somebody tell me what I'm doing wrong, or whether something might be off about the WCS.

import numpy as np
import matplotlib.pyplot as plt
from matplotlib import colors
from astroquery.skyview import SkyView
from astropy.wcs import WCS
import astropy.io.fits as fits
import astropy.units as u
from astropy.coordinates import SkyCoord
from astroquery.mast import Tesscut

test_tic = 113449635
cutout_coord = SkyCoord(113.29538299104,-35.48102008076,unit=u.degree)
dir_ffi = "./temp/"
tess_data = Tesscut.download_cutouts(cutout_coord, size=20, path=dir_ffi)

tess_file = tess_data[2][0]

def extract_wcs(dataheader):
    """
    Generate WCS for a TESS or K2 TPF. dataheader is the header of extension #1
    """
    w5 = WCS(naxis=2)

    w5.wcs.crpix = [dataheader["1CRPX5"],dataheader["2CRPX5"]]
    w5.wcs.cdelt = np.array([dataheader["1CDLT5"],dataheader["2CDLT5"]])
    w5.wcs.crval = [dataheader["1CRVL5"],dataheader["2CRVL5"]]
    w5.wcs.ctype = [dataheader["1CTYP5"],dataheader["2CTYP5"]]
    w5.wcs.pc = [[dataheader["11PC5"],dataheader["12PC5"]],
                 [dataheader["21PC5"],dataheader["22PC5"]]]
    return w5


with fits.open(tess_file) as hdu:

    # not a coadd, but reducing code for the sake of the MWE
    coadd = hdu[1].data["FLUX"][100]
    dataheader = hdu[1].header
    w3 = extract_wcs(dataheader)

# This works just fine
plt.subplot(projection=w3)
plt.imshow(coadd,origin="lower", norm=colors.LogNorm())
plt.xlabel('RA')
plt.ylabel('Dec')
plt.show()

# Get 2MASS K image
twomass_images = SkyView.get_images(position=cutout_coord, survey=['2MASS-K'],pixels=500)
pix_2mass = twomass_images[0][0].data
hdr_2mass = twomass_images[0][0].header
wcs_2mass = WCS(hdr_2mass)

print(w3)
print(wcs_2mass)

# Try to plot the TESS image in the background, with 2MASS overlaid as a contour
# This yields oversized axes with the image in the lower left corner
ax = plt.subplot(projection=w3)
ax.imshow(coadd,origin="lower", norm=colors.LogNorm(),zorder=-10)
ax.contour(pix_2mass, transform=ax.get_transform(wcs_2mass),
           colors='k',zorder=5)
plt.xlabel('RA')
plt.ylabel('Dec')
plt.show()


# Try to plot the 2MASS image in the background, on the TESS WCS
# This yields oversized axes with the image in the lower left corner
ax = plt.subplot(projection=w3)
ax.imshow(pix_2mass, origin="lower", norm=colors.LogNorm(), zorder=-10,
          transform=ax.get_transform(wcs_2mass))
plt.xlabel('RA')
plt.ylabel('Dec')
plt.show()