Sinusoidal signal with variable instantaneous frequency

Below is a simple matlab code that generates a sinusoidal signal with user specified no of cycles. How do i change this such that each cycle has different instantaneous frequency? For example, if i have 10 cycles, then there should be 10 random frequencies in total where each cycle has different frequency.

clc
clear

f1=120;
fs=f1/60; %sampling frequency depends on f1
no_of_cycles=10;
f=no_of_cycles/60; %Frequency of sinusoid
duration=1;%%duration of the sinusoidal signal in minutes
A=1;
t=1/fs:1/fs:duration*no_of_cycles*1/f; %time index
y1=A*sin(2*pi*f*t);%%simulated sinusoidal signal

Answer 1

[John BG]

Hi Sowmia MR

This is John BG ( <mailto:[e-mail] [e-mail]> )

1.

One way to answer your question is, without losing constant amplitude, which is a good quality of FM:

f1_range=120*[1:1/3:4]
f1=f1_range(randi([1 9],1,9))
A=1;
y1=0

for k=1:1:numel(f1)
    y1c=A*sin(2*pi/f1(k)*[1:1:f1(k)]);
    y1=[y1 y1c];
end
plot(y1)

.

2.

Another way would be constraining the amount of samples for each frequency to a constant value, that would be of 1 cycle of f0=120Hz

    f0=120
    f1_range=f0*[.7:.1:1.3]
    f1=f1_range(randi([1 numel(f1_range)],1,numel(f1_range)))
    A=1;
    y1=0

    for k=1:1:numel(f1)
        y1c=A*sin(2*pi*f1(k)/f0*[0:1/f0:1]);
        y1=[y1 y1c];
    end
    plot(y1)

.

3.

Note that now it's not 10 cycles but 7, 7*120 are the 840 samples shown above 2nd graph.

Recommended reading 'Telecommunication Systems' by Bruce Carlson

<http://een.iust.ac.ir/profs/Beheshti/communication%20system/Carlson%20-%20Communication%20Systems%20-%204ed%20.pdf>

4.

Please not that f1_range is on purpose chosen 'not too wide'.

If so alias occurs, data loss, the output signal and plot may be meaningless.

5.

If you want to learn more about FM basics with MATLAB, please check my accepted answer to question in this link:

https://uk.mathworks.com/matlabcentral/answers/318969-fm-modulation-without-using-fmmod

there's an explanation how to FM modulate a signal, with MATLAB and with SIMULINK

So, Sowmia MR

if you find this answer useful would you please be so kind to consider marking my answer as Accepted Answer?

To any other reader, if you find this answer useful please consider clicking on the thumbs-up vote link

thanks in advance

John BG

<mailto:[e-mail] [e-mail]>

[Sowmya MR]

Thank you. This is exactly what i required.

[Sowmya MR]

Hi John, I was playing with the code and realized that each cycle in example 1 has varying number of samples. Is there anyway we can have same number of samples (say 32 samples) for every cycle and still get a clean signal?

[John BG]

Yes, it can be done, you want 120 samples in each cycle, AND amplitude continuity.

One can achieve such amplitude continuity, with a phase carry:

clear all;clc;clf
f0=120
f1_range=120*[1:1/3:4]
f1=f1_range(randi([1 9],1,9))
A=1;
y1=0
dph=0
for k=1:1:numel(f1)
    y1c=A*sin(2*pi*f1(k)/f0*[0:1/f0:f1(k)/f0-1/f0]+dph);
    y1=[y1 y1c([1:1:f0])];
    dph=floor(f0*rem(f1(k),f0)/f0)+dph

end

plot(y1)

the random sequence of frequencies of this example is

f1 =

240.0000  320.0000  160.0000  360.0000  120.0000  200.0000  120.0000  120.0000  400.0000

.

Repeat to test different f1 patterns.

Note that sometimes there is a small gl

itch

It's caused by the rounding of the command floor and the simplified way the carry is implemented.

You may find also of interest the Hull University introduction to GMSK, the modulation used in GSM, document available here:

http://www.emc.york.ac.uk/reports/linkpcp/appD.pdf

GSM uses the modulation GMSK that departs from this basic FM (or FSK) precisely to avoid constellation zero crossings (not to be confused with base band signal zero crossings or modulated signal zero crossings) avoiding sharp transitions of the VCO input signal using a gaussian filter.

GMSK also uses In phase and Quadrature, that is more robust than just using one real carrier.

Regards

John BG

<mailto:[e-mail] [e-mail]>

[John BG]

removing command floor

clear all;clc;clf
f0=120
f1_range=120*[1:1/3:4]
f1=f1_range(randi([1 9],1,9))
A=1;
y1=0
dph=0
for k=1:1:numel(f1)
    y1c=A*sin(2*pi*f1(k)/f0*[0:1/f0:f1(k)/f0-1/f0]+dph);
    y1=[y1 y1c([1:1:f0])];
    dph=rem(f1(k),f0)+dph
end
plot(y1)

fixes the phase glitches shown in the last graph of this answer: no floor better.

John BG

<mailto:[e-mail] [e-mail]>

Answer 2

[Walter Roberson]

Figure out how many samples there are per "cycle" for the purpose of your "10 cycles". The last time for each cycle should end 1 sample before the next full cycle boundary. So for example 10 samples per "cycle" should correspond to relative times 0, 1/10, 2/10, 3/10, 4/10, 5/10, 6/10, 7/10, 8/10, 9/10, and should not go all the way to 1 (which would be 11 samples). If you end up using linspace along the way to generate N points, then use linspace(first_time, last_time, N+1 ) and then throw away the last sample.

Now, generate the list of instantaneous frequencies for each cycle. You should probably limit those to being no more than 1/2 the sampling frequency, to avoid aliasing.

instant_frequency_for_cycle = fs/2 * rand(1, no_of_cycles)

Now replicate each of those by the number of samples per cycle:

f = reshape( repmat( instant_frequency_for_cycle(:), 1, samples_per_cycle ), 1, [])

or if you have a sufficiently new MATLAB,

f = repelem( instant_frequency_for_cycle, samples_per_cycle )

or in any version you can use

f = kron(instant_frequency_for_cycle, ones(1, samples_per_cycle) )

After that it is a simple modification to what you had:

y1 = A * sin(2*pi*f.*t);   %%simulated sinusoidal signal

the f*t got changed to f.*t . With f and t being the same length, this gives an instant frequency to every sample.

[John BG]

Hi this is John BG ( <mailto:[e-mail] [e-mail]> )

1.

Mr Roberson answer only works when

samples_per_cycle=12

for any other value of samples_per_cycle different than 12 the last line returns mismatch error.

2.

Also, not constraining the frequency deviation to a small value compared to the carrier frequency loses some of the good properties of FM, like superior SNR compared to AM, following y1 with Walter's lines

.

.

because with Walter's lines the amplitude is no longer constant, it may leak FM into AM, which may not be desirable if you are trying to conceal the information in the carrier phase only.

I have abridged some of the FM key components explained in

https://en.wikipedia.org/wiki/Frequency_modulation

as follow:

[Sowmya MR]

Hi John, Thanks for the great insight. Could you please help me with changing my code to work for any cycles?

[Walter Roberson]

John, you missed my "With f and t being the same length". My outline does not care about the 120 of the original code. I did not assume that "cycle" for the purpose of "10 cycles" corresponds to "hertz".