# Where is my mistake in this MATLAB code?

1. Oct 25, 2012

### EngWiPy

Hi,

I am trying to plot the AN product sec II.C in the attached paper as:

Code (Text):
clear all
clc

s=0.5;
w=0;
d=5;

f=0:20;

NdB1=(17-30.*log10(f));
N1=10.^(NdB1./10);
NdB2=40+20*(s-0.5)+26.*log10(f)-60.*log10(f+0.03);
N2=10.^(NdB2./10);
NdB3=50+7.5*w^0.5+20.*log10(f)-40.*log10(f+0.4);
N3=10.^(NdB3./10);
NdB4=-15+20.*log10(f);
N4=10.^(NdB4./10);

N=N1+N2+N3+N4;

A=(d^(1.5)).*a.^d;

A_N=1./(A.*N);

A_NdB=10*log10(A_N);

plot(f,A_NdB)
The details of these equations are highlighted in the attached file.

I appreciate if someone points to me my mistake, and the figure should be like Fig 3

Thanks

#### Attached Files:

• ###### MITSG_07-13J.pdf
File size:
265.6 KB
Views:
93
2. Oct 28, 2012

### NemoReally

What problem are you seeing? I'm not that familiar with Matlab, but one potential problem is trying to take log(0) ... how does Matlab handle this (eg, error, NaN)?

I've implemented your equations in Mathcad and the shape of the resultant plot seems to agree with the document's plot, seemingly differing by an offset and/or scaling factor - the peaks seem to agree and the overall size looks OK, so I suspect an offset of about 45 dB.

#### Attached Files:

• ###### phys - 12 10 27 acoustic underwater comms 01.jpg
File size:
74.1 KB
Views:
226
3. Oct 28, 2012

### EngWiPy

4. Oct 28, 2012

### NemoReally

I'm (hopefully) plotting the same range as you, on the assumption that your code was specialized for kHz. Here's an image of the (attached) Mathcad code.

#### Attached Files:

File size:
16.1 KB
Views:
218
• ###### phys - 12 10 27 acoustic underwater comms 02.mcd
File size:
14.8 KB
Views:
62
5. Oct 28, 2012

### EngWiPy

Are you using log10 or natural logarithm in your code?

It looks like very much the same as my code!!

6. Oct 29, 2012

### NemoReally

Yes, it does look pretty much like your code, so that's why I was wondering what problem you're seeing.

Mathcad's "log" is log10 by default; Mathcad uses "ln" for loge. Whereas, it seems, Matlab's "log" is loge. Using loge results in plots that look markedly different from the article's plots - far more peaky.

7. Oct 30, 2012

### uart

There's no real problem with your code, it's just the normalizing constant $A_0$ causing the dB offset. That and the fact that you are plotting at too few values of "f", especially at low frequencies, to get a very accurate picture of what's going on.

Just replace the line

"f = 0:20"

with

"f=[0.1 : 0.1 : 20]" (or even better with f = logspace(-1,log10(20),100)

and you'll get exactly the same graph as figure 3 in the paper (k=1.5, w=0, s = 0.5 and d=5km), except for about a 45 dB offset, presumably due to A0.

BTW. You could have made it a lot easier for people to see what your problem was if you had of actually described the errors or anomalies you were concerned about (or even better posted an image of the plot you were obtaining).

Last edited: Oct 30, 2012
8. Oct 30, 2012

### EngWiPy

I am sorry if I were not clear, and thanks for your reply.

For d=5 Km the curves are very similar except for a certain scale, but for d=100 Km there is a huge difference!! Is that because of A0, too?

Thanks

9. Oct 30, 2012

### uart

Ok I don't know, I only checked the code you posted and it gave the correct results (except for the offset) for d=5. The paper refers to a normalizing constant A0, which is not included in your equations, so certainly that could account for the dB offset. The other problem was that the result is -inf dB at f=0, so that result wasn't being shown on the graph, and plotting at just the integer frequencies 1:20 was not enough to get a good picture of the graph.

10. Oct 30, 2012

### uart

Ok I just checked the result at d=100 and it also agrees exactly with the graph (fig3) in the paper (except for the offset due to whatever normalizing constant they're using).

For d=100 I used the freq range: f = logspace(-2, log10(7), 100) to reproduce the figure in the paper (same shape exactly, but again with an offset of about 45 dB).

11. Oct 30, 2012

### EngWiPy

Thank you so much, it is working now. Just the scaling factor is different which is of minor importance.