[Top][All Lists]

[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

Overlaying map in octave

From: Zilore Mumba
Subject: Overlaying map in octave
Date: Thu, 9 Oct 2014 12:31:42 +0200

I noticed that octave has a possibility of using the m_map utility of matlab to overaly a basemap over contours.
I have not installed m_map but when I do "apropos m_map" on my ubuntu terminal I get "spatialite_osm_map (1) - load OSM-XML maps into a SpatiaLite DB" which I understood to mean that I have m_map.
I then tried example "

5. Oblique Mercator Projection with quiver and contour data" on
with my modifications, of course.
But then I get the error
"error: 'm_proj' undefined near line 67 column 1
error: called from:
error:   /home/zmumba/Optiminterp/example_optiminterp_data.m at line 67, column 1"

Any help or advice will be appreciated.

My full code is below

% load the data file
% the 1st and 2nd column represent the location
x = data(:,1);
y = data(:,2);
% the 3rd column represents the observations to interpolate
obs = data(:,3);

% create a regular mesh from -2 to 2 in x and y direction
% with a grid spacing of 0.1

%[xi,yi] = ndgrid(-2:0.1:2,-2:0.1:2);
[xi,yi] = ndgrid(22:0.1:34,-18:0.1:-8);

% number of influential points

m = 29;

% compute the first guess (or background), here we take the
% mean of all observations
% mean of all observations
mean_obs = mean(obs);

% remove the mean
f = obs - mean_obs;

% var is the error variance ratio between of the observations and background.
% observations are assumed to be 100 times more accurate (in terms of error variance) than the
% background estimation (here the mean)

var = 0.01 * ones(size(f));

% call optiminterp for 2d case with a correlation length of 1.5 in x and y direction

[fi,vari] = optiminterp2(x,y,f,var,1.5,1.5,m,xi,yi);

% put the mean back
obsi = fi + mean_obs;

fid = fopen('mydata.txt','w');
%fprintf(fid, '%5.2f %5.1f %5.1f\n', xi,yi,obsi);

for i=1:numel(xi),
   fprintf(fid, '%5.1f %5.1f %5.1f\n', xi(i), yi(i), obsi(i));

% Define my map grid

% define map projection
m_proj('mercator','lat',[-8 -18],'lon',[22 34],'aspect',.8);

m_coast('patch',[.9 .9 .9],'edgecolor','none');

hold on
xlabel('Simulated something else');


reply via email to

[Prev in Thread] Current Thread [Next in Thread]