pro rdmerlinfits,directory,merdat

;  procedure input the inversion results from MERLIN, save
;  in a common area, including UT

;  INPUTS:	directory	= AZAM data directory
;				merdat	= path name to the Merlin output fits file

;  common area for outputs
common com_merlinresults, $
	fld, psi, azm, eta0, dmp, dop,$
	bzero, b1mu, vmac, cen1, cen2, alpha, $
	delta, cct, chisq, chisq_i, chisq_q, $
	chisq_u, chisq_v, signed_v, niter, occt, $
	rconv, netpol, netpol_full, $
	fld_er, psi_er, azm_er, eta0_er, dmp_er, dop_er, $
	bzero_er, b1mu_er, vmac_er, cen1_er, cen2_er, $
	alpha_er, delta_er,xx,yy,slat,slong, $
	yr,mnth,dy,hr,minit,sec



;  MERLIN data
strr = merlin_fits_read(directory,merdat)


;  setup arrays
fld = strr.field_strength
psi = strr.field_inclination
azm = strr.field_azimuth
eta0 = strr.line_strength
dmp = strr.damping
dop = strr.doppler_width
bzero = strr.source_function
b1mu = strr.source_function_gradient
vmac = strr.macro_turbulence
cen1 = strr.doppler_shift1
cen2 = strr.doppler_shift2
alpha = strr.stray_light_fill_factor
delta = strr.stray_light_shift
cct = strr.continuum_intensity
occt = strr.original_continuum_intensity
chisq = strr.chisq_total
chisq_i = strr.chisq_i
chisq_q = strr.chisq_q
chisq_u = strr.chisq_u
chisq_v = strr.chisq_v
signed_v = strr.stokesv_magnitude
niter = strr.number_of_iterations
rconv = strr.fitting_attribute
netpol = strr.polarization
;mer_netpol_full = dblarr(nmer,nsl)
fld_er = strr.field_strength_error
psi_er = strr.field_inclination_error
azm_er = strr.field_azimuth_error
eta0_er = strr.line_strength_error
dmp_er = strr.damping_error
dop_er = strr.doppler_width_error
bzero_er = strr.source_function_error
b1mu_er = strr.source_function_gradient_error
vmac_er = strr.macro_turbulence_error
cen1_er = strr.doppler_shift1_error
cen2_er = strr.doppler_shift2_error
alpha_er = strr.stray_light_fill_factor_error
delta_er = strr.stray_light_shift_error
xx = strr.x_coordinate
yy = strr.y_coordinate

;  get some dimensions
npix = n_elements(xx)
sz = size(xx)
nxx = sz(1) & nyy = sz(2)

;  get the time information from the FITS file
timez = strr.times.times
yr=fix(strmid(timez,0,4))
mnth=fix(strmid(timez,5,2))
dy=fix(strmid(timez,8,2))
hr=fix(strmid(timez,11,2))
minit=fix(strmid(timez,14,2))
sec=float(strmid(timez,17,6))


;  find the solar latitude,longitude from XX, YY
mdpt = n_elements(hr)/2
timmm = timez(mdpt)
;  get a 2-D string array for time to enable more precise
;  ephemeris calculations throughout maps
time2d = strarr(nxx,nyy)
for jj = 0,nyy-1 do time2d(*,jj) = timez
;  xy2lonlat calculates its own ephemeris values based on times
arcvec = fltarr(2,npix)
arcvec(0,*) = xx
arcvec(1,*) = yy
lonlat = xy2lonlat(arcvec,time2d)
slong = fltarr(nxx,nyy) & slat = fltarr(nxx,nyy)
slong(*,*) = lonlat(0,*)
slat(*,*) = lonlat(1,*)

return
end