;+ ;Purpose: A crib detailing how to plot 2-D slices through THEMIS ESA and SST ; particle distributions using thm_part_dist_array, thm_part_slice2d, ; and plot_part_slice2d ; ; Run "thm_ui_slice2d" on the IDL console to use for the GUI version. ; (Also part of the Analysis menu on the main THEMIS GUI) ; ; ;QUICK DESCRIPTION: Particle distributions from the ESA and SST are used ; to generate 2D slices of the particles’ three-dimensional velocity ; distributions. The slices show contoured particle data (in counts ; or other units) plotted against their Cartesian velocities in the ; given coordinates. If the time range requested covers multiple ; distributions then the data are averaged. ; ; ;TYPE: There are three methods of producing slices ; ; Geomtric (default): ; Exact method showing the bounds of all bins intersecting the slice plane. ; Values are averaged in cases of overlapping bins or when the slice plane ; is coplanar with a bin's boundary. ; ; 2D Interpolation (from thm_esa_slice2d, keyword:/TWO_D_INTERP): ; Datapoints within the specified range are projected onto the slice plane ; and linearly interpolated onto a regular grid. The range may be specified ; by theta value (defualt) or z-axis value (in km/s) with respect to the ; slice plane's coordinates (see ROTATION and COORD keywords). ; ; 3D Interpolation (keyword:/THREE_D_INTERP): ; The entire 3-dimensional distribution is linearly interpolated onto a ; regular 3d grid and a this slice is extracted. This method works best with dense data. ; Note: - Interpolation may occur across data gaps or areas with recorded zeroes ; - Can be slow if slice resolution is high (>250) ; ; ;KEYWORD INFO: ; ; COORD: Available coordinates (DSL by default): ; 'DSL','GSM','GSE' ; ; ROTATION: The rotation keyword can specify the orienation of the slice plane ; within the specified coordinates (DSL by default). Note: the orientation ; of the slice plane can be further specified with the SLICE_ORIENT keyword. ; 'BV': The x axis is parallel to B field; the bulk velocity defines the x-y plane. (DEFAULT) ; 'BE': The x axis is parallel to B field; the B x V(bulk) vector defines the x-y plane. ; 'xy': The x axis is V_x and the y axis is V_y. ; 'xz': The x axis is V_x and the y axis is V_z. ; 'yz': The x axis is V_y and the y axis is V_z. ; 'xvel': The x axis is V_x; the y axis is defined by the bulk velocity. ; 'perp': The x axis is the bulk velocity projected onto the plane normal to the B field; y is B x V(bulk) ; 'perp_xy': Geometric xy coordinates projected onto the plane normal to the B field. ; 'perp_xz': Geometric xz coordinates projected onto the plane normal to the B field. ; 'perp_yz': Geometric yz coordinates projected onto the plane normal to the B field. ; ; UNITS: Available units (DF by default): ; 'Counts', 'DF', 'Rate', 'CRate', 'Flux', 'EFlux' ; ; ;CONTAMINATION REMOVAL: See the documentation in THM_PART_DIST_ARRAY for options on ; contamination removal. The contamination/background removal keywords used with ; THM_PART_GETSPEC and THM_PART_MOMENTS may be used in calls to THM_PART_DIST_ARRAY. ; Also see THM_CRIB_SST_CONTAMINATION and THM_CRIB_ESA_BGND_REMOVE. ; ; ; ;OTHER: A general overview of velocity slices can be found the THEMIS User's guide. ; Detailed information on keywords and usage can be found in the documentation for: ; thm_part_slice2d.pro (main slices routine) ; plot_part_slice2d.pro (slices plotting routine) ; thm_part_dist_array.pro (particle data load routine) ; ; ;NOTE: A second crib detailing how to export multiple slices at once ; is at the bottom of this document. ; ;- PRO thm_crib_part_slice2d ;=========================================================== ; Load Data ;=========================================================== ; Set time range of interest to create array of particle distributions probe = 'b' day = '2008-02-26/' start_time = time_double(day + '04:52:00') end_time = time_double(day + '04:52:30') ; Pad time range to ensure enough data is loaded trange = [start_time - 90, end_time + 90] ; Load magnetic field and l2 velocity data ; *B data is required for all rotations other than geometric (xyz). ; *Velocity data is required for BV, BE, perp, and xvel. If velocity ; data is not loaded explicitly it will be automatically calculated ; from the raw distribution. thm_load_fgm, probe=probe, datatype = 'fgh', level=2, coord='dsl', trange=trange thm_load_esa, probe=probe, datatype = 'peib_velocity_dsl', trange=trange ; Create array of ESA particle distributions ; (specify options for contamination removal here) peib_arr = thm_part_dist_array(probe=probe, type='peib', trange=trange, $ ; pass in tplot names of desired mag/vel data mag_data='thb_fgh_dsl', $ vel_data='thb_peib_velocity_dsl') ; Create array of SST particle distributions ; (specify options for contamination removal here) psif_arr = thm_part_dist_array(probe=probe, type='psif', trange=trange, $ ; set this keyword to use the new SST calibrations ; outlines in thm_crib_sst_calibration.pro ; /sst_cal, $ ; pass in tplot names of desired mag/vel data mag_data='thb_fgh_dsl') ;=========================================================== ; Create Basic Slices ;=========================================================== ; Set the time slice_time = start_time ; time at which to calculate the slice timewin = 30. ; the window over which data will be averaged ; Default geometric ESA slice ; ----------------- ; -cuts along xy plane in DSL thm_part_slice2d, peib_arr, slice_time=slice_time, timewin=timewin, $ part_slice=part_slice, $ ; slice data xgrid=xgrid, $ ; x axis ygrid=ygrid, $ ; y axis slice_info=slice_info ; plotting info ; Plot plot_part_slice2d, part_slice, xgrid, ygrid, slice_info stop ; Lower resolution ESA slice ; ----------------- ; -High resolution (500) is used by default to show well defined bins. ; This example uses the same default algorithm but at lower resolution ; and with smoothing applied to the result. thm_part_slice2d, peib_arr, slice_time=slice_time, timewin=timewin, $ resolution = 51, $ ; default resolution for 2d interpolation smooth = 3, $ ; default smoothing for 2d interpolation part_slice=part_slice, $ xgrid=xgrid, $ ygrid=ygrid, $ slice_info=slice_info plot_part_slice2d, part_slice, xgrid, ygrid, slice_info stop ; ESA slice using 2d interpolation ; ----------------- ; -The 2d interpolation method is (nearly) identical to that ; of thm_esa_slice2d. ; -Many default setting such as slice resolution and smoothing ; are will change when using this method. More documentation on ; defaults is listed below in the "Options Overview" section thm_part_slice2d, peib_arr, slice_time=slice_time, timewin=timewin, $ /two_d_interp, $ ;2d interp keyword part_slice=part_slice, $ xgrid=xgrid, $ ygrid=ygrid, $ slice_info=slice_info ; Plot plot_part_slice2d, part_slice, xgrid, ygrid, slice_info stop ; ESA slice using 3d interpolation ; ----------------- ; -3d interpolation differs from 2d in that the entire ; data set is utilized. ; -see "Options Overview" section for method specific defaults thm_part_slice2d, peib_arr, slice_time=slice_time, timewin=timewin, $ /three_d_interp, $ ;3d interp keyword part_slice=part_slice, $ xgrid=xgrid, $ ygrid=ygrid, $ slice_info=slice_info plot_part_slice2d, part_slice, xgrid, ygrid, slice_info stop ; SST slice in GSM coordinates ; ----------------- ; -slices through the xy plane is GSM coordinates ; -set rotation keyword to use the GSM xz plane ; -note: rotations such as BV will be invarient between ; coordinates, such rotations will essentially override ; the COORD keyword thm_part_slice2d, psif_arr, slice_time=slice_time, timewin=timewin, $ COORD='GSM', $ ; GSM coordinates ; rotation='xz' ; slice along xz plane in GSM part_slice=part_slice, $ xgrid=xgrid, $ ygrid=ygrid, $ slice_info=slice_info plot_part_slice2d, part_slice, xgrid, ygrid, slice_info stop ; Custom aligned SST slice ; ----------------- ; -SLICE_ORIENT keyword specifies the normal for the slice plane ; in the given coordinates/rotation, this example will plot a ; slice normal to <1,0,5> in GSE coordinates ; -note: If called with the ROTATION keyword SLICE_ORIENT will ; specify a normal vector within the coordinates given by ROTATION, ; (e.g. inserting ROTATION='BV' here would yield a slice normal ; to <1,0,5> in BV coordinates) thm_part_slice2d, psif_arr, slice_time=slice_time, timewin=timewin, $ COORD='GSE', $ ; GSE coordinates slice_orient = [1,0,5.], $ ; slie plane normal part_slice=part_slice, $ xgrid=xgrid, $ ygrid=ygrid, $ slice_info=slice_info plot_part_slice2d, part_slice, xgrid, ygrid, slice_info stop ; Plot multiple distributions ; ----------------- ; -Up to four distribution arrays can be used to create a ; single slice using the geometric method. This example ; combines ESA and SST ion distributions. thm_part_slice2d, psif_arr, peib_arr, $ ; pass in up to for distributions ; as normal arguments slice_time=slice_time, timewin=timewin, $ COORD = 'GSE', $ ; GSE coords erange = [0,200000], $ ; limit energy range to exclude ; larger SST energies part_slice=part_slice, $ xgrid=xgrid, $ ygrid=ygrid, $ slice_info=slice_info ; Plot plot_part_slice2d, part_slice, xgrid, ygrid, slice_info stop ;=========================================================== ; Options Overview ;=========================================================== ; Coordinates and Rotation ;------------------------- ; ;coord = 'gsm' ;slice based off GSM coordinates (DSL is default) ;coord = 'gse' ;slice based off GSE coordinates ;rotation = 'xyz' ; no rotation, uses specified coordinates rotation = 'BV' ; the x axis is parallel to B field, the bulk velocity defines the x-y plane. ; Cuts are made along the x-y axis by default, the following ; keyword specifies the slice plane's normal and and offset ; from zero along that normal(in km/s) (Invalid for 2D Interpolation) slice_orient = [0, 0, 1] ; slice plane is the x-y axis ;slice_orient = [0, 1, 0] ; slice plane is the x-z axis ;slice_orient = [1, 0, 0, 400] ; slice plane is the y-z axis, translated 400 km/s from the origin ; Other Options ;-------------- ;erange = [0,5000] ; set this array to the desired energy limits in eV ; data outside this range will not be plotted ;units = 'eflux' ; set this keyword to change the units used ('DF' is default) ;onecount = 1 ; set to 1 to remove any bins that fall below the one-count level after averaging ;smooth = 6 ; smooth data w/ boxcar average (9 point window size) ;center_time = 1 ; set this keyword to make SLICE_TIME the center of the time window ;resolution = 200 ; set this keyword to change the resolution of the final plot ; Method Specific Defaults ;------------------------- ; By default the following options are set for each method ; unless specified otherwise. ; Geometric (default method) ;resolution = 500. ; 2D Interpolation ;resolution = 51. ;smooth = 3 ;thetaRange = [-20,20.] ; 3D Interpolation ;resolution = 150. ; Create slice with listed options ;----------------------------- thm_part_slice2d, peib_arr, slice_time=slice_time, timewin=timewin, $ ; Interpolation methods ; /two_d_interp, $ ; /three_d_interp, $ ; Optional keywords: coord=coord, rotation=rotation, slice_orient=slice_orient, $ onecount=onecount, $ smooth=smooth, $ center_time=center_time, $ ; thetarange=thetarange, $ ;only for 2d interpolation ; Output: part_slice=part_slice, $ xgrid=xgrid, $ ygrid=ygrid, $ slice_info=slice_info, $ ; This varible will contain any error messages fail=fail if keyword_set(fail) then print, fail ; Plot ;------ plot_part_slice2d, part_slice, xgrid, ygrid, slice_info stop ;=========================================================== ; Ploting Options ;=========================================================== ; Create xy ESA slice in DSL thm_part_slice2d, peib_arr, slice_time=slice_time, timewin=timewin, $ part_slice=part_slice, $ xgrid=xgrid, $ ygrid=ygrid, $ slice_info=slice_info ; Optional keywords for plotting: ; ------------------------------- ;olines = 8 ; change the number of contour line levels ; (this is primarily useful when using the interpolation options) range = [1.E-14,1.E-8] ; limit data range for plotting xrange = [1500.,-1500] ; specify x-range to plot yrange = [1500.,-1500] ; specify y-range to plot sundir = 1 ;plot sun direction (may load STATE data) plotaxes = 0 ;turn off dashed lines at x,y=0 ;ecircle = 0 ;turn off plotting of energy limits ;plotbulk = 0 ;stun off plotting bulk velocity projection on the slice plane ; Plot plot_part_slice2d, part_slice, xgrid, ygrid, slice_info, $ ; Optional keywords: range=range, $ xrange=xrange, $ yrange=xrange, $ plotbulk=plotbulk, $ ecircle=ecircle, $ olines=olines, $ sundir=sundir, $ plotaxes=plotaxes stop END ;=========================================================== ; Multiple Slices ;=========================================================== ;+ ; ; A Secondary Crib: Create Multiple Slices and Export: ; ; ;- PRO thm_crib_part_slice2d_multi ; Choose output location: outputfolder='C:\THEMIS\slice_crib\' ;Modify and uncomment this line if you're using Windows ;outputfolder='/THEMIS/slice_crib/' ;Modify and uncomment this line if you're using Mac/Unix/Linux ; Load Data ;----------- ; Set time range day = '2008-02-26/' start_time = time_double(day + '04:50:00') end_time = time_double(day + '04:55:00') ; pad time range to ensure enough data is loaded trange=[start_time - 90, end_time + 90] ; set data types probe = 'b' type = 'pseb' ; This example will use geometric coordinates so skip loading of mag & velocity data ;thm_load_fgm, probe=probe, datatype = 'fgh', level=2, coord='dsl', trange=trange ;thm_load_esa, probe=probe, datatype = 'peeb_velocity_dsl', trange=trange ; Create array of SST particle distributions dist_arr = thm_part_dist_array(probe=probe,type=type, trange=trange, $ ;contamination removal method_sunpulse_clean = 'median', $ mask_remove = .99, $ ; /sst_cal, $ ; mag_data='thb_fgh_dsl', $ ; vel_data='thb_peeb_velocity_dsl') vel_auto=1) ;Set Options ;------------ timewin = 60. ; set the time window for each slice incriment = 30. ; time incriment for next slice's start coord = 'gsm' ; GSM coordinates slice_orient = [0,1,0,0] ; slice along x-z plane ;units = 'counts' ; 'DF' by default range = [2.2e-27, 2.2e-20] ; plot using fixed range this time ; Use loop to create multiple slices and export plots ;---------------------------------------------------- slice_time = start_time while slice_time lt end_time do begin ;Create slice thm_part_slice2d, dist_arr, slice_time=slice_time, timewin=timewin, $ coord=coord, rotation=rotation, slice_orient=slice_orient, $ units=units, $ part_slice=part_slice, $ xgrid=xgrid, $ ygrid=ygrid, $ slice_info=slice_info, $ fail=fail ; the FAIL variable will contain a string message if something goes wrong if keyword_set(fail) then begin print, 'An error occured while creating the slice at '+time_string(slice_time)+':' print, fail endif ;create filename for image outputfile = outputfolder + time_string(format=2,slice_time) + '_th'+probe+'_'+type ;Call plotting procedure plot_part_slice2d, part_slice, xgrid, ygrid, slice_info, $ export=outputfile, $ ;automatically export ; /eps, $ ;set this keyword to export to postscript range=range ;constant range slice_time += incriment ;increment time endwhile END