diff --git a/docs/rfpy.rst b/docs/rfpy.rst index dd38934..53a34d4 100644 --- a/docs/rfpy.rst +++ b/docs/rfpy.rst @@ -48,7 +48,7 @@ where ``RfPy`` can be installed along with some of its dependencies. .. sourcecode:: bash - conda create -n rfpy python=3.10 obspy spectrum -c conda-forge + conda create -n rfpy "python=3.10" "setuptools=60" obspy spectrum -c conda-forge Activate the newly created environment: @@ -62,8 +62,8 @@ Install remaining dependencies using ``pip`` inside the ``rfpy`` environment: pip install stdb -Installing from GitHub master branch ------------------------------------- +Installing from GitHub development branch +----------------------------------------- .. sourcecode:: bash diff --git a/rfpy/plotting.py b/rfpy/plotting.py index 214c19f..4311152 100644 --- a/rfpy/plotting.py +++ b/rfpy/plotting.py @@ -23,8 +23,8 @@ """ Functions to plot single station P-receiver functions as wiggle plots. -Options to plot by receiver functions # vs time/depth, by back-azimuth vs -time/depth or slowness vs time. +Options to plot by receiver functions # vs time/depth, by back-azimuth vs +time/depth or slowness vs time. """ @@ -41,11 +41,11 @@ def wiggle(stream1, stream2=None, sort=None, tmin=0., tmax=30, normalize=True, save=False, title=None, form='png'): """ - Function to plot receiver function traces by index in stream. By default, + Function to plot receiver function traces by index in stream. By default, only one stream is required, which produces a Figure with a single panel. - Optionally, if a second stream is present, the Figure will contain two panels. - The stream(s) can be sorted by stats attribute ``sort``, normalized, and - the Figure can be saved as .eps. + Optionally, if a second stream is present, the Figure will contain two + panels. The stream(s) can be sorted by stats attribute ``sort``, + normalized, and the Figure can be saved as .eps. Parameters ---------- @@ -59,16 +59,16 @@ def wiggle(stream1, stream2=None, sort=None, tmin=0., tmax=30, normalize=True, xmax : float Maximum x-axis value displayed in the Figure. normalize : bool - Whether or not to normalize the traces according to the max amplitude + Whether or not to normalize the traces according to the max amplitude in ``stream1`` save : bool - Whether or not to save the Figure + Whether or not to save the Figure title : str Title of plot """ - if sort: + if sort is not None: try: stream1.traces.sort(key=lambda x: x.stats[sort], reverse=False) except: @@ -150,10 +150,10 @@ def wiggle_bins(stream1, stream2=None, tr1=None, tr2=None, norm=None, save=None, folder=None, show=True): """ Function to plot receiver function according to either baz or - slowness bins. By default, - only one stream is required, which produces a Figure with a single panel. - Optionally, if a second stream is present, the Figure will contain two panels. - If the single trace arguments are present, they will be plotted at the top. + slowness bins. By default, only one stream is required, which produces a + Figure with a single panel. Optionally, if a second stream is present, + the Figure will contain two panels. If the single trace arguments are + present, they will be plotted at the top. Parameters ---------- @@ -172,10 +172,10 @@ def wiggle_bins(stream1, stream2=None, tr1=None, tr2=None, xtyp : str Type of x-axis label (either 'time' or 'depth') norm : float - Normalization value applied to all traces. If not specified, + Normalization value applied to all traces. If not specified, default amplitude ranges will be set. save : str - Filename of figure to be saved, including extension. + Filename of figure to be saved, including extension. folder : str Folder name where figure will be saved. show : bool @@ -183,6 +183,10 @@ def wiggle_bins(stream1, stream2=None, tr1=None, tr2=None, """ + st1 = stream1.copy() + if stream2 is not None: + st2 = stream2.copy() + if not (btyp == 'baz' or btyp == 'slow' or btyp == 'dist'): raise(Exception("Type has to be 'baz' or 'slow' or 'dist'")) if not (xtyp == 'time' or xtyp == 'depth'): @@ -192,10 +196,10 @@ def wiggle_bins(stream1, stream2=None, tr1=None, tr2=None, # Figure out scaling here if norm is not None: - for tr in stream1: + for tr in st1: tr.data /= norm - if stream2: - for tr in stream2: + if stream2 is not None: + for tr in st2: tr.data /= norm if btyp == 'baz': maxval = 10. @@ -218,21 +222,21 @@ def wiggle_bins(stream1, stream2=None, tr1=None, tr2=None, maxvalT = maxval # Time axis - nn = stream1[0].stats.npts - sr = stream1[0].stats.sampling_rate + nn = st1[0].stats.npts + sr = st1[0].stats.sampling_rate time = np.arange(-nn/2, nn/2)/sr # Initialize figure fig = plt.figure() plt.clf() - if stream2 and tr1 and tr2: + if (stream2 is not None) and (tr1 is not None) and (tr2 is not None): # Get more control on subplots ax1 = fig.add_axes([0.1, 0.825, 0.3, 0.05]) ax2 = fig.add_axes([0.1, 0.1, 0.3, 0.7]) ax3 = fig.add_axes([0.45, 0.825, 0.3, 0.05]) ax4 = fig.add_axes([0.45, 0.1, 0.3, 0.7]) - elif stream2 and not tr1 and not tr2: + elif (stream2 is not None) and (tr1 is None) and (tr2 is None): ax1 = None ax2 = fig.add_subplot(121) ax3 = None @@ -258,10 +262,11 @@ def wiggle_bins(stream1, stream2=None, tr1=None, tr2=None, # facecolor='red', facecolor='k', linewidth=0) - ax1.plot(time, tr1.data, + ax1.plot( + time, tr1.data, linewidth=0.25, c='k') ylim = np.max([np.max(np.abs(tr1.data)), np.max(np.abs(tr2.data))]) - if btyp=='slow': + if btyp == 'slow': ylimT = ylim/2. else: ylimT = ylim @@ -272,7 +277,7 @@ def wiggle_bins(stream1, stream2=None, tr1=None, tr2=None, ax1.set_xlim(trange[0], trange[1]) # Plot binned SV traces in back-azimuth on bottom left - for tr in stream1: + for tr in st1: # Define y axis if btyp == 'baz': @@ -295,7 +300,8 @@ def wiggle_bins(stream1, stream2=None, tr1=None, tr2=None, facecolor='k', # facecolor='red', linewidth=0) - ax2.plot(time, y+tr.data*maxval, + ax2.plot( + time, y+tr.data*maxval, linewidth=0.25, c='k') ax2.set_xlim(trange[0], trange[1]) @@ -316,7 +322,7 @@ def wiggle_bins(stream1, stream2=None, tr1=None, tr2=None, elif xtyp == 'depth': ax2.set_xlabel('Depth (km)') ax2.grid(ls=':') - if not ax1: + if (ax1 is None): ax2.set_title('Radial') if ax3: @@ -334,7 +340,8 @@ def wiggle_bins(stream1, stream2=None, tr1=None, tr2=None, facecolor='k', # facecolor='red', linewidth=0) - ax3.plot(time, tr2.data, + ax3.plot( + time, tr2.data, linewidth=0.25, c='k') ax3.set_xlim(trange[0], trange[1]) ax3.set_ylim(-1.*ylimT, ylimT) @@ -344,7 +351,7 @@ def wiggle_bins(stream1, stream2=None, tr1=None, tr2=None, if ax4: # Plot binned SH traces in back-azimuth on bottom right - for tr in stream2: + for tr in st2: # Define y axis if btyp == 'baz': @@ -367,7 +374,8 @@ def wiggle_bins(stream1, stream2=None, tr1=None, tr2=None, # facecolor='red', facecolor='k', linewidth=0) - ax4.plot(time, y+tr.data*maxvalT, + ax4.plot( + time, y+tr.data*maxvalT, linewidth=0.25, c='k') ax4.set_xlim(trange[0], trange[1]) @@ -386,15 +394,15 @@ def wiggle_bins(stream1, stream2=None, tr1=None, tr2=None, ax4.set_xlabel('Depth (km)') ax4.set_yticklabels([]) ax4.grid(ls=':') - if not ax3: + if (ax3 is None): ax4.set_title('Transverse') if save: if folder is not None: - plt.savefig(folder + '/' + stream1[0].stats.station + + plt.savefig(folder + '/' + st1[0].stats.station + '.' + save, format=save.split('.')[-1], dpi=300) else: - plt.savefig('RF_PLOTS/' + stream1[0].stats.station + + plt.savefig('RF_PLOTS/' + st1[0].stats.station + '.' + save, format=save.split('.')[-1], dpi=300) if show: @@ -412,14 +420,16 @@ def wiggle_single_event(rfdata, filt=None, pre_filt=None, trange=None): taxis = np.arange(-nn/2., nn/2.)/sr if pre_filt: - lqtcopy.filter('bandpass', freqmin=pre_filt[0], + lqtcopy.filter( + 'bandpass', freqmin=pre_filt[0], freqmax=pre_filt[1], corners=2, zerophase=True) if filt: - rfcopy.filter('bandpass', freqmin=filt[0], + rfcopy.filter( + 'bandpass', freqmin=filt[0], freqmax=filt[1], corners=2, zerophase=True) - fig, (ax1, ax2, ax3, ax4) = plt.subplots(4,1, figsize=(7,5)) + fig, (ax1, ax2, ax3, ax4) = plt.subplots(4, 1, figsize=(7, 5)) ax1.plot(taxis, lqtcopy[0], label=lqtcopy[0].stats.channel, lw=1) ax2.plot(taxis, lqtcopy[1], label=lqtcopy[0].stats.channel, lw=1) @@ -442,7 +452,7 @@ def wiggle_single_event(rfdata, filt=None, pre_filt=None, trange=None): ax4.legend() plt.suptitle( 'AZ corr: {0:.1f}; BAZ: {1:.1f}\n SNR: {2:.1f}; CC: {3:.1f}'.format( - rfdata.sta.azcorr, rfdata.meta.baz, + rfdata.sta.azcorr, rfdata.meta.baz, rfdata.meta.snr, rfdata.meta.cc)) plt.show() diff --git a/rfpy/scripts/rfpy_plot.py b/rfpy/scripts/rfpy_plot.py index 5937c20..310d782 100644 --- a/rfpy/scripts/rfpy_plot.py +++ b/rfpy/scripts/rfpy_plot.py @@ -37,9 +37,7 @@ def get_plot_arguments(argv=None): """ - Get Options from :class:`~optparse.OptionParser` objects. - - This function is used for data processing on-the-fly (requires web connection) + Get arguments from :class:`~argparse.ArgumentParser` objects. """ @@ -212,7 +210,8 @@ def get_plot_arguments(argv=None): type=str, dest="trange", help="Specify two floats that define the time range (in sec.) for " + - "the x-axis on the RF figure. Negative times are allowed [Default 0., 30.]") + "the x-axis on the RF figure. Negative times are allowed " + + "[Default 0., 30.]") PlotGroup.add_argument( "--save-fig", action="store", @@ -220,17 +219,18 @@ def get_plot_arguments(argv=None): type=str, default=None, help="Specify figure filename if you wish to save the figure. By " + - "default, the station name will be pre-appended to the file name and " + - "saved to 'RF_PLOTS' unless --save-rfs is set. Valid figure formats " + - "are 'png', 'jpg', 'eps', 'pdf'. [Default does not save figure]") + "default, the station name will be pre-appended to the file name " + + "and saved to 'RF_PLOTS' unless --save-rfs is set. Valid figure " + + "formats are 'png', 'jpg', 'eps', 'pdf'. [Default does not save " + + "figure]") PlotGroup.add_argument( "--save-rfs", action="store", dest="rf_folder", type=str, default=None, - help="Specify folder name to save the plotted RFs. Lower case characters " + - "will be capitalized. [Default does not save RFs]") + help="Specify folder name to save the plotted RFs. Lower case " + + "characters will be capitalized. [Default does not save RFs]") PlotGroup.add_argument( "--hide-fig", action="store_true", @@ -238,14 +238,6 @@ def get_plot_arguments(argv=None): default=False, help="Specify if you do not wish to show the figure upon " + "execution. [Default shows the figure]") - # PlotGroup.add_argument( - # "--plot-event-dist", - # action="store_true", - # dest="plot_event_dist", - # default=False, - # help="Plot distribution of events on map. Other Plotting Options " + - # "will be applied to this figure (title, save, etc.). " + - # "[Default no plot]") args = parser.parse_args(argv) @@ -308,14 +300,12 @@ def get_plot_arguments(argv=None): parser.error( "Error: --trange should contain 2 " + "comma-separated floats") - if args.trange[0] < args.trange[1]: - parser.error( - "Error: --trange=arg1,arg2 requires arg1 < arg2") if args.figname is not None: if args.figname.split('.')[-1] not in ['png', 'jpg', 'eps', 'pdf']: parser.error( - "Error: Invalid figure extension. Choose among 'jpg', 'png', 'eps', 'pdf'") + "Error: Invalid figure extension. Choose among 'jpg', " + + "'png', 'eps', 'pdf'") if args.rf_folder is not None: args.rf_folder = args.rf_folder.upper() @@ -592,45 +582,54 @@ def main(): # Save RFs if args.rf_folder is not None: - import json + import json import csv - - # Write the command-line arguments to a hidden JSON file - with open(Path(args.rf_folder) / '.proc.json', 'w') as file: - json.dump(args.__dict__, file) + + # Write the command-line arguments to a hidden JSON file + with open(Path(args.rf_folder) / '.proc.json', 'w') as file: + json.dump(args.__dict__, file) # Write baz and slow to .csv file bb = [tr.stats.baz for tr in rf_tmp[0]] ss = [tr.stats.slow for tr in rf_tmp[0]] with open( Path(args.rf_folder) / 'baz_slow.csv', - mode='w', - newline='') as file: + mode='w', newline='') as file: writer = csv.writer(file) writer.writerow(['# BAZ (deg)', 'SLOW (s/mk)']) for b, s in zip(bb, ss): writer.writerow([f"{b:.0f}", f"{s:.3f}"]) - # Write RFs as SAC files + # Write RFs as SAC files with header if tr1 is not None: - tr1.write(args.rf_folder + '/' + (tr1.stats.station + - '.' + tr1.stats.channel + '.stack.sac'), format='SAC') + tr1.write( + args.rf_folder + '/' + tr1.stats.station + + '.' + tr1.stats.channel + '.stack.sac', + format='SAC') if tr2 is not None: - tr2.write(args.rf_folder + '/' + tr2.stats.station + - '.' + tr2.stats.channel + '.stack.sac', format='SAC') - for tr in rf_tmp[0]: - tr.write(args.rf_folder + '/' + tr.stats.station + - '.' + tr.stats.channel + '.' + f"{tr.stats.baz:.0f}" + - '.' + f"{tr.stats.slow:.3f}".split('.')[-1] + '.sac', format='SAC') - for tr in rf_tmp[1]: - tr.write(args.rf_folder + '/' + tr.stats.station + - '.' + tr.stats.channel + '.' + f"{tr.stats.baz:.0f}" + - '.' + f"{tr.stats.slow:.3f}".split('.')[-1] + '.sac', format='SAC') - - # # Event distribution - # if args.plot_event_dist: - # plotting.event_dist( - # rfRstream, phase=args.phase, save=args.figname) + tr2.write( + args.rf_folder + '/' + tr2.stats.station + + '.' + tr2.stats.channel + '.stack.sac', + format='SAC') + for trR, trT in zip(rf_tmp[0], rf_tmp[1]): + baz = trR.stats.baz + slow = trR.stats.slow + sachdr = { + 'baz': baz, + 'user0': slow, + } + trR.stats.sac = sachdr + trT.stats.sac = sachdr + filestr = '.' + f"{baz:.0f}".zfill(3) + '.' + \ + f"{slow:.3f}".split('.')[-1] + '.sac' + trR.write( + args.rf_folder + '/' + trR.stats.station + + '.' + trR.stats.channel + filestr, + format='SAC') + trT.write( + args.rf_folder + '/' + trT.stats.station + + '.' + trT.stats.channel + filestr, + format='SAC') # Update processed folders procfold.append(stfld)