;+ ; Generate spice product file for a day. ; ; time. A time in UT sec. Only time[0] is used to determine the day. ; probe=. A string of 'a' or 'b'. ; filename=. A string to set the output file name. ; ; To replace rbsp_gen_spice_product. Need to regenerate from 2012-09-25 to 2018-09-18. ;- pro rbsp_read_spice_gen_file, time, probe=probe, $ filename=file, errmsg=errmsg ;---Check inputs. if n_elements(file) eq 0 then begin errmsg = handle_error('No output file ...') return endif if n_elements(probe) eq 0 then begin errmsg = handle_error('No input probe ...') return endif if n_elements(time) eq 0 then begin errmsg = handle_error('No input time ...') return endif ;---Constants and settings. secofday = 86400d deg = 180d/!dpi rad = !dpi/180 re = 6378d & re1 = 1d/re errmsg = '' common_time_step = 1d ; sec. ; Derived settings. date = time[0]-(time[0] mod secofday) time_range = date+[0,secofday] prefix = 'rbsp'+probe+'_' ;---Init file. out_dir = fgetpath(file) if file_test(out_dir,/directory) eq 0 then file_mkdir, out_dir data_file = file if file_test(data_file) eq 1 then file_delete, data_file ; overwrite old files. cdf_id = cdf_create(data_file) ginfo = dictionary($ 'TITLE', 'RBSP position and boom direction from SPICE kernel', $ 'TEXT', 'Generated by Sheng Tian at the University of Minnesota' ) cdf_save_setting, ginfo, filename=cdf_id ;---Save common_times. ; Need to overlapping with next day for interpolation purpose. common_times = make_bins(time_range, common_time_step) utname = 'Epoch' data = stoepoch(common_times,'unix') settings = dictionary($ 'FIELDNAM', 'Epoch', $ 'UNITS', 'ms', $ 'VAR_TYPE', 'support_data' ) cdf_save_var, utname, value=data, filename=cdf_id, cdf_type='CDF_EPOCH' cdf_save_setting, settings, var=utname, filename=cdf_id ;---Load spice kernels for all times. defsysv,'!rbsp_spice', exists=flag if flag eq 0 then rbsp_load_spice_kernel ; Prepare epochs used for q_uvw2gse. scid = strupcase(prefix+'science') cspice_str2et, time_string(common_times[0], tformat='YYYY-MM-DDThh:mm:ss.ffffff'), epoch0 epochs = epoch0+common_times-common_times[0] ; Get r_coord, v_coord, and q_uvw2coord. ;foreach coord, ['gsm','gse'] do begin foreach coord, ['gse'] do begin cap_coord = strupcase(coord) suffix = '_'+coord rbsp_load_spice_state, probe=probe, coord=cap_coord, times=common_times, /no_spice_load ; r_coord. vname = prefix+'r'+suffix tplot_rename, prefix+'state_pos_'+coord, vname data = get_var_data(vname)*re1 store_data, vname, common_times, data settings = dictionary($ 'FIELDNAM', 'R '+cap_coord, $ 'UNITS', 'Re', $ 'VAR_TYPE', 'data', $ 'DEPEND_0', utname ) cdf_save_var, vname, value=data, filename=cdf_id cdf_save_setting, settings, var=vname, filename=cdf_id ; v_coord. vname = prefix+'v'+suffix tplot_rename, prefix+'state_vel_'+coord, vname ; r_coord = get_var_data(prefix+'r'+suffix)*re ; data = r_coord ; for ii=0,2 do data[*,ii] = deriv(common_times, r_coord[*,ii]) data = get_var_data(vname) settings = dictionary($ 'FIELDNAM', 'V '+cap_coord, $ 'UNITS', 'km/s', $ 'VAR_TYPE', 'data', $ 'DEPEND_0', utname ) cdf_save_var, vname, value=data, filename=cdf_id cdf_save_setting, settings, var=vname, filename=cdf_id ; q_uvw2coord. cspice_pxform, scid, cap_coord, epochs, muvw muvw = transpose(muvw) ; [n,3,3]. data = mtoq(muvw) vname = prefix+'q_uvw2'+coord settings = dictionary($ 'FIELDNAM', 'Q UVW2'+cap_coord, $ 'UNITS', '#', $ 'VAR_TYPE', 'data', $ 'DEPEND_0', utname ) cdf_save_var, vname, value=data, filename=cdf_id cdf_save_setting, settings, var=vname, filename=cdf_id ; Save to use for DSC. store_data, prefix+'wsc_'+coord, common_times, muvw[*,*,2] endforeach ; Wsc_gse. vname = prefix+'wsc_gse' data = get_var_data(vname) settings = dictionary($ 'FIELDNAM', 'Wsc in GSE', $ 'UNITS', '#', $ 'VAR_TYPE', 'data', $ 'DEPEND_0', utname ) cdf_save_var, vname, value=data, filename=cdf_id cdf_save_setting, settings, var=vname, filename=cdf_id ; MLT. vname = prefix+'mlt' r_gse = get_var_data(prefix+'r_gse') r_sm = cotran(r_gse, common_times, 'gse2sm') data = atan(r_sm[*,1],r_sm[*,0])*deg/15+12 settings = dictionary($ 'FIELDNAM', 'MLT', $ 'UNITS', 'hr', $ 'VAR_TYPE', 'data', $ 'yrange', [0,24], $ 'DEPEND_0', utname ) cdf_save_var, vname, value=data, filename=cdf_id cdf_save_setting, settings, var=vname, filename=cdf_id ; MLat. vname = prefix+'mlat' mlat = atan(r_sm[*,2],sqrt(r_sm[*,0]^2+r_sm[*,1]^2)) data = mlat*deg settings = dictionary($ 'FIELDNAM', 'MLat', $ 'UNITS', 'deg', $ 'VAR_TYPE', 'data', $ 'DEPEND_0', utname ) cdf_save_var, vname, value=data, filename=cdf_id cdf_save_setting, settings, var=vname, filename=cdf_id ; lshell. vname = prefix+'lshell' dis = snorm(r_gse) data = dis/(cos(mlat)^2) settings = dictionary($ 'FIELDNAM', 'L', $ 'UNITS', '#', $ 'VAR_TYPE', 'data', $ 'DEPEND_0', utname ) cdf_save_var, vname, value=data, filename=cdf_id cdf_save_setting, settings, var=vname, filename=cdf_id ; UVW to DSC, adopted from rbsp_uvw_to_dsc. ; Save sun_phase_[ssha,sshb]. This is the essential info for this rotation. ; From sun_phase (sphase) to rotation, check Line 198-204 in rbsp_uvw_to_dsc.pro. wsc_gse = get_var_data(prefix+'wsc_gse') ncommon_time = n_elements(common_times) ndim = 3 ; The normal to the sun sensor triggering plane, nTsg, is defined as the ; cross product between X_GSE and the spin axis direction. x_gse = [1d,0,0] ## (dblarr(ncommon_time)+1) ntsg = sunitvec(vec_cross(x_gse, wsc_gse)) ; The triggering sun sensor look direction perpendicular to the spin axis ; is defined by the cross product of the spin axis direction and nTsg tsg = sunitvec(vec_cross(wsc_gse,ntsg)) ;foreach the_ssh, ['A','B'] do begin foreach the_ssh, ['A'] do begin ssh_id = 'RBSP'+strupcase(probe)+'_SSH_'+the_ssh cspice_pxform, 'GSE', ssh_id, epochs, gse2ssh tsg_ssh = dblarr(ncommon_time,3) for ii=0,ncommon_time-1 do tsg_ssh[ii,*] = gse2ssh[*,*,ii] ## tsg[ii,*] deg = constant('deg') sphase = atan(tsg_ssh[*,1],tsg_ssh[*,0])*deg index = where(sphase lt 0, count) if count ne 0 then sphase[index] += 360. vname = prefix+'sphase_ssh'+strlowcase(the_ssh) data = sphase settings = dictionary($ 'FIELDNAM', 'sphase of SSH'+the_ssh, $ 'UNITS', 'deg', $ 'VAR_TYPE', 'data', $ 'DEPEND_0', utname ) cdf_save_var, vname, value=data, filename=cdf_id cdf_save_setting, settings, var=vname, filename=cdf_id endforeach ; Spin period and spin phase. timespan, time, 86400d, /seconds rbsp_load_state, probe=probe, dt=common_time_step, /no_spice_load vname = prefix+'spin_period' data = get_var_data(prefix+'spinper') settings = dictionary($ 'FIELDNAM', 'Spin period', $ 'UNITS', 'sec', $ 'VAR_TYPE', 'data', $ 'DEPEND_0', utname ) cdf_save_var, vname, value=data, filename=cdf_id cdf_save_setting, settings, var=vname, filename=cdf_id vname = prefix+'spin_phase' data = get_var_data(prefix+'spinphase') settings = dictionary($ 'FIELDNAM', 'Spin phase', $ 'UNITS', 'deg', $ 'VAR_TYPE', 'data', $ 'DEPEND_0', utname ) cdf_save_var, vname, value=data, filename=cdf_id cdf_save_setting, settings, var=vname, filename=cdf_id cdf_close, cdf_id end ; Last dates of RBSP-A and -B time = time_double(['2019-10-15']) probe = 'a' time = time_double(['2019-07-17']) probe = 'b' ;timespan, time, 86400d, /seconds ;rbsp_read_emfisis, time+[0,86400], probe=probe, id='l2%magnetometer' ;rbsp_uvw_to_dsc, probe, prefix+'b_uvw', /debug ;stop time = time_double(['2013-07-17']) probe = 'a' prefix = 'rbsp'+probe+'_' file = join_path([homedir(),prefix+'spice_products_'+time_string(time,tformat='YYYY_MMDD')+'_v08.cdf']) rbsp_read_spice_gen_file, time, probe=probe, $ filename=file, errmsg=errmsg end