make_n2_files.pro (Kronos IDL library)

;+
; Contains the make_n2_files procedure
;
; :Author:
;    Baptiste Cecconi
;
; :History:
;    2006/02/22: Created
;
;    2006/02/22: Last Edit
;-
;
;+
; fabrication et sauvegarde des fichiers N2 a partir des fichiers K
; du (jour,heure) (aaaajjjd,hd) a (aaaajjjf,hf)
; a partir d'une liste de fichiers 'liste.scratch'
; donnees "analyse" seulement
;
; :Uses:
;    auto_db, cal1, cal, delpath, freq_abc, lit_hfr_1
;
; :Params:
;    aaaajjjd: in, required, type=string
;       date de debut
;    hd:   in, required, type=int
;       heure de debut
;    aaaajjjf: in, required, type=string
;       date de fin
;    hf:       in, required, type=int
;       heure de fin
;-
pro MAKE_N2_FILES, aaaajjjd, hd, aaaajjjf, hf, $
                   verbose=verbose, nodb=nodb

common KRONOSCAL, a1h1,a1h2,ai,dbcal_abc,dbcal_h,dbcal

tABC = [ [[273.,264,242],[364,365.5,356],[399,395,378],[896/2,892/2,864/2]], $
         [[546,530,504],[728,733,737],[798,792,782],[1690/2,1680/2,1660/2]], $
         [[1089,1056,1008],[1456,1464.5,1474],[1596,1580,1560],[3280/2,3260/2,3220/2]], $
         [[2180,2110,2003],[2905,2925,2945],[3192,3160,3120],[6420/2,6400/2,6340/2]] ]
cABC = [ 0.22, 0.09, 0.026 ]
tH1_2 = [[20.1,23.1,24.6,26.2,50.4/2.],[34.7,41.9,41.45,47.5,90.2/2.], $
		[65.5,78.2,84.1,91.2,171.7/2.],[125.9,153.6,163.3,175.7,336.7/2.]]
tH2_2 = [[12.3,13.9,14.75,15.2,29.6/2.],[20.3,23.4,25,26.2,50.3/2.], $
		[35.2,41.4,44.3,47.3,90.4/2.],[66,78,84.5,89.5,172/2.]]
cH = 0.19

print,'% Loading List of Files for the selected period ...'

ptrDataFileList = make_file_list(aaaajjjd, hd, aaaajjjf, hf, level='k', /noemptyfiles)
if (ptr_valid(ptrDataFileList) eq 0) then begin
	nDataFile=0
	if keyword_set(verbose) then print,"# Warning ! No K data file selected..."
endif else begin
	dataFileList=*ptrDataFileList
	nDataFile = n_elements(dataFileList)
	n2files=strarr(nDatafile)
	ptr_free, ptrDataFileList
endelse

for iDataFile = 0,nDataFile-1 do begin
	fichier = dataFileList(iDataFile)
	str_tmp = strsplit(fichier,'k/K',/extract,/regex)
	fileOut= str_tmp(0)+'n2/P'+str_tmp(1)

	print & print,'% fichier : ',fichier
	LIT_HFR_1, fichier, z
	
	if n_elements(z) gt 0 then begin
		data = {data_N2}
		data = replicate(data,n_elements(z))
		data.t97 = AJ_T97((z.an+1900L)*1000L + z.jour) + $
		  reform([3600.d0,60.d0,1.d0,1.d-2]#z.time)/86400.d0
		dt = z.msec*10
		w = where(dt eq 120) & if w(0) ne -1 then dt(w) = 125
		
		nz = n_elements(z) & xf = fltarr(nz) & xdt = xf & xdf = xf
		nfilt = 16 & freq = FREQ_ABC(nfilt) 
		for i = 0L,nz-1L do begin
			if z(i).bande le 2 then begin
				if z(i).filters ne nfilt then begin
					nfilt = z(i).filters & freq = FREQ_ABC(nfilt)
				endif
				xf(i) = freq(z(i).bande*nfilt + z(i).filter)
				xdf(i) = xf(i) *1.6/nfilt
				m = long(alog(dt(i)/125.)/alog(2.)+0.5)
				j = long(alog(nfilt/8.)/alog(2.)+0.5)
				if z(i).df eq 0 or z(i).agc1 eq 255B or z(i).agc2 eq 255B then k=0
				if (z(i).agc1 ne 255B and z(i).agc2 ne 255B) and $
				  (z(i).auto1 ne 255B or z(i).cross1 ne -999) then k=1
				if (z(i).agc1 ne 255B and z(i).agc2 ne 255B) and $
				  (z(i).auto1 ne 255B and z(i).cross1 ne -999) then k=2
				if z(i).df gt 10 then k=3
				xdt(i) = tABC(j,k,m)/3. * cABC(z(i).bande) / xf(i) * $
				  total(freq(z(i).bande*nfilt:z(i).bande*nfilt+nfilt-1)) / nfilt
			endif
			
			if z(i).bande eq 3 then begin
				xf(i) = z(i).channel+25.*(z(i).filter-z(i).filters/2.+0.5)/z(i).filters
				xdf(i) = 25./z(i).filters
				m = long(alog(dt(i)/20.)/alog(2.)+0.5)
				if z(i).df eq 0 or z(i).agc1 eq 255B or z(i).agc2 eq 255B then k=0
				if (z(i).agc1 ne 255B and z(i).agc2 ne 255B) and $
				  (z(i).auto1 eq 255B or z(i).cross1 eq -999) then k=1
				if (z(i).agc1 ne 255B and z(i).agc2 ne 255B) and $
				  (z(i).auto1 ne 255B or z(i).cross1 ne -999) then k=2
				if (z(i).agc1 ne 255B and z(i).agc2 ne 255B) and $
				  (z(i).auto1 ne 255B and z(i).cross1 ne -999) then k=3
				if z(i).df gt 10 then k=4
				xdt(i) = tH1_2(k,m) * cH
			endif
			
			if z(i).bande eq 4 then begin
				xf(i) = z(i).channel+25.*(z(i).filter-z(i).filters/2.+0.5)/z(i).filters
				xdf(i) = 25./z(i).filters
				m = long(alog(dt(i)/10.)/alog(2.)+0.5)
				if z(i).df eq 0 or z(i).agc1 eq 255B or z(i).agc2 eq 255B then k=0
				if (z(i).agc1 ne 255B and z(i).agc2 ne 255B) and $
				  (z(i).auto1 eq 255B or z(i).cross1 eq -999) then k=1
				if (z(i).agc1 ne 255B and z(i).agc2 ne 255B) and $
				  (z(i).auto1 ne 255B or z(i).cross1 ne -999) then k=2
				if (z(i).agc1 ne 255B and z(i).agc2 ne 255B) and $
				  (z(i).auto1 ne 255B and z(i).cross1 ne -999) then k=3
				if z(i).df gt 10 then k=4
				xdt(i) = tH2_2(k,m) * cH
			endif
			
		endfor
		data.f = xf
		data.dt = xdt
		data.df = xdf
		
		ant = z.df mod 10
		w = where(ant ge 1 and ant le 3) & if w(0) ne -1 then ant(w) = ant(w)-1
		data.autoX = CAL1(z.agc1,z.auto1,z.bande,ant,z.channel,z.filters,z.filter)
		w = where(data.autoX ne 0.)
		if w(0) ne -1 then data(w).autoX = 10.^( data(w).autoX /10.)
		dbcalX = fltarr(nz)
		w = where(z.auto1 eq 255b and z.cross1 ne -999)
		if w(0) ne -1 then dbcalX(w) = 10.^(dbcal(w)/10.)
		data.autoZ = $
		  CAL1(z.agc2,z.auto2,z.bande,bytarr(nz)+3b,z.channel,z.filters,z.filter)
		w = where(data.autoZ ne 0.)
		if w(0) ne -1 then data(w).autoZ = 10.^( data(w).autoZ /10.)
		dbcalZ = fltarr(nz)
		w = where(z.auto2 eq 255b and z.cross1 ne -999)
		if w(0) ne -1 then dbcalZ(w) = 10.^(dbcal(w)/10.)
		
		ax = 10.^(AUTO_DB(z.auto1)/10.)
		w = where(dbcalX ne 0.) & if w(0) ne -1 then ax(w) = dbcalX(w)
		az = 10.^(AUTO_DB(z.auto2)/10.)
		w = where(dbcalZ ne 0.) & if w(0) ne -1 then az(w) = dbcalZ(w)
		w = where(z.cross1 ne -999)
		if w(0) ne -1 then begin
			nw = n_elements(w)
			cr = float(z(w).cross1) 
			ci = float(z(w).cross2)
			wr = where(cr ne 0)
			if wr(0) ne -1 then cr(wr) = cr(wr)/abs(cr(wr))*(10.^(AUTO_DB(byte(abs(cr(wr))))/10.))
			wi = where(ci ne 0) 
			if wi(0) ne -1 then ci(wi) = ci(wi)/abs(ci(wi))*(10.^(AUTO_DB(byte(abs(ci(wi))))/10.))
			cr = cr/sqrt(ax(w)*az(w))
			ci = ci/sqrt(ax(w)*az(w))
			cm = sqrt(cr^2+ci^2)
			thph = atan(ci,cr)*!radeg
			restore, getenv('ROOT_RPWS') + '/bin/phases.xdr'
			ph = fltarr(nw)  &  att = intarr(nw) ; hyp. = sans attenuation
			for i=0L,nw-1L do begin
				if z(w(i)).bande le 2 then $
				  ph(i) = phases_abc(z(w(i)).filter,ant(w(i)),z(w(i)).bande, $
									 fix(alog(z(w(i)).filters*1.)/alog(2.)+0.5)-3)
				if z(w(i)).bande eq 3 then begin
					freq = z(w(i)).channel < 4300
					ph(i) = ph_h1(ant(w(i)),fix(freq/25),0) $
					  + phases_h(z(w(i)).filter,fix(alog(z(w(i)).filters*1.)/alog(2.)+0.5))
				endif
				if z(w(i)).bande eq 4 then begin
					freq = z(w(i)).channel < 16075
					ph(i) = ph_h2(ant(w(i)),fix((freq-25)/50),0) + $
					  phases_h(z(w(i)).filter,fix(alog(z(w(i)).filters*1.)/alog(2.)+0.5))
				endif
			endfor
			th = (thph-ph)/!radeg
			cr = cm*cos(th)
			ci = cm*sin(th)
			data(w).crossR = cr
			data(w).crossI = ci
		endif 
		
		data.ant = z.df

		data.ydh = long(strmid(delpath(dataFileList(iDataFile)),1,7))*100l+long(strmid(delpath(dataFileList(iDataFile)),9,2))
	data.num = lindgen(n_elements(data))
		print,'writing ',fileOut,' ...'
		WRITE_DATA_BINARY, fileOut, data
		
	endif

	n2files[iDataFile] = fileOut 
endfor

; Update the Kronos Database with the new files
if ~keyword_set(nodb) and nDatafile ne 0 then begin
	upsert_args = getenv('ROOT_RPWS') + '/pro/kronosdb/upsert.sh'
	if keyword_set(verbose) then upsert_args = [upsert_args, '-v']
	upsert_args = [upsert_args, '-f', n2files]
	spawn, /NOSHELL, upsert_args
end

return
end