creating maps using image FITS having "sin(deg)" as unit of latitude

[Sudha]

Hi,

I was wondering about the concept of "Sine Latitude" based map projection.

The example given in "Loading an HMI Synoptic map" uses the JSOC synoptic series with the FITS header 'cunit2' keyword indicating that the latitude is projected in units of 'Sine Latitude'. As shown in the example, the synoptic map is plotted.

When I try to use other JSOC HMI series like 

(1) Daily synchronic image FITS => hmi.Mrdailysynframe_720s

(2) Time-distance image FITS => hmi.td_fsi_12h,

the FITS header keyword 'cunit2' is 'sin(deg)' and sunpy does not recognize this unit and the following error is displayed:

ValueError: 'sin(deg)' did not parse as unit: At col 0, sin is not a valid unit.  If this is meant to be a custom unit, define it with 'u.def_unit'. To have it recognized inside a file reader or other code, enable it with 'u.add_enabled_units'. For details, see https://docs.astropy.org/en/latest/units/combining_and_defining.html

So I decided to try the following:

sin_deg = u.def_unit('sin_deg')
u.add_enabled_units([sin_deg])

For example, to test:

ang = 30*u.deg

<Quantity 30. deg>

q = np.sin(ang)*sin_deg

Quantity 0.5 sin_deg>

Then, I tried to create a new header using "sunpy.map.make_fitswcs_header" for this daily synchronic map by using the newly enabled unit 'sin_deg'

hdu = fits.open('hmi.Mrdailysynframe_720s.20140918_120000_TAI.data.fits')

data = hdu[0].data

coord = SkyCoord(180*u.arcsec, 0.sin_deg, obstime='2014-09-18', observer='earth', frame="heliographic_carrington")

UnitTypeError: Latitude instances require units equivalent to 'rad', so cannot set it to 'sin_deg'.

Could someone please help!

(1) I would like to know how sunpy recognized the latitude unit 'Sine Latitude' in synoptic file "hmi.Synoptic_Mr.*.fits' and why my newly created unit was not recognized. I might be doing something very wrong.

(2) I really do not know how to modify the header information of the synchronic FITs or time-distance FITs so that the sunpy map can be created. For example, the scale values in the header of daily synchronic file are

CDELT1  = -0.10000000000000001 / [degree/pixel] image scale in the x direction

CDELT2 = 0.0013888888888888889 / [sin(deg)/pixel] image scale in the y direction.

Thanks!

Sudha

[Albert Y. Shih]

Hi, Sudha,

     The short answer is that "Sine Latitude" is not the appropriate unit to use.  `CDELT2` should simply be an angle, because the cylindrical equal-area (CEA) projection will account for how the pixel size changes with latitude:

CUNIT2 = 'deg'

CDELT2 = 0.0013888888888888889 * 180 / np.pi

To see the analogous correction in SunPy, see here: https://github.com/sunpy/sunpy/blob/09d7d286b9fa4487c00493998b9485d6667082c4/sunpy/map/sources/sdo.py#L156-L164

(We don't actually do it this way anymore, but that's irrelevant to this discussion.)

Hope that helps!

Albert

--
You received this message because you are subscribed to the Google Groups "SunPy" group.
To unsubscribe from this group and stop receiving emails from it, send an email to sunpy+un...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/sunpy/b7e9ea2e-c8a1-4c07-a41f-8bcd5278e17cn%40googlegroups.com.

[Sudha]

Hi Albert ,

Thanks so much for your help! Making the map worked well, once the right changes were made. 

I think I've understood...I'll try and make longitude-sin(latitude) CEA maps with other datasets to get a better idea of it.

Thanks and regards,

Sudha.