pro ltf_auto, lun, outfil, dir, infil, dset_arr, filnamdc, $
		imin, imax, qplot, method, $
		xst, xen0, yst, yen0, nx, ny, xstep, ystep, $
		areadnoise2, aavg_off, max_gainc
;
if (lun lt 0) then begin
    ltf, '****', strver
    lun = abs(lun)
    close, lun
    openw, lun, outfil
    printf, lun, 'Program run: ', !stime
    printf, lun, 'Pixels Commanded ' + strtrim(xst,2) + ':' + strtrim(xen0,2) + $
		' , ' + strtrim(yst,2) + ':' + strtrim(yen0,2)
    printf, lun, 'Block Size: ' + strtrim(fix(nx),2) + 'x' + strtrim(fix(ny),2) + $
		'   Step Size: ' + strtrim(xstep,2) + ' & ' + strtrim(ystep,2)
    printf, lun, 'Max gain constant to have histogram generated for = ' + strtrim(max_gainc,2)
    printf, lun, strver
    printf, lun, '   '
    printf, lun, '   '
end
;
ndset = n_elements(dset_arr)
if (ndset le 1) then begin
    print, 'LTF_AUTO: Need at least two images to do LTF
    stop
end
npairs=0
for i=1,ndset-1 do npairs=npairs+i
;
gainc_lin_arr  = fltarr(npairs)
gainc_log_arr  = fltarr(npairs)
gainc_hist_arr = fltarr(npairs)
;
tit1='         Gain Constant             Fit Error       Offset Corr        Linear Fit        Hist Fit    Image
tit2='    Linear     Log     Hist   Linear  Log    Hist   Max    Avg   DC_Off Readnoise  "B"    FWHM BLS_Noise DC_Off
;
if (strtrim(filnamdc,2) ne '') then begin
    ;TODO 
end else begin
    imgdc = 0
end
;
k=0
for iimg1=0, ndset-2 do begin
	rd_sda, dir+infil, dset_arr(iimg1), index1, img1
	filnam1 = infil + ' Dset: ' + strtrim(dset_arr(iimg1),2) + ' ('+fmt_tim(index1.gen) + ')'

	for iimg2=iimg1+1, ndset-1 do begin
            rd_sda, dir+infil, dset_arr(iimg2), index2, img2
	    filnam2 = infil + ' Dset: ' + strtrim(dset_arr(iimg2),2) + ' ('+fmt_tim(index2.gen) + ')'

	    expos = mbe2sec(index2.sxt.explevmode mod 64)
	    ccd = 'CCD: TIJ44C'
	    temp = 'Temperature = ' + string(sxt_ccd_temp( index0.sxt.temp_ccd, format = '(f8.3)')
	    label = ccd + '  ' + temp

	    ltf, filnam1, img1, filnam2, img2, filnamdc, imgdc, label, $
		imin, imax, qplot, method, $
		xst, xen0, yst, yen0, nx, ny, xstep, ystep, $
		areadnoise2, aavg_off, max_gainc, $
		gainc_lin, gainc_log, gainc_hist, $
		avg_off_lin, readnoise2_lin, b_lin, fwhm, $
		bls_noise, avg_off_img, $
		err_linfit, err_logfit, err_histfit, maxoff, avgoff

	   if (k eq 0) then begin
		fmt='(10x, f8.3, 2i6, f7.3, f7.2, 3x, a, 3x, a, 3x, a)'
		printf, lun, expos, imin, imax, areadnoise2, aavg_off, infil, filnamdc, label, format=fmt
		printf, lun, tit1
		printf, lun, tit2
	   end
	   fmt='(1x, 3f9.2, ' + $				;28
		'3f7.2, 2f7.2, ' + $				;35 =  63
		'f9.2, f9.3, f7.3, f7.2, ' + $			;32 =  95
		'2f7.2, 2i3)'					;20 = 115
	   printf, lun, gainc_lin, gainc_log, gainc_hist, $
		err_linfit, err_logfit, err_histfit, maxoff, avgoff, $
		avg_off_lin, readnoise2_lin, b_lin, fwhm, $
		bls_noise, avg_off_img, iimg1, iimg2, format=fmt

	   gainc_lin_arr(k)  = gainc_lin
	   gainc_log_arr(k)  = gainc_log
	   gainc_hist_arr(k) = gainc_hist
	   k=k+1
	end
end
;
avg_lin  = total(gainc_lin_arr) / npairs
avg_log  = total(gainc_log_arr) / npairs
avg_hist = total(gainc_hist_arr) / npairs
fmt='(1x, 3f9.2)'
printf, lun, avg_lin, avg_log, avg_hist, format=fmt
dev_lin  = (max(gainc_lin_arr)-min(gainc_lin_arr)) / max(gainc_lin_arr)/2*100
dev_log  = (max(gainc_log_arr)-min(gainc_log_arr)) / max(gainc_log_arr)/2*100
dev_hist  = (max(gainc_hist_arr)-min(gainc_hist_arr)) / max(gainc_hist_arr)/2*100
!c=0
printf, lun, dev_lin, dev_log, dev_hist, format=fmt
printf, lun, '   '
;
;;lprint,'d'
;
return
end