;-------------------------------------------------------------
;+
; NAME:
;       MAP_SET_SCALE
; PURPOSE:
;       Set map scaling from info embedded in a map image.
; CATEGORY:
; CALLING SEQUENCE:
;       No args.
; INPUTS:
; KEYWORD PARAMETERS:
;       Keywords:
;         NEWPOS=pos2  Give a new position array, overrides
;           existing position.  Good for resizing maps.
;           pos2 must be in device coordinates (not normalized).
;         SCALE_FACTOR=fact Correction factor for map scale.
;         COLOR=clr Color to override border color if any.
;         OUT=out  Scaling info returned as a structure.
;         COMMAND=cmd  Returned map command.
;         /LIST  List values.
;         IMAGE=img  Give image array instead of reading it from
;           the display.
;         INFO=info Give the embedded info byte array instead of
;           reading it from the display or image.
; OUTPUTS:
; COMMON BLOCKS:
; NOTES:
;       Notes: Uses info embedded on bottom image line by
;       map_put_scale, if available.
; MODIFICATION HISTORY:
;       R. Sterner, 1999 Sep 21
;       R. Sterner, 2002 Jan 18 --- complete rewrite to work with map_set2.
;       R. Sterner, 2002 Jan 29 --- Added clip=0 flag.
;       R. Sterner, 2002 Feb 08 --- Added COLOR=clr for border.
;       R. Sterner, 2002 Feb 10 --- Returned pixels/degree in structure.
;       R. Sterner, 2002 Apr 02 --- Return map_set2 command from embedded info.
;       R. Sterner, 2002 Apr 05 --- Fixed for IMAGE=img case.
;       R. Sterner, 2002 Jul 03 --- Added map pixel size to returned structure.
;       R. Sterner, 2002 Jul 31 --- Added map central lat,long to list & struct.
;
; Copyright (C) 1999, Johns Hopkins University/Applied Physics Laboratory
; This software may be used, copied, or redistributed as long as it is not
; sold and this copyright notice is reproduced on each copy made.  This
; routine is provided as is without any express or implied warranties
; whatsoever.  Other limitations apply as described in the file disclaimer.txt.
;-
;-------------------------------------------------------------
	function map_proj_name, proj, help=hlp
 
	if (n_params(0) lt 1) or keyword_set(hlp) then begin
	  print,' Convert map_set projection number to projection name.'
	  print,' name = map_proj_name( proj)'
	  print,'   proj = projection number.  in'
	  print,'   name = projection name.    out'
	  print,' Note: projection number is the number set by the'
	  print,' map_set command and may be found in !map.projection'
	  print,' after map_set is executed.'
	  return,''
	endif
 
	names = ['', $
		'Stereographic', $
		'Orthographic', $
		'Lambert Conic', $
		'Lambert Azimuthal', $
		'Gnomic', $
		'Azimuthal Equidistant', $
		'Satellite', $
		'Cylindrical', $
		'Mercator', $
		'Mollweide', $
		'Sinusoidal', $
		'Aitoff', $
		'Hammer Aitoff', $
		'Albers Equal Area Conic', $
		'Transverse Mercator', $
		'Miller Cylindrical', $
		'Robinson', $
		'Lambert Conic Ellipsoid', $
		'Goodes Homolosine', $
		'']
 
	return, (names([proj]))(0)
 
	end
 
 
;===================================================================
;	map_set_scale = Set map scaling from info embedded in a map image.
;	R. Sterner, 1999 Aug 31
;===================================================================
 
	pro map_set_scale, list=list, image=t0, info=info, color=pclr, $
	  newpos=newpos, scale_factor=fact, out=out, help=hlp, $
	  command=cmdtxt
 
	if keyword_set(hlp) then begin
	  print,' Set map scaling from info embedded in a map image.'
 	  print,'   No args.'
	  print,' Keywords:'
	  print,'   NEWPOS=pos2  Give a new position array, overrides'
	  print,'     existing position.  Good for resizing maps.'
	  print,'     pos2 must be in device coordinates (not normalized).'
	  print,'   SCALE_FACTOR=fact Correction factor for map scale.'
	  print,'   COLOR=clr Color to override border color if any.'
	  print,'   OUT=out  Scaling info returned as a structure.'
	  print,'   COMMAND=cmd  Returned map command.'
	  print,'   /LIST  List values.'
	  print,'   IMAGE=img  Give image array instead of reading it from'
	  print,'     the display.'
	  print,'   INFO=info Give the embedded info byte array instead of'
	  print,'     reading it from the display or image.'
 	  print,' Notes: Uses info embedded on bottom image line by'
	  print,' map_put_scale, if available.'
	  return
	endif
 
	;****************************************************************
	;	Extract embedded info from image and interpret
	;****************************************************************
 
	;=============================================================
	;	Make sure map coordinate info is embedded in image
	;=============================================================
	if n_elements(info) eq 0 then begin	; Info array not given.
	  if n_elements(t0) eq 0 then begin	; Image not given.
	    if !d.x_size lt 160 then return	  ; Image too small.
	    t0 = tvrd(true=3)			  ; Read image.
	    d_x_size = !d.x_size		  ; Window size (=current).
	    d_y_size = !d.y_size
	    wintxt = strtrim(!d.window,2)
	  endif else begin			; Image was given.
	    img_shape, t0, nx=d_x_size, ny=d_y_size	; Window size (=image).
	    wintxt = 'image'
	  endelse
	  img_split, t0, r, g, b		; Split image into RGB.
	  t = b(0:159)			  	; Grab info from B.
	endif else begin			; Info array given.
	  t = info				; Copy info from given array.
	  d_x_size = !d.x_size		  	; Window size (=current).
	  d_y_size = !d.y_size
	  wintxt = strtrim(!d.window,2)
	endelse
	if n_elements(t) lt 160 then return	; Check size of info array.
	m = string(t(0:9))			; Data available flag.
	if m ne '1234567891' then return	; No map scaling info in image.
 
	;=============================================================
	;	Pick info out of byte array
	;=============================================================
	tmp = float(t,10,3) & lat=tmp(0) & lon=tmp(1) & ang=tmp(2)
	flag_noborder = t(22)
	flag_iso = t(23)
	flag_col = t(24)
	if flag_col then clr=long(t,25) else clr=!p.color
	flag_scale = t(29)
	scale = float(t,30)
	flag_azi = t(34)
	azi = float(t,35)
	flag_lim = t(39)
	nlim = long(t,40)
	lim = float(t,44,8)
	lim = lim(0:nlim-1)
	flag_par = t(76)
	std_par = float(t,77,2)
	flag_sat = t(85)
	sat_p = float(t,86,3)
	flag_ell = t(98)
	ellip = float(t,99,3)
	pos = float(t,111,4)
	proj = long(t,127)
	flag_clip = t(131)
 
	;=============================================================
	;	Deal with new POSITION if any
	;=============================================================
	if n_elements(newpos) eq 4 then begin
	  pos = float(newpos)/[!d.x_size,!d.y_size,!d.x_size,!d.y_size]
	endif
 
	;=============================================================
	;	Deal with map scale change.
	;=============================================================
	if n_elements(fact) gt 0 then begin
	  if flag_scale eq 1 then begin
	    scale = fact*scale
	  endif
	endif
 
	;=============================================================
	;	Deal with new color if any
	;=============================================================
	if n_elements(pclr) ne 0 then clr=pclr
 
	;****************************************************************
	;	Use extracted info to set up map_set command
	;****************************************************************
 
	;=========================================
	;	Set up basic map_set command
	;	Must use map_set2 to update
	;	  scaling info.
	;=========================================
	cmd = 'map_set2,/noerase,lat,lon,ang,proj=proj,pos=pos,col=clr'
	if arg_present(cmdtxt) then begin
	  cmdtxt = ['lat = '+strtrim(lat,2)]
	  cmdtxt = ['lon = '+strtrim(lon,2),cmdtxt]
	  cmdtxt = ['ang = '+strtrim(ang,2),cmdtxt]
	  cmdtxt = ['proj = '+strtrim(proj,2)+'  ; '+ $
	    map_proj_name(proj),cmdtxt]
	  cmdtxt = ['pos = ['+commalist(pos)+']',cmdtxt]
	  cmdtxt = ['clr = '+strtrim(clr,2),cmdtxt]
	endif
 
	;=========================================
	;	Deal with /ISO
	;=========================================
	if flag_iso eq 1 then $
	  cmd = cmd + ',/iso'			; Add /ISO.
 
	;=========================================
	;	Deal with /NOBORDER
	;=========================================
	if flag_noborder eq 1 then $
	  cmd = cmd + ',/noborder'		; Add /NOBORDER.
 
	;=========================================
	;	Deal with CLIP=0
	;=========================================
	if flag_clip eq 1 then $
	  cmd = cmd + ',clip=0'			; Add CLIP=0.
 
	;=========================================
	;	Deal with LIMIT=[...]
	;=========================================
	if flag_lim then begin
	  cmd = cmd+',limit=lim'		; Add limit.
	  if arg_present(cmdtxt) then begin
	    cmdtxt = ['lim = ['+commalist(lim)+']',cmdtxt]
	  endif
	endif
 
	;=========================================
	;	Deal with SCALE
	;=========================================
	if flag_scale then begin
	  cmd = cmd+',scale=scale'
	  if arg_present(cmdtxt) then begin
	    cmdtxt = ['scale = '+strtrim(scale,2),cmdtxt]
	  endif
	endif
 
	;=========================================
	;	Deal with CENTRAL_AZI
	;=========================================
	if flag_azi then begin
	  cmd = cmd+',central_azi=azi'
	  if arg_present(cmdtxt) then begin
	    cmdtxt = ['azi = '+strtrim(azi,2),cmdtxt]
	  endif
	endif
 
	;=========================================
	;	Deal with Standard Parallels
	;=========================================
	if flag_par then begin
	  cmd = cmd+',standard_par=std_par'
	  if arg_present(cmdtxt) then begin
	    cmdtxt = ['std_par = ['+commalist(std_par)+']',cmdtxt]
	  endif
	endif
 
	;=========================================
	;	Deal with Satellite
	;=========================================
	if flag_sat then begin
	  cmd = cmd+',sat_p=sat_p'
	  if arg_present(cmdtxt) then begin
	    cmdtxt = ['sat_p = ['+commalist(sat_p)+']',cmdtxt]
	  endif
	endif
 
	;=========================================
	;	Deal with Ellipsoid
	;=========================================
	if flag_ell then begin
	  cmd = cmd+',ellips=ellip'
	  if arg_present(cmdtxt) then begin
	    cmdtxt = ['ellip = ['+commalist(ellip)+']',cmdtxt]
	  endif
	endif
 
	;=========================================
	;	Execute the map_set command
	;=========================================
	err = execute(cmd)
	if err ne 1 then begin
	  print,' Error in map_set_scale: command not executed correctly.'
	  print,' Command was:'
	  print, cmd
	endif else begin
	  if arg_present(cmdtxt) then begin
	    txt = ';-----  Window size: xsize='+strtrim(d_x_size,2)+$
	      ', ysize='+strtrim(d_y_size,2)+'  --------'
	    cmdtxt = [txt,reverse(cmdtxt),cmd]
	  endif
	endelse
 
	;=======================================
        ;       Extract some values
        ;=======================================
	ix1 = pos(0)*d_x_size
	ix2 = pos(2)*d_x_size & idx=ix2+0-ix1
	iy1 = pos(1)*d_y_size
	iy2 = pos(3)*d_y_size & idy=iy2+0-iy1
 
        ;=======================================
	;	Central map scale (pixels/degree)
	;	Also true central long, lat
        ;=======================================
	ixmd=(ix1+ix2)/2. & iymd=(iy1+iy2)/2. & iyy1=iymd-10 & iyy2=iymd+10
	tmp=convert_coord([ixmd,ixmd],[iyy1,iyy2],/dev,/to_data)
	x=tmp(0,*) & y=tmp(1,*)
	d = sphdist(x(0),y(0),x(1),y(1),/deg)
	pix = 20/d
	tmp=convert_coord(ixmd,iymd,/dev,/to_data)	; Long, Lat at
	lon_cen=tmp(0,*) & lat_cen=tmp(1,*)		;   map center.
	lon_cen=lon_cen(0) & lat_cen=lat_cen(0)
 
	;=============================================================
	;	Pack up info in a structure
	;=============================================================
	out = { projection:map_proj_name(proj), projection_code:proj, $
	  lat:lat, lon:lon, ang:ang, $
	  lon_cen:lon_cen, lat_cen:lat_cen, $
	  iso_used:flag_iso, clip_used:flag_clip, $
	  noborder_used:flag_noborder, color_used:flag_col, color:clr, $
	  scale_used:flag_scale, scale:scale, central_azi_used:flag_azi, $
	  azi:azi, limit_used:flag_lim, limit:lim, $
	  stand_par_used:flag_par, parallels:std_par, $
	  satellite_par_used:flag_sat, sat_p:sat_p, $
	  ellipsoid_used:flag_ell, ellipsoid:ellip, $
	  position:pos, pix_deg:pix, idx:fix(idx), idy:fix(idy), $
	  ix:fix(ix1), iy:fix(iy1) } 
 
	;=======================================
	;	List
	;=======================================
	if keyword_set(list) then begin
	  print,' '
          print,' Values set from embedded scaling (window = '+wintxt+'):'
          print,' '
	  print,' Map projection = ',map_proj_name(proj),' ('+$
	    strtrim(proj,2)+')'
	  print,' Reference lat, long = ',strtrim(lat,2),', ',strtrim(lon,2), $
	  '    Map angle = ',strtrim(ang,2)
	  print,' Map center lat, long = ',strtrim(lat_cen,2),$
	    ', ',strtrim(lon_cen,2)
	  txt = ''
	  if flag_noborder eq 0 then begin
	    txt=' Map Border of color = '+strtrim(clr,2)
	  endif else begin
	    txt=' No Map Border'
	  endelse
	  if flag_iso eq 1 then txt=txt+',   /ISO used' else $
	    txt=txt+',   No /ISO used'
	  if flag_clip eq 1 then txt=txt+',   CLIP=0 used' else $
	    txt=txt+',   No CLIP=0'
	  print,txt
	  if max(abs(lim)) eq 0 then begin
	    print,' No LIMIT used'
	  endif else begin
	    print,' LIMIT = ',strtrim(lim,2)
	  endelse
	  if flag_scale then print,' Scale = ',strtrim(scale,2)
	  if flag_azi then print,' Central_azi = ',strtrim(azi,2)
	  if flag_par then begin
	    print,' Standard parallels = ',strtrim(stdpar,2)
	  endif
	  if flag_sat then begin
	    print,' Satellite parameters = ',strtrim(sat_p,2)
	  endif
	  ;------  Screen position  -----------
	  six1 = strtrim(fix(ix1),2)
	  six2 = strtrim(fix(ix2),2) & sidx=strtrim(fix(idx),2)
	  siy1 = strtrim(fix(iy1),2)
	  siy2 = strtrim(fix(iy2),2) & sidy=strtrim(fix(idy),2)
	  print,' Screen window: '+$
	    'ix1, ix2, iy1, iy2: '+six1+','+ six2+','+ siy1+','+siy2
	  print,'   position format: '+$
            'pos=['+six1+','+siy1+','+ six2+','+siy2+'],/dev'
          print,'   tvrd format:     a=tvrd('+$
            six1+','+siy1+','+ sidx+','+sidy+')'
          print,'   plots format:    plots,['+$
            six1+','+six2+','+six2+','+ six1+','+six1+'],['+$
            siy1+','+siy1+','+siy2+','+ siy2+','+siy1+'],/dev' 
          print,'   crop image:      a=tvrd('+$
            six1+','+siy1+','+ sidx+','+sidy+')'
          print,'                    swindow,xs='+sidx+',ys='+sidy+' & tv,a'
	  ;------  Central map scale (pixels/degree)  -------
	  print,' Map scale at center = '+strtrim(pix,2)+' pixels/degree'
	endif
 
	return
	end