Re: onde em raddrizzate. Esistono ?

From: JTS <pireddag_at_hotmail.com>
Date: Mon, 11 Dec 2017 21:42:48 +0100

Am 01.12.2017 um 22:56 schrieb JTS:

>
> Dovresti quindi nuovamente ricevere un segnale a media nulla, e di forma
> diversa da quello che avevi inviato. Ho provato un piccolo calcolo
> numerico e mi pare che dia risultati ragionevoli. Lo metto in ordine e
> poi lo posto.


Ecco il codice per il calcolo. Lo ho scritto per Matlab, dovrebbe
funzionare anche con Octave (che e' gratis). Se non dovesse funzionare
su Octave sono disponibile a metterlo a posto (se interessa a qualcuno).

Lo copio qui sotto, lo ho scritto rispettando il limite delle 72
colonne, perche' cosi' vengono mantenute intiere le linee

%% Script test_dc_filtering_electric_dipole_radiation

% This script takes an oscillatory signal, rectifies it, considers the
% rectified signal as the time-dependent amplitude of an electric dipole
% and calculates the corresponding radiated field. One can observe that
% the time dependence of the radiated field is different from the time
% dependence of the dipole, as propagation suppresses low frequencies
% and enchances high frequencies.


%% Parameters

% It is possible to consider an oscillating signal limited in time (by
% setting gaussian_evenlope to true) or a purely cosinusoidal
% oscillation (gaussian_envelope = false)
gaussian_envelope = true;
% gaussian_envelope = false;

w0 = 0.625; % [Hz] carrier
w1 = 0.125; % [Hz] signal

A = 10;
m = 1;
% case of signal limited in time: width of the signal in the time domain
% (width of gaussian)
delta_t = 20;

% A parameter for smoothing the rectification
rectification_constant = 1/6;

%%

t = linspace(-96, 96, 4096*4);

% Set up an axis for the FFT (see Numerical Recipes and Matlab help)
omega= 2*pi*linspace ...
     (-1/(2*(t(2)-t(1))), ...
     1/(2*(t(2)-t(1))) - 1/(length(t)*(t(2)-t(1))),...
     length (t));
% Note: not tested enough and I fear that there may be a mistake, so
% that the zero frequency element of the fft does not correspond to the
% zero in the omega array; but the error should be small

%% Calculations

% A single-sideband signal and its fft; I consider only the upper band
% and I do not include the carrier.

% As far as I can see, multiplying the original signa times a gaussian
% in the time domain reduces the high frequencies of the rectified
% signal, but I do not understand why. A way to look for the reason is
% to compare the following two procedures: first multiple times a
% gaussian, then rectify, or first rectify, then multiply times a
% gaussian. We know that in the frequency domain they should give rise
% to different spectra, but it is interesting to look at the differences
% in the time domain, as they could give a hint to the reason for the
% differences.

if gaussian_envelope
     signal_ub = (m*A/2)*cos((w0+w1)*t).*exp(-t.^2/delta_t^2);
else
     signal_ub = (m*A/2)*cos((w0+w1)*t);
end
signal_ub_ft = fft_w_shift(signal_ub);

% Rectify the signal, eliminating the negative parts.

% I want smooth transitions in order to limit the maximum frequency of
% the rectified signal. To achieve this I keep a small portion of the
% negative part of the signal.

% I proceed in this way. Call f the absolute value of the signal. I
% multiply the negative part of the signal times 1/(1+(1/r)*f/max(f)).
% In this way the negative part of the function is limited smoothly to
% max(f)/r

% signal_ub_rectified = signal_ub.*(signal_ub > 0);
signal_ub_rectified_positive = signal_ub.*(signal_ub > 0);
signal_ub_rectified_negative = signal_ub.*(signal_ub < 0).* ...
     (1./(1+abs(signal_ub)/(max(signal_ub)*rectification_constant)));
signal_ub_rectified = ...
     signal_ub_rectified_positive + signal_ub_rectified_negative;

% Propagate, using the omega.^2 propagator in the frequency domain (see
% the formula for the radiation of a dipole).

% Use rhe Matlab fftshift and ifftshift functions to bring the zro
% frequencies to the center of the arrays
signal_ub_rectified_ft = ...
     fftshift(fft(ifftshift(signal_ub_rectified)));

signal_ub_propagated_ft = signal_ub_rectified_ft .* omega.^2;

signal_ub_propagated = ...
     ifftshift(ifft(fftshift(signal_ub_propagated_ft)));

%% Graphics

figure(101)
clf
set(gcf, 'Name', 'Signal in frequency domain and propagator')

[AX,H1,H2] = plotyy(omega,abs(signal_ub_rectified_ft), omega, omega.^2);
title('Signal in frequency domain and propagator')
xlabel('frequency (Hz)')

set(get(AX(1),'Ylabel'),'String','signal amplitude (arb. units)')
set(get(AX(2),'Ylabel'),'String','propagator (arb. units)')

xlim(AX(1),[-15 15])
xlim(AX(2),[-15 15])
ylim(AX(1),[-1000 5000])
ylim(AX(2),[-40 200])
set(AX(2),'ytick', [0 50 100 150])
  % Turn off box of axis 1, which removes its right-hand ticks
set(AX(1),'Box','off')
% Turn off box of axis 2, which removes its left-hand ticks
set(AX(2),'Box','off')
% Set x-axis on top for ax(2) and remove the ticks: this completes the
% box
set(AX(2),'XAxisLocation','top','XTickLabel',[]);

figure(102)
clf
set(gcf, 'Name', 'Signal in frequency domain')

subplot(3,1,1)
plot(omega, abs(signal_ub_ft))
title('Signal in frequency domain')
set(gca,'yTick',[0 4000 8000])
xlim([-15 15])
ylim([-2000 8000])

subplot(3,1,2)
% Zoom on y-axis to see high frequency components resulting from
% rectification
plot(omega,abs(signal_ub_rectified_ft))
title('Rectified signal')
ylabel('signal amplitude (arb. units)')
set(gca,'yTick',[0 4000 8000])
xlim([-15 15])
ylim([-2000 8000])

subplot(3,1,3)
% High frequency components resulting from rectification are enhanced by
% propagation relative to low-frequency components. DC is not propagated
% at all (the propagator is zero)
plot(omega, abs(signal_ub_propagated_ft))
title('Propagated signal')
xlabel('frequency (Hz)')
set(gca,'yTick',[0 4000 8000])
xlim([-15 15])
ylim([-2000 8000])

figure(103)
clf
set(gcf, 'Name', 'Signal in time domain')

subplot(3,1,1)
plot(t, signal_ub)
title('Signal in time domain')
subplot(3,1,2)
plot(t, signal_ub_rectified)
title('Rectified signal')
ylabel('signal amplitude (arb. units)')
subplot(3,1,3)
plot(t, real(signal_ub_propagated))
% Here I take the real part of the propagated signal; the very small
% imaginary part comes from numerical errors
title('Propagated signal (note different scale of y-axis)')
xlabel('time (s)')
Received on Mon Dec 11 2017 - 21:42:48 CET