Possible incorrect AACGMv2 usage

[RemingtonRohel]

Discussion

Topic: Possible incorrect usage of AACGMv2 library (not using correct flags for calculations).

Description

AACGMv2 is used in invmag.c, called by the make_grid binary. AACGM_v2_Convert can do conversions to/from AACGM coordinates, using either geodetic (default) or geocentric as the coordinate being converted from/to. In invmag.c, the usage of AACGM_v2_Convert uses the geodetic option (code & 16 == 1 indicating geocentric, called in this case with code=0), however it is clear from the code that geocentric coordinates are passed in.

Category

• User-friendly
• Software Design
• Data related
• Capabilities
• Clarity
• Workflow

Code

rst/codebase/superdarn/src.lib/tk/rpos.1.7/src/invmag.c

Lines 394 to 420 in a224011

	/* Convert range/beam position from geocentric latitude/longitude (flat,flon)
	* at virtual height (tmp_ht) to AACGM magnetic latitude/longitude
	* coordinates (mlat,mlon) */
	if (magflg == 2) {
	s=geod2ecdip(gdlat,gdlon,tmp_ht,out);
	*mlat = out[0];
	*mlon = out[1];
	} else if (magflg == 1) s=AACGMConvert(flat,flon,tmp_ht,mlat,mlon,&dummy,0);
	else s=AACGM_v2_Convert(flat,flon,tmp_ht,mlat,mlon,&dummy,0);
	if (s==-1) return -1;
	/* Calculate pointing direction latitude/longitude (xlat,xlon) given
	* distance (rsep) and bearing (azc) from the radar position (flat,flon)
	* at the field point radius (frho) */
	fldpnt_sph(frho,flat,flon,azc,rsep,&xlat,&xlon);
	/* Convert pointing direction position from geocentric latitude/longitude
	* (xlat,xlon) at virtual height (tmp_height) to AACGM magnetic
	* latitude/longitude coordinates (nlat,nlon) */
	if (magflg == 2) {
	geodtgc(-1,&gdlat,&gdlon,&gdrho,&xlat,&xlon,&dummy);
	s=geod2ecdip(gdlat,gdlon,tmp_ht,out);
	nlat = out[0];
	nlon = out[1];
	} else if (magflg == 1) s=AACGMConvert(xlat,xlon,tmp_ht,&nlat,&nlon,&dummy,0);
	else s=AACGM_v2_Convert(xlat,xlon,tmp_ht,&nlat,&nlon,&dummy,0);
	if (s==-1) return -1;

Documentation

AACGM_v2_Convert docstring:

/*-----------------------------------------------------------------------------
;
; NAME:
;       AACGM_v2_Convert
;
; PURPOSE:
;       Main function called by many SD plotting routines that are written
;       in C.
;
; CALLING SEQUENCE:
;       err = AACGM_v2_Convert(in_lat, in_lon, height,
;                 out_lat, out_lon, r, code);
;     
;     Input Arguments:  
;       int_lat       - double precision input latitude in degrees
;       int_lon       - double precision input longitude in degrees
;       height        - altitude in km
;       code          - bitwise code for passing options into converter
;                       G2A         - geographic (geodetic) to AACGM-v2
;                       A2G         - AACGM-v2 to geographic (geodetic)
;                       TRACE       - use field-line tracing, not coefficients
;                       ALLOWTRACE  - use trace only above 2000 km
;                       BADIDEA     - use coefficients above 2000 km
;                       GEOCENTRIC  - assume inputs are geocentric w/ RE=6371.2
;
;     Output Arguments:  
;       out_lat       - output latitude in degrees
;       out_lon       - output longitude in degrees
;       r             - geocentric radial distance in Re
;
;     Return Value:
;       error code
;
;
; NOTES:
;
;       All AACGM-v2 conversions are done in geocentric coordinates using a
;           value of 6371.2 km for the Earth radius.
;
;       For G2A conversion inputs are geographic latitude, longitude and
;           height (glat,glon,height), specified as either geocentric or
;           geodetic (default). For geodetic inputs a conversion to geocentric
;           coordinates is performed, which changes the values of
;           glat,glon,height. The output is AACGM-v2 latitude, longitude and
;           the geocentric radius (mlat,mlon,r) using the geocentric height
;           in units of RE.
;
;        For A2G conversion inputs are AACGM-v2 latitude, longitude and the
;            geocentric height (mlat,mlon,height). The latter can be obtained
;            from the r output of the G2A conversion. The output is geographic
;            latitude, longitude and height (glat,glon,height). If the
;            gedodetic option is desired (default) a conversion of the outputs
;            is performed, which changes the values of glat,glon,height.
;
;+-----------------------------------------------------------------------------
*/

[egthomas]

@RemingtonRohel thanks for identifying this! We took a quick look this afternoon, and agree the code value should probably be 16 when calling AACGM_v2_Convert (to indicate geocentric input) instead of 0. I'll try looking at the other places that function is called in the RST to see if a similar adjustment is needed.

In terms of scientific impact, my guess is the difference in calculated latitude / longitude is very small, especially relative to the ~110 km grid cell size.

[RemingtonRohel]

Good to know, I figured the difference would probably be small but I thought since I saw something I should say something.