;---------------------------
;IDL MAKEFILE
;---------------------------
;
;---------------------------
;Central file of used IDL analysis
;---------------------------
pro main ;has to correspond to the name specified in call command in makefile
;---------------------------
;DAS PARAMETERS INPUT
;---------------------------
dt_slow = 0.5D-5
dt_fast = 1d-6
polarity = [1d,1d,-1d,-1d,1d,-1d,1d,1d,1d,-1d,-1d,1d,1d,1d,0d,-1d] ;[-]
normalisation = [1d,1d,1d,1d,1d,2d,2d,2d,2d,2d,2d,2d,2d,2d,0d,1d] ;[-]
uloop_cal = 5.5d ;[-]
bt_cal = 70.42d ;[T/Wb]
;---------------------------
;INITIALISATION
;---------------------------
path_fast = '../../../../DAS/0311NIturbo.ON/'+['NIturbo_01','NIturbo_02','NIturbo_03','NIturbo_04','NIturbo_05','NIturbo_06','NIturbo_07','NIturbo_08','NIturbo_09','NIturbo_10','NIturbo_11','NIturbo_12','NIturbo_13','NIturbo_14','NIturbo_15','NIturbo_16'] ;location of central data with respect to folder containing this program
path_slow = '../../0411ShotHomepage.ONN/usbscopes' ;location of basic parameters data
read = size(file_reading(4,path_slow),/dimensions)
N_slow = read[1]
read = size(file_reading(2,path_fast[0]), /dimensions)
N_fast = read[1]
shot_no = long(file_reading(1, '../../../../ShotNo'))
time_slow = dindgen(N_slow)*dt_slow*1d3 ;[ms]
time_fast = dindgen(N_fast)*dt_fast*1d3 ;[ms]
theta = 360d/16d*(dindgen(15)-1d) ;[deg] poloidal angle
plasmadetect = long(file_reading(1,'../../0411ShotHomepage.ONN/Plasma'))
if plasmadetect eq 1L then begin
plasmastart = file_reading(1,'../../0411ShotHomepage.ONN/PlasmaStart')
plasmaend = file_reading(1,'../../0411ShotHomepage.ONN/PlasmaEnd')
plasmastart = plasmastart*1d3 ;[mus]
plasmaend = plasmaend*1d3 ;[mus]
time_low = plasmastart[0]
time_up = plasmaend[0]
endif else begin
time_low = 8d3
time_up = 16d3
endelse
N_fft = 400d ;[mus] length of fft window
fft_jmp_rat = 0.10d ;[-] ratio of time_len by what the fft windows are moved forward
N_fft_time = long((time_up-time_low)/N_fft/fft_jmp_rat+1d)
a_w = 2d ;[mus] coeff. in STFT for Gaussian windowing
N_vid_time = 200d
vid_jmp_rat = 0.1d
N_vid = long((time_up-time_low)/N_vid_time/vid_jmp_rat + 1)
;---------------------------
;BASIC PARAMETERS
;---------------------------
;FILE READING
data_slow = file_reading(4,path_slow)
;CALCULATION
offset_1 = total(data_slow[0,0:499])/5d2
offset_2 = total(data_slow[1,0:499])/5d2
offset_3 = total(data_slow[2,0:499])/5d2
data_slow[0,*] = uloop_cal*(data_slow[0,*]-offset_1)
data_slow[1,*] = bt_cal*total(data_slow[1,*]-offset_2,/cumulative)*dt_slow
data_slow[2,*] = 1d-3*i_plasma(data_slow[0,*],data_slow[2,*]-offset_3,dt_slow)
q_est = 90.3125d*data_slow[1,*]/data_slow[2,*]
;PLOTTING
make_plot = plot_multi(time_slow,[data_slow[0:2,*],q_est],['Uloop [V]','Bt [T]','Ip [kA]','q [-]'],'Basic parameters of shot no. '+strtrim(shot_no,2),'../outputs/plot_'+strtrim(shot_no,2)+'_basic',[-1d,-5d-2,-2d-1,0d],[15d,0.5d,6d,8d])
;---------------------------
;MHD ANALYSIS
;---------------------------
;FILE READING AND DATA PREPARATION
data_fast = dblarr(16,N_fast)
for i=0,15 do begin
data_temp = file_reading(2,path_fast[i])
data_fast[i,*] = data_temp[1,*]
endfor
for i=0,15d do begin
data_fast[i,*] = polarity[i]*normalisation[i]*data_fast[i,*]
endfor
for i=0,15 do begin
data_fast[i,*] = smooth(data_fast[i,*],5)
endfor
;FFT ANALYSIS
FFT_arr = dblarr(N_fft_time, N_fft/2, 16)
fft_time = time_low + dindgen(N_fft_time)/double(N_fft_time)*(time_up-time_low) ;[mus]
fft_freq = dindgen(N_fft/2)/double(N_fft/2)/dt_fast/2d
data_fft = data_fast
w_t = exp(-0.5d*(a_w*(dindgen(N_fft)-0.5d*double(N_fft))/(0.5d*double(N_fft)))^2) ;windowing function (Gaussian)
for i=0,15 do begin
for j=0, N_fft_time-1 do begin
data_temp = abs(fft(w_t*data_fft[i,time_low+double(j)*fft_jmp_rat*N_fft:time_low+(double(j)*fft_jmp_rat+1d)*N_fft-1], /double))
FFT_arr[j,*,i] = smooth(data_temp[0:N_fft/2-1],3)
endfor
endfor
for i=0,8 do begin
make_plot = plot_fft(fft_time, 1d-3*fft_freq[0:N_fft/7-1], fft_arr[*,0:N_fft/7-1,i],'Time [ms]','f [kHz]','Frequencies at coil no. '+strtrim(i+1,2)+' of shot no. '+strtrim(shot_no,2), '../outputs/plot_'+strtrim(shot_no,2)+'_fft_0'+strtrim(i+1,2))
endfor
for i=9,15 do begin
make_plot = plot_fft(fft_time, 1d-3*fft_freq[0:N_fft/7-1], fft_arr[*,0:N_fft/7-1,i],'Time [ms]','f [kHz]','Frequencies at coil no. '+strtrim(i+1,2)+' of shot no. '+strtrim(shot_no,2), '../outputs/plot_'+strtrim(shot_no,2)+'_fft_'+strtrim(i+1,2))
endfor
;VIDEO PLOTS
data = dblarr(15,N_fast) ;this is all due to coil 15 being off
data[0,*] = data_fast[15,*]
data[1:14,*] = data_fast[0:13,*]
if long(N_vid) le 99L then begin
for i=0,9 do begin
make_plot = plot_profile(time_fast[time_low+double(i)*vid_jmp_rat*N_vid_time:time_low+(double(i)*vid_jmp_rat+1d)*N_vid_time-1], theta, transpose(data[*,time_low+double(i)*vid_jmp_rat*N_vid_time:time_low+(double(i)*vid_jmp_rat+1d)*N_vid_time-1]),'Time [ms]','Theta [deg]','Perturbations of shot no. '+strtrim(shot_no,2), '../outputs/video/plot_pert_000'+strtrim(i,2),0.5d)
endfor
for i=10,N_vid-1 do begin
make_plot = plot_profile(time_fast[time_low+double(i)*vid_jmp_rat*N_vid_time:time_low+(double(i)*vid_jmp_rat+1d)*N_vid_time-1], theta, transpose(data[*,time_low+double(i)*vid_jmp_rat*N_vid_time:time_low+(double(i)*vid_jmp_rat+1d)*N_vid_time-1]),'Time [ms]','Theta [deg]','Perturbations of shot no. '+strtrim(shot_no,2), '../outputs/video/plot_pert_00'+strtrim(i,2),0.5d)
endfor
endif
if (long(N_vid) gt 99L) and (long(N_vid) le 999L) then begin
for i=0,9 do begin
make_plot = plot_profile(time_fast[time_low+double(i)*vid_jmp_rat*N_vid_time:time_low+(double(i)*vid_jmp_rat+1d)*N_vid_time-1], theta, transpose(data[*,time_low+double(i)*vid_jmp_rat*N_vid_time:time_low+(double(i)*vid_jmp_rat+1d)*N_vid_time-1]),'Time [ms]','Theta [deg]','Perturbations of shot no. '+strtrim(shot_no,2), '../outputs/video/plot_pert_000'+strtrim(i,2),0.5d)
endfor
for i=10,99 do begin
make_plot = plot_profile(time_fast[time_low+double(i)*vid_jmp_rat*N_vid_time:time_low+(double(i)*vid_jmp_rat+1d)*N_vid_time-1], theta, transpose(data[*,time_low+double(i)*vid_jmp_rat*N_vid_time:time_low+(double(i)*vid_jmp_rat+1d)*N_vid_time-1]),'Time [ms]','Theta [deg]','Perturbations of shot no. '+strtrim(shot_no,2), '../outputs/video/plot_pert_00'+strtrim(i,2),0.5d)
endfor
for i=100,N_vid-1 do begin
make_plot = plot_profile(time_fast[time_low+double(i)*vid_jmp_rat*N_vid_time:time_low+(double(i)*vid_jmp_rat+1d)*N_vid_time-1], theta, transpose(data[*,time_low+double(i)*vid_jmp_rat*N_vid_time:time_low+(double(i)*vid_jmp_rat+1d)*N_vid_time-1]),'Time [ms]','Theta [deg]','Perturbations of shot no. '+strtrim(shot_no,2), '../outputs/video/plot_pert_0'+strtrim(i,2),0.5d)
endfor
endif
if (long(N_vid) gt 999L) then begin
for i=0,9 do begin
make_plot = plot_profile(time_fast[time_low+double(i)*vid_jmp_rat*N_vid_time:time_low+(double(i)*vid_jmp_rat+1d)*N_vid_time-1], theta, transpose(data[*,time_low+double(i)*vid_jmp_rat*N_vid_time:time_low+(double(i)*vid_jmp_rat+1d)*N_vid_time-1]),'Time [ms]','Theta [deg]','Perturbations of shot no. '+strtrim(shot_no,2), '../outputs/video/plot_pert_000'+strtrim(i,2),0.5d)
endfor
for i=10,99 do begin
make_plot = plot_profile(time_fast[time_low+double(i)*vid_jmp_rat*N_vid_time:time_low+(double(i)*vid_jmp_rat+1d)*N_vid_time-1], theta, transpose(data[*,time_low+double(i)*vid_jmp_rat*N_vid_time:time_low+(double(i)*vid_jmp_rat+1d)*N_vid_time-1]),'Time [ms]','Theta [deg]','Perturbations of shot no. '+strtrim(shot_no,2), '../outputs/video/plot_pert_00'+strtrim(i,2),0.5d)
endfor
for i=100,999 do begin
make_plot = plot_profile(time_fast[time_low+double(i)*vid_jmp_rat*N_vid_time:time_low+(double(i)*vid_jmp_rat+1d)*N_vid_time-1], theta, transpose(data[*,time_low+double(i)*vid_jmp_rat*N_vid_time:time_low+(double(i)*vid_jmp_rat+1d)*N_vid_time-1]),'Time [ms]','Theta [deg]','Perturbations of shot no. '+strtrim(shot_no,2), '../outputs/video/plot_pert_0'+strtrim(i,2),0.5d)
endfor
for i=1000,N_vid-1 do begin
make_plot = plot_profile(time_fast[time_low+double(i)*vid_jmp_rat*N_vid_time:time_low+(double(i)*vid_jmp_rat+1d)*N_vid_time-1], theta, transpose(data[*,time_low+double(i)*vid_jmp_rat*N_vid_time:time_low+(double(i)*vid_jmp_rat+1d)*N_vid_time-1]),'Time [ms]','Theta [deg]','Perturbations of shot no. '+strtrim(shot_no,2), '../outputs/video/plot_pert_'+strtrim(i,2),0.5d)
endfor
endif
end
