绘制地图：图形化显示海地地震危机数据

可以参见

https://gist.github.com/dalguji/10979733

从这里找到的解决方法

地图

第八章258页

x, y =m(cat_data.LONGITUDE,cat_data.LATITUDE)

运行中出错

错误代码

Traceback (most recent call last):

File “C:/Python27/base.py”, line 72, in

x, y =m(cat_data.LONGITUDE,cat_data.LATITUDE)

File “C:\Python27\lib\site-packages\mpl_toolkits\basemap__init__.py”, line 1148, in call

xout,yout = self.projtran(x,y,inverse=inverse)

File “C:\Python27\lib\site-packages\mpl_toolkits\basemap\proj.py”, line 286, in call

outx,outy = self._proj4(x, y, inverse=inverse)

File “C:\Python27\lib\site-packages\mpl_toolkits\basemap\pyproj.py”, line 388, in call

_proj.Proj._fwd(self, inx, iny, radians=radians, errcheck=errcheck)

File “_proj.pyx”, line 122, in _proj.Proj._fwd (src/_proj.c:1571)

RuntimeError

错误已经修改了

把x, y =m(cat_data.LONGITUDE,cat_data.LATITUDE)

改为

x, y = m(cat_data.LONGITUDE.values, cat_data.LATITUDE.values)

# -*- coding: utf-8 -*-
"""
Created on Mon Nov 14 11:11:40 2016

@author: Administrator
"""
import pandas as pd
from pandas import DataFrame
import numpy as np
path = r'F:\python-code\python-for-data-analysis\pydata-book-master\ch08\Haiti.csv'
data = pd.read_csv(path)
#data.describe()
data = data[(data.LATITUDE>18)&(data.LATITUDE<20)&(data.LONGITUDE>-75)&(data.LONGITUDE<-70)&data.CATEGORY.notnull()]
#data2.head()
def to_cat_list(catstr):
stripped = (x.strip() for x in catstr.split(','))
return [x for x in stripped if x]
def get_all_categories(cat_series):
cat_sets = (set(to_cat_list(x)) for x in cat_series)
return sorted(set.union(*cat_sets))

def get_english(cat):
code, names = cat.split('.')
if '|' in names:
names = names.split('|')[1]
return code ,names.strip()
#get_english('2. Urgences logistiques | Vital Lines')

all_cats = get_all_categories(data.CATEGORY)

english_mapping = dict(get_english(x) for x in all_cats)
english_mapping['2a']

english_mapping['6c']
def get_code(seq):
return [x.split('.')[0] for x in seq if x]
all_codes = get_code(all_cats)
code_index = pd.Index(np.unique(all_codes))
dummy_frame = DataFrame(np.zeros((len(data),len(code_index))),index=data.index,columns=code_index)
#dummy_frame.ix[:,:6]
for row, cat in zip(data.index,data.CATEGORY):
codes=get_code(to_cat_list(cat))
dummy_frame.ix[row,codes]=1
data = data.join(dummy_frame.add_prefix('category_'))
data.ix[:,10:15]

from mpl_toolkits.basemap import Basemap

import matplotlib.pyplot as plt
def basic_haiti_map(ax=None, lllat = 17.25,urlat = 20.25,lllon = -75,urlon=-71):
m = Basemap(ax=ax,projection = 'stere',
lon_0=(urlon+lllon)/2,
lat_0 = (urlat+lllat)/2,
llcrnrlat = lllat,urcrnrlat=urlat,
llcrnrlon = lllon,urcrnrlon = urlon,
resolution = 'f')
m.drawcoastlines()
m.drawstates()
m.drawcountries()
return m

fig,axes = plt.subplots(nrows = 2,ncols =2,figsize=(12,10))
fig.subplots_adjust(hspace =0.05,wspace=0.05)
to_plot = ['2a','1','3c','7a']
lllat =17.25
urlat=20.25
lllon =-75
urlon=-71
for code,ax in zip(to_plot,axes.flat):
m = basic_haiti_map(ax,lllat=lllat,urlat=urlat,lllon=lllon,urlon=urlon)
cat_data=data[data['category_%s' % code] ==1]
#x, y =m(cat_data.LONGITUDE,cat_data.LATITUDE)
x, y = m(cat_data.LONGITUDE.values, cat_data.LATITUDE.values)
m.plot(x,y,'k.',alpha = 0.5)
ax.set_title('%s:%s' % (code, english_mapping[code]))

结果