@C:\Users\robgillies\SuperDARN\code\time.pro
@C:\Users\robgillies\SuperDARN\code\rprm.pro
@C:\Users\robgillies\SuperDARN\code\oldfit.pro
@C:\Users\robgillies\SuperDARN\code\genlib.pro
@C:\Users\robgillies\SuperDARN\code\fit.pro
@C:\Users\robgillies\SuperDARN\code\datamap.pro
@C:\Users\robgillies\SuperDARN\code\fit_reduce.pro
@C:\Users\robgillies\SuperDARN\code\fit_clipped.pro
;@C:\Users\robgillies\SuperDARN\code\read_fit_data.pro
;---------------------------------------------------------------------------------------------------------------
;READ_FIT_DATA
;
;Written by Pasha PONOMARENKO 26 January 2011.
;
;This is the universal reading routine to extract fitted data 
;from both old (.fit) and new (.fitacf) types of files.
;
;It automatically detects the data type from the file extension.
;
;This routine can read data for all beams or for a single beam.
;(see BEAM keyword). 
;
;
;----------------------------------------------------------------------------------------------------------------


pro read_fit_data, fname=fname, prm_=prm_, fit_=fit_, beam=beam, simple_stat=simple_stat, max_range=max_range, plt=plt

;Keywords:
;FNAME - filename
;PRM_ - parameter structure
;FIT_ - fitted data structure
;BEAM - required beam number (0-15); if BEAM is set to the negative number then the data from all available beams are read
;MAX_RANGE - maximum number of range gates. If this keyword is not set, then the data will be read from all used range gates
;PLT -- setting this to 1 will result in plotting a range-time velocity map in a separate window.
;SIMPLE_STAT - this keyword allows to save computer resources by reading only basic data variables (see FIT_CLIPPED.PRO fro details).
;
;If no filename supplied then open the search window
if ~keyword_set(fname) then fname=dialog_pickfile(filter='*.fit*', get_path=g_path, $
path='/data/fitred/', /read)

;Fixing the problem with beam 0 -- NULL value fro a keyword means no keyword set!
;Without this line the program would crash complaining about absence of the BEAM variable.
if ~keyword_set(beam) then beam=0;  



;Creating output array with n_rec records 
;n_rec=long(60000)
n_rec=long(120000) ;mod by R. Gillies 2011/06/15
;(this should be enough even for a 24-hr high-resolution mode file)


;Detecting the file type from the filename extension (old -- 'fit', new -- 'acf')
type=STRMID(fname, 2 , 3, /REVERSE_OFFSET)   ;;;;; change this to (fname, 2,3,..)
;type=STRMID(fname, 5 , 3, /REVERSE_OFFSET)
;print, 'type=', type
; Applying respective reading procedures for old (fit) and new (fitacf) data formats.
;To dsitinguish between these two we use the last three characters 
;from the filename extension as determined in the previous line 


;Creating empty data structures to fill.
;This is a much faster process then reading the data record by record and appending to the
;initial array becaus concatenation (i.e. c=[c,a]) is a very slow procedure.
;Arrays for the new dta format:
if (type eq 'acf' or type eq 'ACF') then begin
fitfp=FitOpen(fname,/read);
;fitfp=FitOpen(fname);
		dumb=fitRead(fitfp, prm, fit)   
		if ~keyword_set(max_range) then max_range=prm[0].nrang      
;Creating a "clipped" version of the data structure (saving only required variables, see FIT_CLIPPED)
		fit_clipped, fit=fit1, max_range=max_range, simple_stat=simple_stat    
		fit_=replicate(fit1, n_rec)
		prm_=replicate(prm, n_rec)
close, fitfp
endif

;Arrays for the old data format
if (type eq 'fit' or type eq 'fit') then begin
fitfp=oldFitOpen(fname);
		dumb=oldfitRead(fitfp, prm, fit)
		if ~keyword_set(max_range) then max_range=prm[0].nrang
		fit_clipped, fit=fit1, max_range=max_range, simple_stat=simple_stat
		fit_=replicate(fit1, n_rec)
		prm_=replicate(prm, n_rec)
s=oldfitClose(fitfp)
endif

;Record counter
count=long(0)


;Reading the data from the file record by record and filling the data strucutres:
;New data format:
if (type eq 'acf' or type eq 'ACF') then begin
fitfp=fitOpen(fname,/read)
;fitfp=fitOpen(fname)
;Reading and filling the arrays until there is no more records
while fitRead(fitfp,prm,fit) ne -1 do begin
	if beam lt 0  or prm.bmnum eq beam then begin
			struct_assign, fit, fit1
			fit_(count)=fit1
			prm_(count)=prm
			count=count+1
	endif
endwhile
close, fitfp
endif

;Old data format
if (type eq 'fit' or type eq 'fit') then begin
fitfp=oldfitOpen(fname)
while oldfitRead(fitfp,prm,fit) ne -1 do begin
	if beam lt 0  or prm.bmnum eq beam then begin
			struct_assign, fit, fit1
			fit_(count)=fit1
			prm_(count)=prm
			count=count+1
	endif
endwhile


s=oldFitClose(fitfp)
endif

;Removing remaining empty record fields from the data structures 
if count gt 0 then begin
	fit_=fit_(0:count-1)
	prm_=prm_(0:count-1)
endif else begin
;If there is no data read then we just kill the whole arrays
;This helps to save the computer's memory
	fit_=fit_[0]
	prm_=prm_[0]
endelse

;Plotting the range-time velocity field:
if keyword_set(plt) then begin
  
  v=transpose(fit_.v)
  sdarn_title,prm=prm_, title=title, t_jul=t_jul
	if beam lt 0 then beam_title=' , all beams' else beam_title=' , beam '+string(beam, '(i2)')
	polygone_plot, transpose(fit_.v), t_jul, findgen(max_range), ctab=4, zrange=[-1000,1000], $
	julian=1, xtit='UT', tit=title+beam_title, ytit='Range gate', unit='m/s'
endif

close, /all
end