Improve the calculation of rotation angle at sub-km grid scales

[mgduda]

This PR improves the calculation of rotation angle at sub-km grid scales.

Previously, the sine and cosine of the rotation angle (SINALPHA and COSALPHA) were computed using finite differences of the latitude and longitude fields. This calculation generally worked well at grid scales larger than ~1 km and had the benefit of being independent of map projection. However, at smaller grid scales, the limited precision of the latitude and longitude fields led to numerical errors in the computed rotation angle.

This PR adds projection-specific code to compute the rotation angle. For the Mercator and cylindrical equidistant projections, the rotation angle is simply zero everywhere, and for the Lambert conformal conic and polar stereographic projections, the rotation angle is derived from the difference between a grid cell's longitude and the standard longitude. For the Cassini projection, a formula for the rotation angle is difficult to derive, and so with the expectation that this projection will find limited use at high resolution, the same finite differencing scheme is still employed.

[mgduda]

Using the following &geogrid namelist

&geogrid
 e_we = 300,
 e_sn = 300,
 dx = 500,
 dy = 500,
 map_proj  = 'lambert',
 ref_lat   =  40.00,
 ref_lon   = -105.25,
 truelat1  =  30.0,
 truelat2  =  60.0,
 stand_lon = -110.0,
/

here's an example of the SINALPHA field before and after the changes in this PR:

[mgduda]

@dudhia Can you think of anyone else who might be added as a reviewer for this PR? Do you recall who originally reported the issue with rotation angles at fine grid scales?

[dudhia]

add @weiwangncar who is also familiar with this issue. It came from a RAL team and we could pass it to them for testing.