Color plot Y2 as a function of X,Y1

Hi guys,

I have salinity measurements (y2) in different depths (y1) from different distances (x) from the shore. I would like to colormap the salinity values and get something like this:

What I've done so far is that I combined the depth values to get a common y-axis and filled up the missing values with NaN. Then I combined it with the distance values (X) to get a grid. At this point I tried several functions but I couldnt solve the issue. Any help is greatly appreciated.

Answer 1

[Chad Greene]

I'd probably concatenate the x position and depth of every salinity measurement, then grid them up with gridfit. Here's an example that you can run for four CTD casts.

% x locations along a transect:
x = [0 10 12 18]; 

% 30 or so random depths from 0 to 600ish dbar:
d1 = sort(500*rand(30,1)); 
d2 = sort(540*rand(25,1)); 
d3 = sort(585*rand(35,1)); 
d4 = sort(510*rand(28,1)); 

% Salinity values at each measurement: 
s1 = 33 + (.002*d1).^3; 
s2 = 33.35 + (.002*d2).^3; 
s3 = 33.3 + (.002*d3).^3; 
s4 = 33.2 + (.002*d4).^3; 

% Plot the salinity measurements: 
scatter(x(1)*ones(size(s1)),d1,40,s1,'filled'); 
hold on
scatter(x(2)*ones(size(s2)),d2,40,s2,'filled'); 
scatter(x(3)*ones(size(s3)),d3,40,s3,'filled'); 
scatter(x(4)*ones(size(s4)),d4,40,s4,'filled'); 

% Make a colorbar:
cb = colorbar; 

% Set colormap: 
cmocean haline

% Label axes: 
xlabel('distance along transect (km)') 
ylabel('pressure (dbar)') 
ylabel(cb,'salinity') 

% reverse order of y axes: 
axis ij tight

Note, instead of the hideous default ODV rainbow color scheme I'm using cmocean.

Now to fill in the space between each measurement the way ODV does we need to interpolate to get a value at each pixel on the screen. There are different ways of interpolating and each method can have drawbacks. With interp1 you can interpolate each profile to common depth values and then use interp1 again to interpolate those interpolated values to distances along the transect, but that will leave some big gaps in the data (and such a method require more manual steps on your part than necessary).

Alternatively, you can use the built in functions griddata or scatteredInterpolant, but both tend to produce nonphysical artifacts like striping or voronoi type patterns. Thus, I prefer the simple, elegant method of using John D'Errico's gridfit like this:

% x values corresponding to each measurement: 
xi = [repmat(x(1),size(d1));repmat(x(2),size(d2));repmat(x(3),size(d3));repmat(x(4),size(d4))]; 

% All depths: 
di = [d1;d2;d3;d4]; 

% All salinity measurments: 
si = [s1;s2;s3;s4]; 

% Define some common distances along transect: 
X = 0:0.1:20; 

% Define some common depths: 
D = 0:5:600; 
S = gridfit(xi,di,si,X,D); 

p = pcolor(X,D,S) ; 
shading interp
uistack(p,'bottom') 

% Clarify the location of each measurement: 
plot(xi,di,'.','color',0.5*[1 1 1])

% Plot bathymetry: 
p = patch([x 20 0],[max(d1) max(d2) max(d3) max(d4) 600 600],'k');
set(p,'facecolor',0.5*[1 1 1])