[Matplotlib-users] Basemap/ orthographic projection plot doesn't respect globe boundary

[Will H. <wh...@le...>]

Hi forum/ mailing list,

When I plot in the orthographic projection I'm getting the large artefact
shown below extending away from the north east of the globe.

I'm not finding the same problem when plotting in a full globe projection so
I'm presuming the problem is with the way I'm projecting everything rather
than my data itself.

I've included my plotting code below, if anyone is able to spot some glaring
omissions/ errors I'd be most grateful (I've been using python/ matplotlib
for only a couple of weeks now!).


http://old.nabble.com/file/p28117655/binploterr.png 


#!/usr/local/bin/python2.6

import numpy as np

import matplotlib.pyplot as plt

from mpl_toolkits.basemap import Basemap

import sys, glob


#input must be 3 col file of lons lats and data

#bins input values into half degree grid, ignores negative values


plts = glob.glob('*.plt')

x = np.arange(-180, 180, 0.5); y = np.arange(-90, 90, 0.5)

grid_lon, grid_lat = np.meshgrid(x,y) #regularly spaced 2D grid

n_vals = np.zeros((360,720)) #mean divisor

dat = np.zeros((360,720)) #2D grid of zeros 


for pt in plts:

    

    in_file = pt

    data = np.loadtxt(in_file, comments = ';')

    fname = in_file.split('.')[0]


    lon = data[:,0] #original 1D list

    lat = data[:,1] #original 1D list

    slcol = data[:,2] #z data


    lon = (np.around(lon*2))/2 #round to nearest .0 or 0.5

    lat = (np.around(lat*2))/2 #round to nearest .0 or 0.5


    ##keep the below between files


    j=0


    for i in slcol:

        if lon[j] < 0:

            grid_lon_ind = 360+(lon[j]*2)

            grid_lat_ind = 180+(lat[j]*2)

        else:

            grid_lon_ind = 360-(lon[j]*2)

            grid_lat_ind = 180+(lat[j]*2)

            

        if i > 0:

            dat[grid_lat_ind, grid_lon_ind] += i #add i'th value

            n_vals[grid_lat_ind, grid_lon_ind] += 1 #increase cell counter
by 1 for each extra value

        j+=1


dat = np.nan_to_num(dat/n_vals)


#create map object

fig = plt.figure()

m = Basemap(projection='ortho', lon_0=lon[(len(lon)/2)], lat_0=0,
resolution='l', area_thresh=10000.)

#m = Basemap(projection='moll',lon_0=0,resolution='c', area_thresh=10000.)


X,Y = m(grid_lon, grid_lat)


#pass all 2d arrays to pcolor

im = m.pcolormesh(X,Y,dat)


#add coastlines, globe boundary and colourbar

m.drawcoastlines()

m.drawmapboundary()

m.drawparallels(np.arange(-90, 90,30))

m.drawmeridians(np.arange(-180,180,30))


fig.colorbar(im)

plt.title('CH20 and ting')

plt.savefig('binplot.png')



Thanks for your help,


Will.
-- 
View this message in context: http://old.nabble.com/Basemap--orthographic-projection-plot-doesn%27t-respect-globe-boundary-tp28117655p28117655.html
Sent from the matplotlib - users mailing list archive at Nabble.com.