How to plot a line on top of the continuous wavelet transform (CWT) output?

Hello. According to the attached figure, I want to add a horizontal green line (at a scale of 60) on top of the continuous wavelet transform (CWT) output. However, the line is positioned below the CWT output. Your assistance with this request would be greatly appreciated.

Answer 1

[dpb]

data = xlsread('signal.xlsx','1');
% Extract the time and signal columns
t = data(:, 1);
signal = data(:, 2);
dt = t(3) - t(2); % sampling time
fs = 1/dt; % freq
wname = 'bior6.8';
N = length(signal);
tv = (0:N-1)/fs;
% Compute the CWT coefficients
scales = 1:200;
C = cwt(signal,scales,wname);
% Create a meshgrid fot the time and scales
[T,S] = meshgrid(tv,scales);
% Plot the CWT coefficients with time and scales
figure;
hS=surf(T,S,abs(C),'EdgeColor','none');
% leave in 3D view for now for demo purposes
%view(0,90); % Set the view to 2D
xlabel('Time (sec)');
ylabel('Scale');
xlim([0, T(end)]);
grid off
colormap (pink(240));
set(gca,'fontname','Times New Roman','FontSize',10)
hold on
SCA = 60;
% left following in temporarily for demo purposes...
plot(xlim,[1 1]*(SCA),'-r','LineWidth',2)

Illustrates the problem/issue; a 2D line parallel to the x axis with its z height equal to 0 is obscured by the surface and the color map chosen is dark enough to prevent the line from showing through. Making the surface less opaque lets one see where the line is located, but it took going to such a low value for .FaceAlpha that the surface itself is then difficult to see for the lighter shades.

So, how to fix? The problem is one needs to draw the line at or above the surface in the 3D axes space even though the desired presentation is from the Z axis direction directly towards the x,y plane. Hence, yline or its workalikes aren't the correct tool for the job, we need plot3 instead.

Z=max(abs(C),[],'all');         % max Z so can be at or above any point in surface
hL3=plot3(xlim,[SCA SCA],[Z Z],'g-','linewidth',2);
pause(10)                       % for demonstration -- remove when going to final use
view(0,90);                     % Set the view to 2D after seeing original for a while

I left in the initial line on the Z==0 plane and commented on it and later...to illustrate, the above code as is will show the surface and then the red line in 3D view and add the green line at Z==max(C), then pause for 10 seeconds so can see the result before then switch the viewpoint to the 2D plane view at which time the green line will appear on top and obscure the red line. Remove those extra pieces to use "in anger", obviously.

Sorry it took me so long to think of the cause of the problem after which the solution is obvious...

[Navid]

Dear @dpb,

Thank you for your time and support. When I ran the code, I encountered the following error:

"Error using max

Invalid option. Option must be 'omitnan' or 'includenan'."

It seems that the 'all' option is not supported in my version. Instead, I tried the following:

[~, idx] = max(abs(C));
Z=max(abs(C(:)));

This change finally worked. Once again, I appreciate your valuable guidance. Thank you.

[dpb]

Glad to help...as noted, I seemed to have had a short-circuit in not realizing the issue clearly from the beginning...as this illustrated, posting sufficient information/data to let folks here reproduce an issue is key to getting accurate and timely answers.

I forgot you were on an earlier release; don't recall and didn't look to see when the 'all' option was introduced.

The max(X(:)) idiom is the best solution failing having it; you'll see a lot of code written with something like

Z=max(max(X));

which works for arrays up to 2D, but breaks after that.

[~, idx] = max(abs(C));

isn't needed here; as the above note refers to, idx will be a row vector of the maximum positions for each column in C as max operates by column if the optional dim argument isn't given (which is the null [] placeholder in the original code).

[Navid]

I want to express my gratitude once again. Thank you sincerely for your valuable time.

[dpb]

You're welcome...glad to help.

Answer 2

[dpb]

From one of the examples it appears that other than (rudely) executing a clf that clears any current figure when called, it looks like cwt simply creates an ordinary axis that can be accessed in the normal HG2 ways...

load kobe
cwt(kobe,1)
hAx=gca;
hold on
yline(60,'g','linewidth',2)

Attach a minimum working example that produces your example plot above if want/need something more specific,

[Navid]

Thank you for your response. The "yline" syntax is not defined in my version (R2016b).

[dpb]

Use line or plot then...

[dpb]

A "poor man's" yline substitute for 2D axes could be as simple as

function hL=yLine(yval,w,c,s)
% rudimentary yline substitute for 2D axes, not vectorized
% USAGE
%  hL = yLine(yValue,Color,LineStyle,LineWidth)

  % set default values if not provided
  if nargin<4, s='-'; end       % solid line
  if nargin<3, c='k'; end      % black color
  if nargin<2, w=1;   end      % unity width
      
  hAx=gca;
  hold on
  hAx.YAxis.LimitMode='manual'; % keep from autoranging on us
  h=line(xlim,[y y],'color',c,'linestyle',s,'linewidth',w)
  if nargout>0, hL=h; end
end

[Walter Roberson]

A few improvements...

function hL=yLine(yval,w,c,s)
% rudimentary yline substitute for 2D axes, not vectorized
% USAGE
%  hL = yLine(yValue,Color,LineStyle,LineWidth)

  % set default values if not provided
  if nargin<4, s='-';       % solid line
  if nargin<3, c='k';       % black color
  if nargin<2, w=1;         % unity width
      
  hAx=gca;
  CurrentHoldState = get(hAx, 'NextPlot');
  hold(hAx, 'on')
  hAx.YAxis.LimitMode='manual'; % keep from autoranging on us
  h = line(hAx, xlim, [y y], 'color', c, 'linestyle', s, 'linewidth', w);
  set(hAx, 'NextPlot', CurrentHoldState);
  if nargout>0, hL=h; end
end

[dpb]

A bug fix I see I forgot to include, Walter...

function hL=yLine(yval,w,c,s)
% rudimentary yline substitute for 2D axes, not vectorized
% USAGE
%  hL = yLine(yValue,Color,LineStyle,LineWidth)

  % set default values if not provided
  if nargin<4, s='-'; end      % solid line
  if nargin<3, c='k'; end      % black color
  if nargin<2, w=1;   end      % unity width
      
  hAx=gca;
  CurrentHoldState = get(hAx, 'NextPlot');
  hold(hAx, 'on')
  hAx.YAxis.LimitMode='manual'; % keep from autoranging on us
  h = line(hAx, xlim, [y y], 'color', c, 'linestyle', s, 'linewidth', w);
  set(hAx, 'NextPlot', CurrentHoldState);
  if nargout>0, hL=h; end
end

The "end"s were missing on the default if clauses...

On retrospect I think probably the autorange is not a good idea; yline will if use a value outside the current axes ylim values.

[Navid]

Dear Walter Roberson,

As always, thank you for your detailed response. Please let me restate my question. I am following up on my previous inquiry about adding a horizontal line to the output of the wavelet transform.

Wavelet transform output has been obtained from:

C = cwt (signal, scales, wname);
[T,S] = meshgrid (tv, scales);
% where tv = (0:N-1)/fs, N = signal length & fs = sampling rate
figure;
surf (T, S, abs(C), 'EdgeColor', 'none')
view (0, 90);

Unfortunately, after adding the desired horizontal line, it appears below the wavelet transform figure (see attached image). I want to bring this line forward to make it more visible. Your attention to this matter is greatly appreciated.

Best regards,

Navid

Attached image:

[dpb]

It was my code @Walter Roberson added an enhancement to, just to be picky, @Navid.

C = cwt (signal, scales, wname);

Unrecognized function or variable 'signal'.

[T,S] = meshgrid (tv, scales);
% where tv = (0:N-1)/fs, N = signal length & fs = sampling rate
figure;
surf (T, S, abs(C), 'EdgeColor', 'none')
view (0, 90);

Again, we can't duplicate your issue without the data with which to generate the figure and ALL the code needed to create a minimum working example.

Post everything needed to make the shown figure such that somebody here can duplicate the issue.

Your posting doesn't even show the horizontal line code you used at all...

ADDENDUM

I suppose perhaps could get by if you would save the output figure to a .fig file and attach it...then could at least load the figure to mung on...

[dpb]

Lacking real data, we'll use one of the surf examples to 'spearmint with...

[X,Y,Z] = peaks(25);
CO(:,:,1) = zeros(25); % red
CO(:,:,2) = ones(25).*linspace(0.5,0.6,25); % green
CO(:,:,3) = ones(25).*linspace(0,1,25); % blue
hS=surf(X,Y,Z,CO,'edgecolor','none');
view(0,90)
xlabel('X'),ylabel('Y')
yline(0,'r')

I thought maybe the peaks would hide the line but they didn't seem to...we simply must have the actual code.

[Navid]

Dear @dpb,

Thank you for your guidance and time. I am sharing my code for further review. Upon running the code, I noticed that the horizontal line appears behind the figure.

I sincerely appreciate your help and look forward to any input or suggestions you may have. Thank you for your time and support.

Best regards,

Navid

clc
clear
close all

data = xlsread('signal.xlsx','1');
% Extract the time and signal columns
t = data(:, 1);
signal = data(:, 2);
dt = t(3) - t(2); % sampling time
fs = 1/dt; % freq
wname = 'bior6.8';
N = length(signal);
tv = (0:N-1)/fs;
% Compute the CWT coefficients
scales = 1:200;
C = cwt(signal,scales,wname);
% Create a meshgrid fot the time and scales
[T,S] = meshgrid(tv,scales);
% Plot the CWT coefficients with time and scales
figure;
surf(T,S,abs(C),'EdgeColor','none')
view(0,90); % Set the view to 2D
xlabel('Time (sec)');
ylabel('Scale');
xlim([0, T(end)]);
grid off
colormap (pink(240));
set(gca,'fontname','Times New Roman','FontSize',10)
hold on
SCA = 60;
plot(xlim,[1 1]*(SCA),'-r','LineWidth',2)
hold off