import astropy.units as u
import matplotlib.pyplot as plt
import numpy as np
from astropy.coordinates import SkyCoord
from reproject import reproject_interp


import sunpy.map
from sunpy.data.sample import HMI_LOS_IMAGE


hmi = sunpy.map.Map(HMI_LOS_IMAGE)

new_data = np.empty((720, 1440))
new_frame = SkyCoord(0*u.deg, 0*u.deg, obstime=hmi.date, frame="heliographic_stonyhurst")
new_header = sunpy.map.make_fitswcs_header(new_data, new_frame,
                                           scale=[0.25, 0.25]*u.deg/u.pix,
                                           instrument="HMI",
                                           projection_code="CAR")
new_header['DATE-OBS'] = hmi.date.fits

outmap = sunpy.map.Map(new_data, new_header)

# Update outmap with the reprojected version of the hmi data
reproject_interp(hmi, outmap)

hmi.plot_settings['cmap'] = 'hmimag'
hmi.plot_settings['norm'] = plt.Normalize(-1500, 1500)
outmap.plot_settings = hmi.plot_settings


fig = plt.figure(figsize=(16,8))

ax1 = fig.add_subplot(1,2,1, projection=hmi)
hmi.plot(axes=ax1)
hmi.draw_grid()

ax2 = fig.add_subplot(1,2,2, projection=outmap)

outmap.plot(axes=ax2)
ax2.coords.grid(color='k')
plt.show()