% Symbol rate processing for rx side (demodulator) -------------------------------------------------------
-function [rx_symb rx_bits rx_symb_linear amp_ phi_ EsNo_ cohpsk] = qpsk_symbols_to_bits(cohpsk, ct_symb_buf)
+function [rx_symb rx_bits rx_symb_linear amp_ phi_ sig_rms noise_rms cohpsk] = qpsk_symbols_to_bits(cohpsk, ct_symb_buf)
framesize = cohpsk.framesize;
Nsymb = cohpsk.Nsymb;
Nsymbrow = cohpsk.Nsymbrow;
% now correct phase of data symbols
rx_symb = zeros(Nsymbrow, Nc);
- rx_symb_linear = zeros(1, Nsymbrow*Nc);
+ rx_symb_linear = zeros(1, Nsymbrow*Nc*Nd);
rx_bits = zeros(1, framesize);
for c=1:Nc*Nd
for r=1:Nsymbrow
rx_symb(r,c) = ct_symb_buf(2+r,c);
end
i = (c-1)*Nsymbrow + r;
+ rx_symb_linear(i) = rx_symb(r,c);
%printf("phi_ %d %d %f %f\n", r,c,real(exp(-j*phi_(r,c))), imag(exp(-j*phi_(r,c))));
end
end
for d=1:Nd-1
div_symb += rx_symb(r,c + Nc*d);
end
- rx_symb_linear(i) = div_symb;
rx_bits((2*(i-1)+1):(2*i)) = qpsk_demod(div_symb);
end
end
% position. However this is complicated by fading, which means the
% amplitude of the symbols is constantly changing.
- % Now the scatter diagram in a fading channel is a X shape. The
- % noise can be resolved into two components at right angles to
- % each other. The component along the the "thickness" of the arms
- % is proportional to the noise power and not affected by fading.
+ % Now the scatter diagram in a fading channel is a X or cross
+ % shape. The noise can be resolved into two components at right
+ % angles to each other. The component along the the "thickness"
+ % of the arms is proportional to the noise power and not affected
+ % by fading. We only use points further along the real axis than
+ % the mean amplitude so we keep out of the central nosiey blob
- v = zeros(1,Nsymb);
- for i=1:Nsymb
+ sig_rms = mean(abs(rx_symb_linear));
+
+ sum_x = 0;
+ sum_xx = 0;
+ n = 0;
+ for i=1:Nsymb*Nd
s = rx_symb_linear(i);
- if abs(real(s)) > abs(imag(s))
- v(i) = imag(s);
- else
- v(i) = real(s);
+ if abs(real(s)) > sig_rms
+ sum_x += imag(s);
+ sum_xx += imag(s)*imag(s);
+ n++;
end
- %printf("s: %f %f v: %f\n", real(s), imag(s), v(i));
end
-
- % Note we are only measuring variance in one axis, as other axis is obscured by fading. We assume
- % that total noise power is shared between both axis so multiply by sqrt(2) to get an estimate of
- % total noise pwr. Small constant prevents divide by zero errors on start up.
-
- No_ = var(v)*sqrt(2) + 1E-6;
-
- % Estimate signal power
-
- Es_ = mean(amp_ .^ 2);
-
- EsNo_ = Es_/No_;
- %printf("Es_: %f No_: %f Es/No: %f Es/No dB: %f\n", Es_, No_, Es_/No_, 10*log10(EsNo_));
+
+ noise_var = 0;
+ if n > 1
+ noise_var = (n*sum_xx - sum_x*sum_x)/(n*(n-1));
+ end
+ noise_rms = sqrt(noise_var);
endfunction
% Plot Octave outputs from cohpsk_demod, c2dec, to visualise whats going on
% when errors hit the system
-% $ ./c2enc 700 ../../raw/ve9qrp_10s.raw - | ./cohpsk_mod - - | ./cohpsk_ch - - -60 50 1 1 | ./cohpsk_demod - - | ./c2dec 700 - - --dump ve9qrp | play -t raw -r 8000 -s -2 - -q
+% $ ./c2enc 700 ../../raw/ve9qrp_10s.raw - | ./cohpsk_mod - - | ./cohpsk_ch - - -60 50 1 1 | ./cohpsk_demod - - cohpsk_demod.txt | ./c2dec 700 - - --dump ve9qrp | play -t raw -r 8000 -s -2 - -q
-% ./c2enc 700 ../../raw/ve9qrp_10s.raw - | ./cohpsk_mod - - | ./cohpsk_ch - - -30 50 1 1 | ./cohpsk_demod - - | ./c2dec 700 - - --dump ve9qrp_snr3 | play -t raw -r 8000 -s -2 - -q
+% ./c2enc 700 ../../raw/ve9qrp_10s.raw - | ./cohpsk_mod - - | ./cohpsk_ch - - -30 50 1 1 | ./cohpsk_demod - - cohpsk_demod.txt | ./c2dec 700 - - --dump ve9qrp_snr3 | play -t raw -r 8000 -s -2 - -q
graphics_toolkit ("gnuplot");
figure(8)
clf;
-f1=fft(ve9qrp_ak_(st:en,:)',128);
-f2=fft(ve9qrp_snr3_ak_(1:en-st+1,:)',128);
-d = (20*log10(abs(f1)) - 20*log10(abs(f2)));
+f1=1./fft(ve9qrp_ak_(st:en,:)',128);
+f2=1./fft(ve9qrp_snr3_ak_(1:en-st+1,:)',128);
+%d = (20*log10(abs(f1)) - 20*log10(abs(f2)));
+d = 20*log10(abs(f2));
sdsq = mean(d.^2);
plot(sdsq)
title('spectral distortion clean and channel SNR=3dB')
y = 1:en-st+1;
x = 1:40;
%mesh(y,x,-20*log10(abs(f2(1:40,:))));
-mesh(y,x,d(1:40,:));
+mesh(y,x,d(x,:));
grid
title('Synthesis filter difference between clean and channel SNR=3dB');
xlabel('Time (codec frames)')
ylabel('Frequency 0 to 2500Hz');
zlabel('Difference (dB)');
+
+% map soft decn information to LSPs
+
+mel1 = ve9qrp_snr3_softdec(:,10:12);
+mel2 = ve9qrp_snr3_softdec(:,13:14);
+mel3 = ve9qrp_snr3_softdec(:,15:18);
+mel4 = ve9qrp_snr3_softdec(:,19:21);
+mel5 = ve9qrp_snr3_softdec(:,22:24);
+mel6 = ve9qrp_snr3_softdec(:,25:26);
+softdec_mel = [sum(mel1'.^2); sum(mel2'.^2); sum(mel3'.^2); sum(mel4'.^2); sum(mel5'.^2); sum(mel6'.^2)];
+
+figure(10)
+clf;
+y = 1:en-st+1;
+x = 1:6;
+%mesh(y,x,-20*log10(abs(f2(1:40,:))));
+mesh(y, x, softdec_mel(:,y));
+grid
+xlabel('Codec frame')
+ylabel('LSP')
+zlabel('Power')
+%axis([1 (en-st+1) 1 6 -10 5])
+
+% plot symbol energy against SD
+
+figure(11)
+semilogx(mean(softdec_mel(:,1:en-st+1)), sdsq,'+')
+grid
% + Only a small change in fading perf with filter on/off
% + however other filters may have other effects, should test this,
% e.g. scatter plots, some sort of BER metric?
-% [ ] pilot based EsNo estimation
+% [X] EsNo estimation
% [ ] filter reqd with compression?
% + make sure not too much noise passed into noise floor
% [X] different diversity combination
% + ssb filter
% + compression
% + make sure it's flat with many errors
-% [ ] linear tracking of ampl ests
-% + seem to only change every
graphics_toolkit ("gnuplot");
more off;
% select which test ----------------------------------------------------------
-%test = 'compare to c';
-test = 'awgn';
+test = 'compare to c';
+%test = 'awgn';
%test = 'fading';
% some parameters that can be over ridden, e.g. to disable parts of modem
% should be BER around 0.015 to 0.02
if strcmp(test, 'awgn')
- frames = 10;
+ frames = 100;
foff = 0;
dfoff = -0/Fs;
- EsNodB = 80;
+ EsNodB = 8;
fading_en = 0;
hf_delay_ms = 2;
compare_with_c = 0;
% Similar to AWGN - should be BER around 0.015 to 0.02
if strcmp(test, 'fading');
- frames = 200;
+ frames = 100;
foff = -10.5;
dfoff = 0.0/Fs;
- EsNodB = 12;
+ EsNodB = 9;
fading_en = 1;
hf_delay_ms = 2;
compare_with_c = 0;
ratio_log = [];
foff_log = [];
f_est_log = [];
+sig_rms_log = [];
+noise_rms_log = [];
+noise_rms_filt_log = [];
% Channel modeling and BER measurement ----------------------------------------
h = freqz(b,a,(600:2600)/(Fs/(2*pi)));
filt_gain = (2600-600)/sum(abs(h) .^ 2); % ensures power after filter == before filter
+noise_rms_filt = 0;
+
% main loop --------------------------------------------------------------------
% run mod and channel as aseparate loop so we can resample to simulate sample rate differences
% if we are in sync complete demodulation with symbol rate processing
if (next_sync == 1) || (sync == 1)
- [rx_symb rx_bits rx_symb_linear amp_ phi_ EsNo_] = qpsk_symbols_to_bits(acohpsk, acohpsk.ct_symb_ff_buf);
+ [rx_symb rx_bits rx_symb_linear amp_ phi_ sig_rms noise_rms] = qpsk_symbols_to_bits(acohpsk, acohpsk.ct_symb_ff_buf);
rx_symb_log = [rx_symb_log; rx_symb];
rx_amp_log = [rx_amp_log; amp_];
rx_phi_log = [rx_phi_log; phi_];
tx_bits_prev_log = [tx_bits_prev_log prev_tx_bits2];
ratio_log = [ratio_log acohpsk.ratio];
ct_symb_ff_log = [ct_symb_ff_log; acohpsk.ct_symb_ff_buf(1:acohpsk.Nsymbrowpilot,:)];
+ sig_rms_log = [sig_rms_log sig_rms];
+ noise_rms_log = [noise_rms_log noise_rms];
+ noise_rms_filt = 0.9*noise_rms_filt + 0.1*noise_rms;
+ noise_rms_filt_log = [noise_rms_filt_log noise_rms_filt];
% BER stats
check(rx_timing_log, rx_timing_log_c, 'rx_timing',0.001);
check(rx_bits_log, rx_bits_log_c, 'rx_bits');
check(f_est_log, f_est_log_c, 'f_est');
+ check(sig_rms_log, sig_rms_log_c, 'sig_rms');
+ check(noise_rms_log, noise_rms_log_c, 'noise_rms');
else
figure(6)
clf;
+ subplot(211)
stem(nerr_log)
title('Bit Errors');
+ subplot(212)
+ plot(noise_rms_filt_log,'r', sig_rms_log,'g');
+ title('Est rms signal and noise')
figure(7);
clf;
plot(error_positions_hist)
title('histogram of bit errors')
+
end
ber = 0.0;
burst_length = burst_period = 0.0;
burst_timer = 0.0;
- dump = natural = 0;
+ dump = natural = softdec = 0;
codec2 = codec2_create(mode);
nsam = codec2_samples_per_frame(codec2);
buf = (short*)malloc(nsam*sizeof(short));
nbyte = (nbit + 7) / 8;
bits = (unsigned char*)malloc(nbyte*sizeof(char));
- softdec_bits = (int*)malloc(nbit*sizeof(float));
+ softdec_bits = (float*)malloc(nbit*sizeof(float));
frames = bit_errors = bits_proc = 0;
nstart_bit = 0;
nend_bit = nbit-1;
codec2_set_natural_or_gray(codec2, !natural);
//printf("%d %d\n", nstart_bit, nend_bit);
+ //fprintf(stderr, "softdec: %d natural: %d\n", softdec, natural);
if (softdec)
ret = (fread(softdec_bits, sizeof(float), nbit, fin) == (size_t)nbit);
else
byte++;
}
}
- }
+ codec2_set_softdec(codec2, softdec_bits);
+ }
codec2_decode_ber(codec2, buf, bits, ber_est);
fwrite(buf, sizeof(short), nsam, fout);
FILE *fout;
short *buf;
unsigned char *bits;
- int nsam, nbit, nbyte, gray, unpacked;
- int *unpacked_bits;
+ int nsam, nbit, nbyte, gray, softdec;
+ float *unpacked_bits;
int bit, byte,i;
if (argc < 4) {
- printf("usage: c2enc 3200|2400|1600|1400|1300|1200|700 InputRawspeechFile OutputBitFile [--natural] [--unpacked]\n");
+ printf("usage: c2enc 3200|2400|1600|1400|1300|1200|700 InputRawspeechFile OutputBitFile [--natural] [--softdec]\n");
printf("e.g c2enc 1400 ../raw/hts1a.raw hts1a.c2\n");
printf("e.g c2enc 1300 ../raw/hts1a.raw hts1a.c2 --natural\n");
exit(1);
nbyte = (nbit + 7) / 8;
bits = (unsigned char*)malloc(nbyte*sizeof(char));
- unpacked_bits = (int*)malloc(nbit*sizeof(int));
+ unpacked_bits = (float*)malloc(nbit*sizeof(float));
- if (argc == 5) {
- if (strcmp(argv[4], "--natural") == 0)
+ gray = 1;
+ softdec = 0;
+ for (i=4; i<argc; i++) {
+ if (strcmp(argv[i], "--natural") == 0) {
gray = 0;
- else
- gray = 1;
- codec2_set_natural_or_gray(codec2, gray);
- if (strcmp(argv[4], "--unpacked") == 0) {
- unpacked = 1;
}
- else
- unpacked = 0;
+ if (strcmp(argv[i], "--softdec") == 0) {
+ softdec = 1;
+ }
}
+ codec2_set_natural_or_gray(codec2, gray);
+ //fprintf(stderr,"gray: %d softdec: %d\n", gray, softdec);
while(fread(buf, sizeof(short), nsam, fin) == (size_t)nsam) {
codec2_encode(codec2, bits, buf);
- if (unpacked) {
- /* unpack bits, MSB first */
+ if (softdec) {
+ /* unpack bits, MSB first, send as soft decision float */
bit = 7; byte = 0;
for(i=0; i<nbit; i++) {
- unpacked_bits[i] = (bits[byte] >> bit) & 0x1;
+ unpacked_bits[i] = 1.0 - 2.0*((bits[byte] >> bit) & 0x1);
bit--;
if (bit < 0) {
bit = 7;
byte++;
}
}
- fwrite(unpacked_bits, sizeof(int), nbit, fout);
+ fwrite(unpacked_bits, sizeof(float), nbit, fout);
}
else
fwrite(bits, sizeof(char), nbyte, fout);
\*---------------------------------------------------------------------------*/
-struct CODEC2 * CODEC2_WIN32SUPPORT codec2_create(int mode)
+struct CODEC2 * codec2_create(int mode)
{
struct CODEC2 *c2;
int i,l;
\*---------------------------------------------------------------------------*/
-void CODEC2_WIN32SUPPORT codec2_destroy(struct CODEC2 *c2)
+void codec2_destroy(struct CODEC2 *c2)
{
assert(c2 != NULL);
free(c2->bpf_buf);
\*---------------------------------------------------------------------------*/
-int CODEC2_WIN32SUPPORT codec2_bits_per_frame(struct CODEC2 *c2) {
+int codec2_bits_per_frame(struct CODEC2 *c2) {
if (c2->mode == CODEC2_MODE_3200)
return 64;
if (c2->mode == CODEC2_MODE_2400)
\*---------------------------------------------------------------------------*/
-int CODEC2_WIN32SUPPORT codec2_samples_per_frame(struct CODEC2 *c2) {
+int codec2_samples_per_frame(struct CODEC2 *c2) {
if (c2->mode == CODEC2_MODE_3200)
return 160;
if (c2->mode == CODEC2_MODE_2400)
return 0; /* shouldnt get here */
}
-void CODEC2_WIN32SUPPORT codec2_encode(struct CODEC2 *c2, unsigned char *bits, short speech[])
+void codec2_encode(struct CODEC2 *c2, unsigned char *bits, short speech[])
{
assert(c2 != NULL);
assert(
codec2_encode_700(c2, bits, speech);
}
-void CODEC2_WIN32SUPPORT codec2_decode(struct CODEC2 *c2, short speech[], const unsigned char *bits)
+void codec2_decode(struct CODEC2 *c2, short speech[], const unsigned char *bits)
{
codec2_decode_ber(c2, speech, bits, 0.0);
}
-void CODEC2_WIN32SUPPORT codec2_decode_ber(struct CODEC2 *c2, short speech[], const unsigned char *bits, float ber_est)
+void codec2_decode_ber(struct CODEC2 *c2, short speech[], const unsigned char *bits, float ber_est)
{
assert(c2 != NULL);
assert(
encode_mels_scalar(indexes, mel, LPC_ORD_LOW);
for(i=0; i<LPC_ORD_LOW; i++) {
- pack(bits, &nbit, indexes[i], mel_bits(i));
+ pack_natural_or_gray(bits, &nbit, indexes[i], mel_bits(i), c2->gray);
}
-
+
pack_natural_or_gray(bits, &nbit, spare, 2, c2->gray);
assert(nbit == (unsigned)codec2_bits_per_frame(c2));
e[3] = decode_energy(e_index, 3);
for(i=0; i<LPC_ORD_LOW; i++) {
- indexes[i] = unpack(bits, &nbit, mel_bits(i));
+ indexes[i] = unpack_natural_or_gray(bits, &nbit, mel_bits(i), c2->gray);
}
+
decode_mels_scalar(mel, indexes, LPC_ORD_LOW);
for(i=0; i<LPC_ORD_LOW; i++) {
f_ = 700.0*( pow(10.0, (float)mel[i]/2595.0) - 1.0);
check_lsp_order(&lsps[3][0], LPC_ORD_LOW);
bw_expand_lsps(&lsps[3][0], LPC_ORD_LOW, 50.0, 100.0);
-
+
+ #ifdef MASK_NOT_FOR_NOW
+ /* first pass at soft decn error masking, needs further work */
+ /* If soft dec info available expand further for low power frames */
+
+ if (c2->softdec) {
+ float e = 0.0;
+ for(i=9; i<9+17; i++)
+ e += c2->softdec[i]*c2->softdec[i];
+ e /= 6.0;
+ //fprintf(stderr, "e: %f\n", e);
+ //if (e < 0.3)
+ // bw_expand_lsps(&lsps[3][0], LPC_ORD_LOW, 150.0, 300.0);
+ }
+ #endif
+
/* interpolate ------------------------------------------------*/
/* LSPs, Wo, and energy are sampled every 40ms so we interpolate
}
}
-void CODEC2_WIN32SUPPORT codec2_set_lpc_post_filter(struct CODEC2 *c2, int enable, int bass_boost, float beta, float gamma)
+void codec2_set_lpc_post_filter(struct CODEC2 *c2, int enable, int bass_boost, float beta, float gamma)
{
assert((beta >= 0.0) && (beta <= 1.0));
assert((gamma >= 0.0) && (gamma <= 1.0));
Experimental method of sending voice/data frames for FreeDV.
*/
-int CODEC2_WIN32SUPPORT codec2_get_spare_bit_index(struct CODEC2 *c2)
+int codec2_get_spare_bit_index(struct CODEC2 *c2)
{
assert(c2 != NULL);
for convenience.
*/
-int CODEC2_WIN32SUPPORT codec2_rebuild_spare_bit(struct CODEC2 *c2, int unpacked_bits[])
+int codec2_rebuild_spare_bit(struct CODEC2 *c2, int unpacked_bits[])
{
int v1,v3;
return -1;
}
-void CODEC2_WIN32SUPPORT codec2_set_natural_or_gray(struct CODEC2 *c2, int gray)
+void codec2_set_natural_or_gray(struct CODEC2 *c2, int gray)
{
assert(c2 != NULL);
c2->gray = gray;
}
+void codec2_set_softdec(struct CODEC2 *c2, float *softdec)
+{
+ assert(c2 != NULL);
+ c2->softdec = softdec;
+}
+
#ifndef __CODEC2__
#define __CODEC2__
-/* set up the calling convention for DLL function import/export for
- WIN32 cross compiling */
-
-#ifdef __CODEC2_WIN32__
-#ifdef __CODEC2_BUILDING_DLL__
-#define CODEC2_WIN32SUPPORT __declspec(dllexport) __stdcall
-#else
-#define CODEC2_WIN32SUPPORT __declspec(dllimport) __stdcall
-#endif
-#else
-#define CODEC2_WIN32SUPPORT
-#endif
-
#define CODEC2_MODE_3200 0
#define CODEC2_MODE_2400 1
#define CODEC2_MODE_1600 2
struct CODEC2;
-struct CODEC2 * CODEC2_WIN32SUPPORT codec2_create(int mode);
-void CODEC2_WIN32SUPPORT codec2_destroy(struct CODEC2 *codec2_state);
-void CODEC2_WIN32SUPPORT codec2_encode(struct CODEC2 *codec2_state, unsigned char * bits, short speech_in[]);
-void CODEC2_WIN32SUPPORT codec2_decode(struct CODEC2 *codec2_state, short speech_out[], const unsigned char *bits);
-void CODEC2_WIN32SUPPORT codec2_decode_ber(struct CODEC2 *codec2_state, short speech_out[], const unsigned char *bits, float ber_est);
-int CODEC2_WIN32SUPPORT codec2_samples_per_frame(struct CODEC2 *codec2_state);
-int CODEC2_WIN32SUPPORT codec2_bits_per_frame(struct CODEC2 *codec2_state);
-
-void CODEC2_WIN32SUPPORT codec2_set_lpc_post_filter(struct CODEC2 *codec2_state, int enable, int bass_boost, float beta, float gamma);
-int CODEC2_WIN32SUPPORT codec2_get_spare_bit_index(struct CODEC2 *codec2_state);
-int CODEC2_WIN32SUPPORT codec2_rebuild_spare_bit(struct CODEC2 *codec2_state, int unpacked_bits[]);
-void CODEC2_WIN32SUPPORT codec2_set_natural_or_gray(struct CODEC2 *codec2_state, int gray);
+struct CODEC2 * codec2_create(int mode);
+void codec2_destroy(struct CODEC2 *codec2_state);
+void codec2_encode(struct CODEC2 *codec2_state, unsigned char * bits, short speech_in[]);
+void codec2_decode(struct CODEC2 *codec2_state, short speech_out[], const unsigned char *bits);
+void codec2_decode_ber(struct CODEC2 *codec2_state, short speech_out[], const unsigned char *bits, float ber_est);
+int codec2_samples_per_frame(struct CODEC2 *codec2_state);
+int codec2_bits_per_frame(struct CODEC2 *codec2_state);
+
+void codec2_set_lpc_post_filter(struct CODEC2 *codec2_state, int enable, int bass_boost, float beta, float gamma);
+int codec2_get_spare_bit_index(struct CODEC2 *codec2_state);
+int codec2_rebuild_spare_bit(struct CODEC2 *codec2_state, int unpacked_bits[]);
+void codec2_set_natural_or_gray(struct CODEC2 *codec2_state, int gray);
+void codec2_set_softdec(struct CODEC2 *c2, float *softdec);
#endif
coh->fdmdv = fdmdv;
+ coh->sig_rms = coh->noise_rms = 0.0;
+
/* disable optional logging by default */
coh->rx_baseband_log = NULL;
void qpsk_symbols_to_bits(struct COHPSK *coh, float rx_bits[], COMP ct_symb_buf[][COHPSK_NC*ND])
{
- int p, r, c, i, pc, d;
+ int p, r, c, i, pc, d, n;
float x[NPILOTSFRAME+2], x1;
COMP y[NPILOTSFRAME+2], yfit;
+ COMP rx_symb_linear[NSYMROW*COHPSK_NC*ND];
COMP m, b;
COMP __attribute__((unused)) corr, rot, pi_on_4, phi_rect, div_symb;
float mag, __attribute__((unused)) phi_, __attribute__((unused)) amp_;
+ float sum_x, sum_xx, noise_var;
+ COMP s;
pi_on_4.real = cosf(M_PI/4); pi_on_4.imag = sinf(M_PI/4);
for(c=0; c<COHPSK_NC*ND; c++) {
-//#define AVERAGE
-#ifdef AVERAGE
- /* Average pilots to get phase and amplitude estimates using
- two pilots at the start of each frame and two at the end */
-
- corr.real = 0.0; corr.imag = 0.0; mag = 0.0;
- for(p=0; p<NPILOTSFRAME+2; p++) {
- corr = cadd(corr, fcmult(coh->pilot2[p][c], ct_symb_buf[sampling_points[p]][c]));
- mag += cabsolute(ct_symb_buf[sampling_points[p]][c]);
- }
-
- phi_ = atan2f(corr.imag, corr.real);
- amp_ = mag/(NPILOTSFRAME+2);
- for(r=0; r<NSYMROW; r++) {
- coh->phi_[r][c] = phi_;
- coh->amp_[r][c] = amp_;
- }
-#else
/* Set up lin reg model and interpolate phase. Works better than average for channels with
quickly changing phase, like HF. */
x[p] = sampling_points[p];
pc = c % COHPSK_NC;
y[p] = fcmult(coh->pilot2[p][pc], ct_symb_buf[sampling_points[p]][c]);
- //printf("%f %f\n", y[p].real, y[p].imag);
}
- //printf("\n");
+
linreg(&m, &b, x, y, NPILOTSFRAME+2);
for(r=0; r<NSYMROW; r++) {
x1 = (float)(r+NPILOTSFRAME);
yfit = cadd(fcmult(x1,m),b);
coh->phi_[r][c] = atan2(yfit.imag, yfit.real);
- //printf(" %f", coh->phi_[r][c]);
}
- //printf("\n");
+
/* amplitude estimation */
mag = 0.0;
for(r=0; r<NSYMROW; r++) {
coh->amp_[r][c] = amp_;
}
-#endif
}
- //exit(0);
+
/* now correct phase of data symbols */
for(c=0; c<COHPSK_NC*ND; c++) {
- //rot.real = 1.0; rot.imag = 0.0;
for (r=0; r<NSYMROW; r++) {
phi_rect.real = cosf(coh->phi_[r][c]); phi_rect.imag = -sinf(coh->phi_[r][c]);
coh->rx_symb[r][c] = cmult(ct_symb_buf[NPILOTSFRAME + r][c], phi_rect);
- //printf("%d %d %f %f\n", r,c, coh->rx_symb[r][c].real, coh->rx_symb[r][c].imag);
- //printf("phi_ %d %d %f %f\n", r,c, ct_symb_buf[NPILOTSFRAME + r][c].real, ct_symb_buf[NPILOTSFRAME + r][c].imag);
+ i = c*NSYMROW + r;
+ rx_symb_linear[i] = coh->rx_symb[r][c];
}
}
- /* and finally optional diversity combination and make decn on bits */
+ /* and finally optional diversity combination, note output is soft decm a "1" is < 0 */
for(c=0; c<COHPSK_NC; c++) {
for(r=0; r<NSYMROW; r++) {
rx_bits[2*i] = rot.imag;
}
}
+
+
+ /* estimate RMS signal and noise */
+
+ mag = 0.0;
+ for(i=0; i<NSYMROW*COHPSK_NC*ND; i++)
+ mag += cabsolute(rx_symb_linear[i]);
+ coh->sig_rms = mag/(NSYMROW*COHPSK_NC*ND);
+
+ sum_x = 0;
+ sum_xx = 0;
+ n = 0;
+ for (i=0; i<NSYMROW*COHPSK_NC*ND; i++) {
+ s = rx_symb_linear[i];
+ if (fabsf(s.real) > coh->sig_rms) {
+ sum_x += s.imag;
+ sum_xx += s.imag*s.imag;
+ n++;
+ }
+ }
+
+ noise_var = 0;
+ if (n > 1) {
+ noise_var = (n*sum_xx - sum_x*sum_x)/(n*(n-1));
+ }
+ coh->noise_rms = sqrtf(noise_var);
+
}
int frame;
float ratio;
+ float sig_rms;
+ float noise_rms;
+
struct FDMDV *fdmdv;
+
/* optional log variables used for testing Octave to C port */
COMP *rx_baseband_log;
AUTHOR......: David Rowe
DATE CREATED: April 5 2015
- Given an input file of bits (note one bit per int, not compressed),
+ Given an input file of bits (note one bit per float, soft decision format),
outputs a raw file (8kHz, 16 bit shorts) of COHPSK modem samples
ready to send over a HF radio channel.
{
FILE *fin, *fout;
struct COHPSK *cohpsk;
+ float tx_bits_sd[COHPSK_BITS_PER_FRAME];
int tx_bits[COHPSK_BITS_PER_FRAME];
COMP tx_fdm[COHPSK_SAMPLES_PER_FRAME];
short tx_fdm_scaled[COHPSK_SAMPLES_PER_FRAME];
frames = 0;
- while(fread(tx_bits, sizeof(int), COHPSK_BITS_PER_FRAME, fin) == COHPSK_BITS_PER_FRAME) {
+ while(fread(tx_bits_sd, sizeof(float), COHPSK_BITS_PER_FRAME, fin) == COHPSK_BITS_PER_FRAME) {
frames++;
+ for(i=0; i<COHPSK_BITS_PER_FRAME; i++)
+ tx_bits[i] = tx_bits_sd[i] < 0.0;
cohpsk_mod(cohpsk, tx_fdm, tx_bits);
cohpsk_clip(tx_fdm);
int rx_bits_log[COHPSK_BITS_PER_FRAME*FRAMES];
FILE *fout;
- int f, r, c, log_r, log_data_r, noise_r, ff_log_r;
+ int f, r, c, log_r, log_data_r, noise_r, ff_log_r, i;
int *ptest_bits_coh, *ptest_bits_coh_end;
double foff;
COMP foff_rect, phase_ch;
//COMP rx_onesym[COHPSK_NC*ND];
//int rx_baseband_log_col_index = 0;
//COMP rx_baseband_log[COHPSK_NC*ND][(M+M/P)*NSYMROWPILOT*FRAMES];
- float f_est_log[FRAMES];
+ float f_est_log[FRAMES], sig_rms_log[FRAMES], noise_rms_log[FRAMES];
int f_est_samples;
int log_bits;
tmp = nin_frame;
cohpsk_demod(coh, rx_bits_sd, &reliable_sync_bit, &ch_fdm_frame_log_out[ch_fdm_frame_log_index], &nin_frame);
for(i=0; i<COHPSK_BITS_PER_FRAME; i++)
- rx_bits[i] = rx_bits_sd[i] > 0.0;
+ rx_bits[i] = rx_bits_sd[i] < 0.0;
ch_fdm_frame_log_index += tmp;
}
memcpy(&rx_bits_log[COHPSK_BITS_PER_FRAME*log_bits], rx_bits, sizeof(int)*COHPSK_BITS_PER_FRAME);
log_bits++;
- f_est_log[f_est_samples++] = coh->f_est;
+ f_est_log[f_est_samples] = coh->f_est;
+ sig_rms_log[f_est_samples] = coh->sig_rms;
+ noise_rms_log[f_est_samples] = coh->noise_rms;
+ f_est_samples++;;
}
assert(log_r <= NSYMROWPILOT*FRAMES);
octave_save_complex(fout, "rx_symb_log_c", (COMP*)rx_symb_log, log_data_r, COHPSK_NC*ND, COHPSK_NC*ND);
octave_save_int(fout, "rx_bits_log_c", rx_bits_log, 1, COHPSK_BITS_PER_FRAME*log_bits);
octave_save_float(fout, "f_est_log_c", &f_est_log[1], 1, f_est_samples-1, f_est_samples-1);
+ octave_save_float(fout, "sig_rms_log_c", sig_rms_log, 1, f_est_samples, f_est_samples-1);
+ octave_save_float(fout, "noise_rms_log_c", noise_rms_log, 1, f_est_samples, f_est_samples);
#ifdef XX
#endif
fclose(fout);
tcohpsk first (any variant) to load the function into Octave, e.g.:
octave:17> tcohpsk
- octave:18> write_noise_file("../raw/fading_samples.float", 7500, 7500*60)
+ octave:18> write_noise_file("../raw/fast_fading_samples.float", 7500, 7500*60)
+ octave:19> write_noise_file("../raw/slow_fading_samples.float", 75000, 7500*60)
*/
#define FADING_FILE_NAME "../../raw/fading_samples.float"
int ch_buf_n;
float tx_pwr, rx_pwr, noise_pwr;
int error_positions_hist[COHPSK_BITS_PER_FRAME];
- int log_data_r, c, j, tmp, ret;
+ int log_data_r, c, j, tmp;
for(i=0; i<COHPSK_BITS_PER_FRAME; i++)
error_positions_hist[i] = 0;
projects / freetel-svn-tracking.git / commitdiff