1

Let me introduce a bit of context. I have a binary file where I have to extract data. This is a huge matrix which depth of the pixel for each longitude/latitude. The first part of the code is a little code in order to read and store each data (in a lon/lat interval predefined). However, I would like to draw these data on a map of the world and I'm blocked. I can plot data on a blank figure with grids but I don't know how to insert a map of the world (using Basemap).

# -*- coding: utf-8 -*-
import numpy as np
from matplotlib import pyplot as plt
import matplotlib.pyplot as plt
from mpl_toolkits.basemap import Basemap
from matplotlib.colors import LinearSegmentedColormap

# Definition of the grid
GLDB_res = 1./120.
GLDB_lon = np.linspace(-180.+0.5*GLDB_res, 180.-0.5*GLDB_res, 43200)
GLDB_lat = np.linspace(  90.-0.5*GLDB_res, -90.+0.5*GLDB_res, 21600)

j_Asia = pl.find(pl.logical_and(GLDB_lon>=0.,GLDB_lon<=40.))
i_Asia= pl.find(pl.logical_and(GLDB_lat>=50.,GLDB_lat<=90.))

# Read depth

GLDB_depth_file = 'fileURL/GlobalLakeDepth.dat'
GLDB_depth = np.fromfile(GLDB_depth_file,'<i2')
GLDB_depth_r = GLDB_depth.reshape((21600,43200))

GLDB_depth_Asia = GLDB_depth_r[i_Asia[0]:i_Asia[-1]+1,j_Asia[0]:j_Asia[-1]+1].astype(float)
GLDB_depth_Asia[np.where(GLDB_depth_Asia==0.)] = np.nan

# Création de la Basemap

map=Basemap(projection='cyl', llcrnrlat=50,urcrnrlat=90,\
        llcrnrlon=0,urcrnrlon=180, resolution='l')
map.drawcoastlines()
map.fillcontinents(color='#cc9955')

#Colormap
cmap_depth  = LinearSegmentedColormap.from_list("my_colormap", ((200./255,200/255.,1.), (20./255,20./255,165./255)), N=255, gamma=1.0)

#Graphique

plt.close()
fig=plt.figure(1,facecolor='w',figsize=(11,7))
ax=plt.subplot(1,2,1)
x,y=map(j_Asia,i_Asia) #associe un x et un y au lon et lat 
plt.scatter(x,y, c='my_colormap') #trace les données lacs avec Basemap

c=plt.colorbar(orientation='horizontal')
for i in np.arange(0,len(i_Asia),120/2): #res:1/120th deg
map=plt.plot(np.array([0,len(i_Asia)]),np.array([i,i]),'-',c=(150./255,150./255,150./255))
map=plt.plot(np.array([i,i]),np.array([0,len(i_Asia)]),'-',c=(150./255,150./255,150./255))

plt.show()

In summary, I want to insert my array Global_Depth_Asia which represents all the depth on a map of the region. The problem is that I have a problem of invalid rgb argument!

Thank you

  • 1
    Related/helpful? https://stackoverflow.com/questions/29590365/scatter-plot-data-does-not-appear-on-continents-in-hammer-basemap – DavidG Mar 26 '18 at 15:48
  • could you please elaborate on your problem part? – Aditya Mar 26 '18 at 16:01
  • In my case the depth have a specific colormap which is defined before the creation of the plot ("my_colormap"). However it seems that I have a conflict with the color between Basemap and the color linked with my array (in the plt.scatter) – Thib Ödlaniug Mar 27 '18 at 07:13

1 Answers1

1

I obviously cannot run the code to test, but there are some strange things in there.

  • plt.close() closes your basemap plot. Why would you do that?
  • plt.figure creates a new figure without the basemap plot. This is probably undesired. Also remove the plt.subplot if you want to use the basemap plot to plot to.
  • Use map.scatter() to plot the data to the map.
  • Provide the data according to which the color should be chosen to the scatter's c argument
  • Provide the colormap to the cmap argument.

    m.scatter(x,y, c=data, cmap=cmap_depth)
    
ImportanceOfBeingErnest
  • 321,279
  • 53
  • 665
  • 712
  • Thank you for your question -For the `plt.close()`, It erases the older figure when I create a new one. -I thought it would be important to create another figure in order to import and plot my array. – Thib Ödlaniug Mar 27 '18 at 07:28
  • Yes, you can use `plt.close()` **before** `Basemap()` to close previous plots if needed. But if you use if afterwards it closes the Basemap plot you want to show. – ImportanceOfBeingErnest Mar 27 '18 at 07:38
  • it seems obvious... thanks. Actually I still have an error with my colors, it seems that my colormap is empty ("length of rgba sequence should be either 3 or 4) – Thib Ödlaniug Mar 27 '18 at 07:44
  • The colormap seems to be created correctly, it consists of 255 colors. Your problem might be `c=data` where I do not know what `data` is. – ImportanceOfBeingErnest Mar 27 '18 at 07:49
  • Actually I had a problem with my data, it was not the right shape. Now I can plot it. The only remaining issue is that I have the Basemap over my data. Do you have an idea on how to inverse this? (an other way than adding an alpha in my Basemap) – Thib Ödlaniug Mar 27 '18 at 08:49
  • Sorry, I'm getting a bit tired of having to guess all the time. If you have a problem with code, make it reproducible for others as in [mcve]. – ImportanceOfBeingErnest Mar 27 '18 at 10:51
  • I totally understand your point but I can't provide the data corresponding to the issue. Thank you a lot for your help – Thib Ödlaniug Mar 27 '18 at 15:35
  • It will sure be possible for you to create some dataset programmatically to make this reproducible, right? – ImportanceOfBeingErnest Mar 27 '18 at 15:53