Azimuthal Equidistant projection: IDL result reading the HDFLook product geotiff file

Overlap: Red and Green drawings from HDFLook, white drawing from IDL
Click on the image to enlarge.

Azimuthal Equidistant projection IDL shell

; IDL SHELL READ A GEOTIFF FILE ( AZIMUTHAL PROJECTION ) AND ADD A MAP
;
image = READ_TIFF('~/IDL/DATA_HDFLOOK/GEOTIFF_AZIMUT.tif', GEOTIFF = geotag,  /VERBOSE)

sz=reverse(SIZE(image, /DIM))

limit_meters=DBLARR(8)
limit_meters[0:1]=[geotag.modeltiepointtag[3], geotag.modeltiepointtag[4]]
scalex=geotag.modelpixelscaletag[0]
scaley=geotag.modelpixelscaletag[1]
 
                                                                                                            
limit_meters[2]=limit_meters[0]+sz[0]*scalex
limit_meters[3]=limit_meters[1]
                                                                                                                                           
limit_meters[4]=limit_meters[0]+sz[0]*scalex
limit_meters[5]=limit_meters[1]-sz[1]*scaley
                                                                                                                                           
limit_meters[6]=limit_meters[0]
limit_meters[7]=limit_meters[1]-sz[1]*scaley
                                                                                                                                           

window, xsize=500, ysize=500
tv, reverse(image, 3), true=1

print, 'Center Longitude', geotag.PROJCENTERLONGGEOKEY

res = map_proj_init('azimuthal', /gctp,  center_latitude = geotag.PROJCENTERLATGEOKEY,  center_longitude = geotag.PROJCENTERLONGGEOKEY,  sphere_radius = geotag.GEOGSEMIMAJORAXISGEOKEY)
limit = map_proj_inverse(limit_meters, map_structure=res)

limit=reform(limit,2,4)
limit=REVERSE(limit)

map_set, geotag.PROJCENTERLATGEOKEY, geotag.PROJCENTERLONGGEOKEY, pos=[0,0,1,1], /azimuthal,  limit=limit, /hires, /noborder, /noerase
map_grid, londel=2, latdel=2
map_continents, /hires, /coast