pro training_20180328_advanced ; ***************** ; settings for the training ; ***************** erg_init, remote_data_dir='https://ergsc.isee.nagoya-u.ac.jp/data/ergsc_training/nagoya_201803/' uname = '***************' pass = '****************' ; ***************** ; set time span ; ***************** timespan, '2017-03-28', 1, /day ; ***************** ; load orbit CDF & set var label ; ***************** set_erg_var_label ; ***************** ; load OFA-MATRIX ; ***************** erg_load_pwe_ofa, datatype='matrix', uname=uname, pass=pass pr_matrix = 'erg_pwe_ofa_matrix_l2_' ; ***************** ; load MGF L2 CDF ; ***************** erg_load_mgf, datatype='64hz', coord='sgi', uname=uname, pass=pass erg_load_mgf, datatype='8sec', uname=uname, pass=pass ; ***************** ; analyze OFA-MATRIX ; ***************** get_data, pr_matrix+'Bx_Bx_132', data=s00 get_data, pr_matrix+'Bx_By_132', data=s01 get_data, pr_matrix+'Bx_Bz_132', data=s02 get_data, pr_matrix+'By_Bx_132', data=s10 get_data, pr_matrix+'By_By_132', data=s11 get_data, pr_matrix+'By_Bz_132', data=s12 get_data, pr_matrix+'Bz_Bx_132', data=s20 get_data, pr_matrix+'Bz_By_132', data=s21 get_data, pr_matrix+'Bz_Bz_132', data=s22 rr=dblarr(3,3,n_elements(s00.x),n_elements(s00.v),2) rr[0,0,*,*,*]=s00.y & rr[1,0,*,*,*]=s01.y & rr[2,0,*,*,*]=s02.y rr[0,1,*,*,*]=s10.y & rr[1,1,*,*,*]=s11.y & rr[2,1,*,*,*]=s12.y rr[0,2,*,*,*]=s20.y & rr[1,2,*,*,*]=s21.y & rr[2,2,*,*,*]=s22.y ; ***************** ; analyze MGF data ; ***************** split_vec, 'erg_mgf_l2_mag_64hz_sgi' tinterpol_mxn, 'erg_mgf_l2_mag_64hz_sgi_x', pr_matrix+'Bx_Bx_132' tinterpol_mxn, 'erg_mgf_l2_mag_64hz_sgi_y', pr_matrix+'Bx_Bx_132' tinterpol_mxn, 'erg_mgf_l2_mag_64hz_sgi_z', pr_matrix+'Bx_Bx_132' get_data, 'erg_mgf_l2_mag_64hz_sgi_x_interp', data=data_x, dlim=dlim_x, lim=lim_x get_data, 'erg_mgf_l2_mag_64hz_sgi_y_interp', data=data_y, dlim=dlim_y, lim=lim_y get_data, 'erg_mgf_l2_mag_64hz_sgi_z_interp', data=data_z, dlim=dlim_z, lim=lim_z rotmat=dblarr(3,3,n_elements(data_x.x)) rotmat_t=dblarr(3,3,n_elements(data_x.x)) ; rotation matrix for i=0, n_elements(data_x.x)-1 do begin bvec=[data_x.y[i],data_y.y[i],data_z.y[i]] zz=[0.,0.,1.] yhat=crossp(zz,bvec) xhat=crossp(yhat,bvec) zhat=bvec xhat=xhat/sqrt(xhat[0]^2+xhat[1]^2+xhat[2]^2) yhat=yhat/sqrt(yhat[0]^2+yhat[1]^2+yhat[2]^2) zhat=zhat/sqrt(zhat[0]^2+zhat[1]^2+zhat[2]^2) rotmat[*,*,i]=([[xhat],[yhat],[zhat]]) rotmat_t[*,*,i]=transpose([[xhat],[yhat],[zhat]]) endfor ; MGF rotation data_rot=dblarr(n_elements(data_x.x), 3) for i=0, n_elements(data_x.x)-1 do begin data=[[data_x.y[i]], [data_y.y[i]],[data_z.y[i]]] data_rot[i,*] = rotmat[*,*,i] ## data endfor store_data, 'erg_mgf_l2_mag_64hz_sgi_rot_x', data={x:data_x.x, y:data_rot[*,0]}, dlim=dlim_x, lim=lim_x store_data, 'erg_mgf_l2_mag_64hz_sgi_rot_y', data={x:data_y.x, y:data_rot[*,1]}, dlim=dlim_y, lim=lim_y store_data, 'erg_mgf_l2_mag_64hz_sgi_rot_z', data={x:data_z.x, y:data_rot[*,2]}, dlim=dlim_z, lim=lim_z ylim, 'erg_mgf_l2_mag_64hz_sgi_rot_?', -4E4, 4E4, 0 tplot, [ 'erg_mgf_l2_mag_64hz_sgi_' + ['x', 'y', 'z'], 'erg_mgf_l2_mag_64hz_sgi_rot_' + ['x', 'y', 'z'] ] stop ; OFA-MATRIX rotation for i=0, n_elements(s00.x)-1 do begin for j=0, n_elements(s00.v)-1 do begin for k=0, 1 do begin rr[*,*,i,j,k] = (rotmat[*,*,i] ## rr[*,*,i,j,k]) ## rotmat_t[*,*,i] endfor endfor endfor ; ************************************ ; auto-spectra ; ************************************ get_data, pr_matrix + 'Bx_Bx_132',data=data, dlim=dlim, lim=lim store_data, pr_matrix + 'Bx_Bx_re', data={x:data.x,y:data.y[*,*,0],v:data.v}, dlim=dlim, lim=lim store_data, pr_matrix + 'Bx_Bx_im', data={x:data.x,y:data.y[*,*,1],v:data.v}, dlim=dlim, lim=lim ; get_data, pr_matrix + 'By_By_132',data=data, dlim=dlim, lim=lim store_data, pr_matrix + 'By_By_re', data={x:data.x,y:data.y[*,*,0],v:data.v}, dlim=dlim, lim=lim store_data, pr_matrix + 'By_By_im', data={x:data.x,y:data.y[*,*,1],v:data.v}, dlim=dlim, lim=lim ; get_data, pr_matrix + 'Bz_Bz_132',data=data, dlim=dlim, lim=lim store_data, pr_matrix + 'Bz_Bz_re', data={x:data.x,y:data.y[*,*,0],v:data.v}, dlim=dlim, lim=lim store_data, pr_matrix + 'Bz_Bz_im', data={x:data.x,y:data.y[*,*,1],v:data.v}, dlim=dlim, lim=lim options, pr_matrix + 'Bx_Bx_re', ytitle='Re(BxBx*)' & options, pr_matrix + 'Bx_Bx_im', ytitle='Im(BxBx*)' options, pr_matrix + 'By_By_re', ytitle='Re(ByBy*)' & options, pr_matrix + 'By_By_im', ytitle='Im(ByBy*)' options, pr_matrix + 'Bz_Bz_re', ytitle='Re(BzBz*)' & options, pr_matrix + 'Bz_Bz_im', ytitle='Im(BzBz*)' options, pr_matrix + 'B?_B?_??', ysubtitle = 'Frequency [kHz]' options, pr_matrix + 'B?_B?_??', ztitle = '[pT!U2!N/Hz]' ylim, pr_matrix + 'B?_B?_??', 0.064, 20, 1 ; kHz zlim, pr_matrix + 'B?_B?_??', 1E-3,1E2, 1 ; nT tplot, pr_matrix + ['Bx_Bx_re', 'Bx_Bx_im', 'By_By_re', 'By_By_im', 'Bz_Bz_re', 'Bz_Bz_im'] stop ; ************************************ ; wave normal angle ; ************************************ wna=dblarr(n_elements(s00.x),n_elements(s00.v)) for i=0, n_elements(s00.x)-1 do begin for j=0, n_elements(s00.v)-1 do begin ; wave normal wna_x = rr[2,1,i,j,1] / sqrt(rr[2,1,i,j,1]^2 + rr[0,2,i,j,1]^2 + rr[1,0,i,j,1]^2) wna_y = rr[0,2,i,j,1] / sqrt(rr[2,1,i,j,1]^2 + rr[0,2,i,j,1]^2 + rr[1,0,i,j,1]^2) wna_z = rr[1,0,i,j,1] / sqrt(rr[2,1,i,j,1]^2 + rr[0,2,i,j,1]^2 + rr[1,0,i,j,1]^2) wna[i,j] = abs(atan(sqrt(wna_x^2+wna_y^2)/wna_z)/!dtor) endfor endfor store_data, 'kvec',data={x:s00.x, y:wna, v:s00.v} options, 'kvec', ytitle='wave normal angle', ztitle='[degree]', ysubtitle='Frequency [kHz]', spec = 1 ylim, 'kvec', 0.064, 20, 1 ; kHz zlim, 'kvec', 0., 90, 0 ; degree tplot,[pr_matrix+'Bx_Bx_re', pr_matrix+'By_By_re', pr_matrix+'Bz_Bz_re', 'kvec'] stop ; ************************************ ; polarization ; ************************************ polarization=dblarr(n_elements(s00.x),n_elements(s00.v)) for i=0, n_elements(s00.x)-1 do begin for j=0, n_elements(s00.v)-1 do begin ; polarization polarization[i,j] = (2 * rr[1,0,i,j,1]) / (rr[0,0,i,j,0] + rr[1,1,i,j,0]) endfor endfor store_data, 'polarization',data={x:s00.x, y:polarization, v:s00.v} options, 'polarization', ytitle='polarization', ysubtitle='Frequency [kHz]', spec = 1 ylim, 'polarization', 0.064, 20, 1 ; kHz zlim, 'polarization', -1., 1., 0 ; tplot,[pr_matrix+'Bx_Bx_re', pr_matrix+'By_By_re', pr_matrix+'Bz_Bz_re', 'polarization'] stop ; ************************************ ; mask ; ************************************ get_data, pr_matrix+'Btotal_132', data=data_ref ; kvec get_data, 'kvec', data=data, dlim=dlim, lim=lim data.y[where(data_ref.y LT 1E-2)] = 'NaN' store_data, 'kvec_mask', data={x:data.x, y:data.y, v:data.v}, dlim=dlim, lim=lim ; polarization get_data, 'polarization', data=data, dlim=dlim, lim=lim data.y[where(data_ref.y LT 1E-2)] = 'NaN' store_data, 'polarization_mask', data={x:data.x, y:data.y, v:data.v}, dlim=dlim, lim=lim tplot,[pr_matrix+'Bx_Bx_re', pr_matrix+'By_By_re', pr_matrix+'Bz_Bz_re', 'kvec_mask', 'polarization_mask'] stop ; ************************************ ; overplot fce ; ************************************ ; calculate fce get_data, 'erg_mgf_l2_magt_8sec', data=data fce = data.y/ 10^(9.) * 1.6 * 10^(-19.) / (9.1093D * 10^(-31.)) / 2. / !pi / 1000. fce_half = fce / 2. store_data, 'fce', data={x:data.x, y:fce} store_data, 'fce_half', data={x:data.x, y:fce_half} options, 'fce', colors=fsc_color('yellow'), thick=2 options, 'fce_half', colors=fsc_color('magenta'), thick=2 store_data, pr_matrix+'Btotal_132_gyro', data=[pr_matrix+'Btotal_132', 'fce', 'fce_half'] store_data, 'kvec_mask_gyro', data=['kvec_mask', 'fce', 'fce_half'] store_data, 'polarization_mask_gyro', data=['polarization_mask', 'fce', 'fce_half'] ylim, '*_gyro', 0.064, 20, 1 ; kHz zlim, pr_matrix+'Btotal_132_gyro', 1E-2, 1E2, 1 ; pT^2/Hz options, pr_matrix+'Btotal_132_gyro', ytitle='B total', ysubtitle='Frequency [kHz]', ztitle='[pT!U2!N/Hz]' tplot, [pr_matrix+'Btotal_132', 'kvec_mask', 'polarization_mask']+'_gyro' stop ; ***************** ; load OFA-COMPLEX ; ***************** erg_load_pwe_ofa, datatype='complex', uname=uname, pass=pass pr_complex = 'erg_pwe_ofa_complex_l2_' get_data, pr_complex+'Ex_132', data=data_ex, dlim=dlim, lim=lim c_ex = dcomplex(data_ex.y[*,*,0]*1E-3, data_ex.y[*,*,1]*1E-3) get_data, pr_complex+'Ey_132', data=data, dlim=dlim, lim=lim c_ey = dcomplex(data.y[*,*,0]*1E-3, data.y[*,*,1]*1E-3) get_data, pr_complex+'Bx_132', data=data_bx, dlim=dlim, lim=lim c_bx = dcomplex(data_bx.y[*,*,0]*1E-12, data_bx.y[*,*,1]*1E-12) get_data, pr_complex+'By_132', data=data, dlim=dlim, lim=lim c_by = dcomplex(data.y[*,*,0]*1E-12, data.y[*,*,1]*1E-12) get_data, pr_complex+'Bz_132', data=data, dlim=dlim, lim=lim c_bz = dcomplex(data.y[*,*,0]*1E-12, data.y[*,*,1]*1E-12) ; ***************** ; calc Ez ; ***************** c_ez = dcindgen(n_elements(data_ex.x), n_elements(data_ex.v)) idx = nn(data_bx.x, data_ex.x) for i=0, n_elements(data_ex.x)-1 do begin c_ez[i,*] = (-c_ex[i,*]*c_bx[idx[i],*]-c_ey[i,*]*c_by[idx[i],*]) / c_bz[idx[i],*] endfor Sx=dindgen(n_elements(data_bx.x),n_elements(data_bx.v)) Sy=dindgen(n_elements(data_bx.x),n_elements(data_bx.v)) Sz=dindgen(n_elements(data_bx.x),n_elements(data_bx.v)) ; calc Poynting flux idx = nn(data_ex.x, data_bx.x) for i=0, n_elements(data_bx.x)-1 do begin for j=0, n_elements(data_bx.v)-1 do begin Sx[i,j]= double(c_ey[idx[i],j]*conj(c_bz[i,j])-c_ez[idx[i],j]*conj(c_by[i,j])) Sy[i,j]=-double(c_ex[idx[i],j]*conj(c_bz[i,j])-c_ez[idx[i],j]*conj(c_bx[i,j])) Sz[i,j]= double(c_ex[idx[i],j]*conj(c_by[i,j])-c_ey[idx[i],j]*conj(c_bx[i,j])) endfor endfor ; ***************** ; analyze MGF data ; ***************** split_vec, 'erg_mgf_l2_mag_64hz_sgi' ; interpolation tinterpol_mxn, 'erg_mgf_l2_mag_64hz_sgi_x',pr_complex+'Bx_132',newname='erg_mgf_l2_mag_64hz_sgi_x_interp' tinterpol_mxn, 'erg_mgf_l2_mag_64hz_sgi_y',pr_complex+'Bx_132',newname='erg_mgf_l2_mag_64hz_sgi_y_interp' tinterpol_mxn, 'erg_mgf_l2_mag_64hz_sgi_z',pr_complex+'Bx_132',newname='erg_mgf_l2_mag_64hz_sgi_z_interp' get_data,'erg_mgf_l2_mag_64hz_sgi_x_interp',data=data_x get_data,'erg_mgf_l2_mag_64hz_sgi_y_interp',data=data_y get_data,'erg_mgf_l2_mag_64hz_sgi_z_interp',data=data_z ; rotation matrix rotmat=dblarr(3,3,n_elements(data_x.x)) rotmat_t=dblarr(3,3,n_elements(data_x.x)) for i=0, n_elements(data_x.x)-1 do begin bvec=[data_x.y[i],data_y.y[i],data_z.y[i]] zz=[0.,0.,1.] yhat=crossp(zz,bvec) xhat=crossp(yhat,bvec) zhat=bvec xhat=xhat/sqrt(xhat[0]^2+xhat[1]^2+xhat[2]^2) yhat=yhat/sqrt(yhat[0]^2+yhat[1]^2+yhat[2]^2) zhat=zhat/sqrt(zhat[0]^2+zhat[1]^2+zhat[2]^2) rotmat[*,*,i]=([[xhat],[yhat],[zhat]]) rotmat_t[*,*,i]=transpose([[xhat],[yhat],[zhat]]) endfor ; ***************** ; Poynting vector calculation ; ***************** S=dindgen(n_elements(data_bx.x),n_elements(data_bx.v),3) for i=0, n_elements(data_bx.x)-1 do begin for j=0, n_elements(data_bx.v)-1 do begin S[i,j,*] = rotmat[*,*,i] ## [Sx[i,j],Sy[i,j],Sz[i,j]] endfor endfor theta=acos(S[*,*,2]/sqrt(S[*,*,0]^2+S[*,*,1]^2+S[*,*,2]^2))/!dtor store_data, 'S', data={x:data_bx.x, y:theta, v:data_bx.v}, dlim=dlim, lim=lim ylim, [pr_complex+'Etotal_132', pr_complex+'Btotal_132', 'S'], 0.064, 20, 1 zlim, pr_complex+'Etotal_132', 1E-7, 1E0, 1 ; mV^2/m^2/Hz zlim, pr_complex+'Btotal_132', 1E-2, 1E2, 1 ; pT^2/Hz zlim, 'S', 0., 180., 0 ; degree options, [pr_complex+'Etotal_132', pr_complex+'Btotal_132', 'S'], ysubtitle='Frequency [kHz]' options, pr_complex+'Etotal_132', ytitle='E total', ztitle='[mV!U2!N/m!U2!N/Hz]' options, pr_complex+'Btotal_132', ytitle='B total', ztitle='[pT!U2!N/Hz]' options, 'S', ytitle='Poynting vector', ztitle='[degree]' tplot, [pr_complex+'Etotal_132', pr_complex+'Btotal_132', 'S'] stop ; ************************************ ; mask ; ************************************ get_data, pr_complex+'Btotal_132', data=data_ref ; S get_data, 'S', data=data, dlim=dlim, lim=lim data.y[where(data_ref.y LT 1E-2)] = 'NaN' store_data, 'S_mask', data={x:data.x, y:data.y, v:data.v}, dlim=dlim, lim=lim tplot, [pr_complex+'Etotal_132', pr_complex+'Btotal_132', 'S_mask'] stop ; ************************************ ; overplot fce ; ************************************ store_data, pr_complex+'Etotal_132_gyro', $ data=[pr_complex+'Etotal_132', 'fce', 'fce_half'] store_data, pr_complex+'Btotal_132_gyro', data=[pr_complex+'Btotal_132', 'fce', 'fce_half'] store_data, 'S_mask_gyro', data=['S_mask', 'fce', 'fce_half'] ylim, '*_gyro', 0.064, 20, 1 ; kHz zlim, pr_complex+'Etotal_132_gyro', 1E-7, 1E0, 1 ; mV^2/m^2/Hz zlim, pr_complex+'Btotal_132_gyro', 1E-2, 1E2, 1 ; pT^2/Hz options, pr_complex+'Etotal_132_gyro', ytitle='E total', ysubtitle='Frequency [kHz]', ztitle='[mV!U2!N/m!U2!N/Hz]' options, pr_complex+'Btotal_132_gyro', ytitle='B total', ysubtitle='Frequency [kHz]', ztitle='[pT!U2!N/Hz]' tplot, [pr_complex+'Etotal_132', pr_complex+'Btotal_132', 'kvec_mask', 'polarization_mask', 'S_mask'] + '_gyro' stop end