分类： 学习

import numpy as np

​import matplotlib.pyplot as plt 

 from mpl_toolkits.basemap import Basemap 

 from obspy import read_events 

from obspy.imaging.beachball import beach 

 event = read_events( 'https://earthquake.usgs.gov/archive/product/moment-tensor/' 'us_20005ysu_mww/us/1470868224040/quakeml.xml', format='QUAKEML')[0] 

origin = event.preferred_origin() or event.origins[0] 

ocmec = event.preferred_focal_mechanism() or event.focal_mechanisms[0] 

tensor = focmec.moment_tensor.tensor 

import time
a='1987-06-05 12:31:28'
tm=time.strptime(a,'%Y-%m-%d %H:%M:%S')
epo_time=time.mktime(tm)

import matplotlib.pyplot as plt
import matplotlib.dates as dates

ax.xaxis.set_minor_locator(dates.YearLocator())
ax.xaxis.set_major_locator(dates.YearLocator(5))

可以自己写个简单函数，也可以使用math模块下的modf函数
import math
a,b=math.modf(32.005)
print a,b

a=0.005000000000002558
b=32.0

*****注意a，b都是浮点型

已经有了pandas数据框格式的地震目录，需要获取日期的datetime数值，并存储到新的一列dt中，可以这么操作
cata['dt']=cata.apply(lambda x: datetime(x['yr'],x['mn'],x['dy']),axis=1)
这个是按行操作的，很慢。
或者用pandas的to_datetime方法：
pd.to_datetime(df['year'].astype(str) + '-'+ df['month'].astype(str)+ '-1')

而批量操作不行，比如
datetime(x['yr'],x['mn'],x['dy'])
因为datetime不支持pandas的批量操作。

https://docs.scipy.org/doc/numpy/reference/generated/numpy.c_.html



Examples

>>> np.c_[np.array([1,2,3]), np.array([4,5,6])]

array([[1, 4], [2, 5], [3, 6]])



​>>> np.c_[np.array([[1,2,3]]), 0, 0, np.array([[4,5,6]])]

​array([[1, 2, 3, 0, 0, 4, 5, 6]])

来自于这里：https://www.zhihu.com/question/25369412

比如这个项目：https://github.com/rveciana/geoexamples/tree/master/python/shaded_relief

将tree/master改成trunk，然后用svn命令

即

svn checkout https://github.com/rveciana/geoexamples/trunk/python/shaded_relief​

就可以了。​

pandas的date_range函数可以方便的生成日期序列

from datetime import datetime​

import pandas as pd

t_series=pd.date_range(start=datetime(2000,1,1),end=datetime(2000,12,31),freq='D')

这里生成从start到end按天的日期序列，如果freq是M，则是每个自然月的月底​

需要注意的是，现在的t_series并不是datetime类型的序列，而是​timestamp类型的，因此在使用的时候还要注意转换，采用map和匿名函数可以进行方便的转换

t_series_new=t_series.map(lambda t:t.date())​

此时的t_series_new就是date类型的了，还需要注意datetime和date类型是不可以进行比较的，要转成一致才可以。​

import matplotlib.pyplot as plt

import matplotlib.dates as mdates

import datetime



years=mdates.YearLocator()

months=mdates.MonthLocator()

yearsFmt=mdates.DateFormatter('%Y')

mag_min=np.zeros(len(mag_array))

fig,ax=plt.subplots(nrows=2,ncols=1,figsize=(6,8))

ax[0].vlines(date_array,mag_min,mag_array,lw=0.1) ### 绘制m-t图

ax[0].set_xlim(datetime.datetime(2009,1,1),datetime.datetime(2018,2,28))

ax[0].set_ylim(0,np.ceil(max(mag_array)))

ax[0].xaxis.set_major_locator(years)  ### x坐标主标记刻度

ax[0].xaxis.set_major_formatter(yearsFmt)  ### x坐标主标记格式

ax[0].xaxis.set_minor_locator(months) ### x坐标次刻度

ax[0].set_xlabel('time')

ax[0].set_ylabel('magnitude')

​

ax[1].bar(date_markers,date_freq,width=35) ### 绘制n-t图，用bar绘制柱状图，这里data_markers是datetime格式的月份，因此width要30天左右

ax[1].set_xlim(datetime.datetime(2009,1,1),datetime.datetime(2018,2,28))

#ax[1].set_ylim(0,np.ceil(max(date_freq)))

ax[1].set_xlabel('time')

ax[1].set_ylabel('monthly freq')

plt.show()