phase model uses fractional pitch, error measure matches Octave plamp.m for consistency

git-svn-id: https://svn.code.sf.net/p/freetel/code@52 01035d8c-6547-0410-b346-abe4f91aad63

diff --git a/codec2/octave/phase.m b/codec2/octave/phase.m
index 537e193b1db15f0b1c22a9c2e2fb237752216971..fc13b47bf54a19c30e8ae2fd6b3ad142e05a1359 100644 (file)
--- a/codec2/octave/phase.m
+++ b/codec2/octave/phase.m
@@ -10,7 +10,7 @@ function phase(samname, F0, randphase, png)
   A = 10000/L;
   phi = zeros(1,L);
   for m=1:L
-    phi(m) = m*Wo*0.3*m;
+    phi(m) = m*Wo*50.125;
   end
   if (randphase == 1) 
     rand("seed",0);

diff --git a/codec2/octave/plamp.m b/codec2/octave/plamp.m
index cf7618cef12e6725a3ac594f132b935e1913f1c4..cdf6176ff99318da572c87ec316ec52d881a5e41 100644 (file)
--- a/codec2/octave/plamp.m
+++ b/codec2/octave/plamp.m
@@ -45,6 +45,11 @@ function plamp(samname, f)
     phase_ = load(phase_name_);
   endif
 
+  snr_name = strcat(samname,"_snr.txt");
+  if (file_in_path(".",snr_name))
+    snr = load(snr_name);
+  endif
+
   do 
     figure(1);
     clg;
@@ -74,16 +79,21 @@ function plamp(samname, f)
       plot((1:L)*Wo*4000/pi, 20*log10(Amq) - 20*log10(Am), Am_err_label);
     endif
 
-    % phase model - determine SNR and error spectrum
+    if (file_in_path(".",snr_name))
+      snr_label = sprintf(";phase SNR %4.2f dB;",snr(f));
+      plot(1,1,snr_label);
+    endif
+
+    % phase model - determine SNR and error spectrum for phase model 1
 
     if (file_in_path(".",phase_name_))
       orig  = Am.*exp(j*phase(f,1:L));
       synth = Am.*exp(j*phase_(f,1:L));
-      signal = orig * orig'
-      noise = (orig-synth) * (orig-synth)'
-      snr = 10*log10(signal/noise);
+      signal = orig * orig';
+      noise = (orig-synth) * (orig-synth)';
+      snr_phase = 10*log10(signal/noise);
 
-      phase_err_label = sprintf(";phase_err SNR %4.2f dB;",snr);
+      phase_err_label = sprintf(";phase_err SNR %4.2f dB;",snr_phase);
       plot((1:L)*Wo*4000/pi, 20*log10(orig-synth), phase_err_label);
     endif
 
@@ -97,16 +107,28 @@ function plamp(samname, f)
 
     if (file_in_path(".",phase_name))
       figure(3);
-      plot((1:L)*Wo*4000/pi, unwrap(phase(f,1:L)), ";phase;");
+      plot((1:L)*Wo*4000/pi, phase(f,1:L), ";phase;");
       axis;
       if (file_in_path(".",phase_name_))
         hold on;
-        plot((1:L)*Wo*4000/pi, unwrap(phase_(f,1:L)), ";phase_;");
+        plot((1:L)*Wo*4000/pi, phase_(f,1:L), ";phase_;");
        hold off;
       endif
       figure(2);
     endif
 
+    % autocorrelation function to research voicing est
+    
+    %M = length(s);
+    %sw = s .* hanning(M)';
+    %for k=0:159
+    %  R(k+1) = sw(1:320-k) * sw(1+k:320)';
+    %endfor
+    %figure(4);
+    %R_label = sprintf(";R(k) %3.2f;",max(R(20:159))/R(1));
+    %plot(R/R(1),R_label);
+    %grid
+
     % interactive menu
 
     printf("\rframe: %d  menu: n-next  b-back  p-png  q-quit ", f);

diff --git a/codec2/src/dump.c b/codec2/src/dump.c
index 02a354d57419187babf8523ccfaffbaab5341204..c52d574b955537b9bc47d7c20a9e2aeea8c2a735 100644 (file)
--- a/codec2/src/dump.c
+++ b/codec2/src/dump.c
@@ -45,6 +45,7 @@ static FILE *ffw = NULL;
 static FILE *fe = NULL;
 static FILE *fsq = NULL;
 static FILE *fdec = NULL;
+static FILE *fsnr = NULL;
 
 static char  prefix[MAX_STR];
 
@@ -80,6 +81,8 @@ void dump_off(){
        fclose(fsq);
     if (fdec != NULL)
        fclose(fdec);
+    if (fsnr != NULL)
+       fclose(fsnr);
 }
 
 void dump_Sn(float Sn[]) {
@@ -223,6 +226,20 @@ void dump_phase_(float phase_[]) {
     fprintf(fphase_,"\n");    
 }
 
+void dump_snr(float snr) {
+    char s[MAX_STR];
+
+    if (!dumpon) return;
+
+    if (fphase_ == NULL) {
+       sprintf(s,"%s_snr.txt", prefix);
+       fsnr = fopen(s, "wt");
+       assert(fsnr != NULL);
+    }
+
+    fprintf(fsnr,"%f\n",snr);
+}
+
 void dump_Pw(COMP Pw[]) {
     int i;
     char s[MAX_STR];

diff --git a/codec2/src/dump.h b/codec2/src/dump.h
index 031445af6e7eac6918fc54d089cc43508a7e961e..2d3bbc25b20539e4ab440745c5ce0fc566e4ef7a 100644 (file)
--- a/codec2/src/dump.h
+++ b/codec2/src/dump.h
@@ -38,6 +38,10 @@ void dump_model(MODEL *m);
 void dump_quantised_model(MODEL *m);
 void dump_Pw(COMP Pw[]);
 void dump_lsp(float lsp[]);
+
+/* phase modelling */
+
+void dump_snr(float snr);
 void dump_phase(float phase[]);
 void dump_phase_(float phase[]);
 

diff --git a/codec2/src/phase.c b/codec2/src/phase.c
index 42b1c658285d98fc2c1312f22a482aaa19392d7b..ffb132af76b501648f50aad44a25c6d34fbd05ea 100644 (file)
--- a/codec2/src/phase.c
+++ b/codec2/src/phase.c
@@ -126,10 +126,10 @@ int    order;
 \*---------------------------------------------------------------------------*/
 
 float phase_model_first_order(
-  float aks[],                  /* LPC coeffs for this frame      */
-  COMP  H[],                   /* LPC filter freq doamin samples */
-  int  *i_min,                  /* pulse position for min error   */ 
-  COMP *minAm                   /* complex gain for min error     */
+  float  aks[],                  /* LPC coeffs for this frame      */
+  COMP   H[],                   /* LPC filter freq doamin samples */
+  float *n_min,                  /* pulse position for min error   */ 
+  COMP  *minAm                   /* complex gain for min error     */
 ) 
 {
   float G;                     /* LPC gain */
@@ -144,11 +144,15 @@ float phase_model_first_order(
   COMP Em;                     /* error for m-th band */
   float den;                   /* energy of synthesised */
   float snr;                   /* snr of each excitation source */
+  int   Lmax;
+  float n;
+
+  Lmax = model.L;
 
   /* Construct target vector */
 
   sig = 0.0;
-  for(m=1; m<=model.L; m++) {
+  for(m=1; m<=Lmax; m++) {
     A[m].real = model.A[m]*cos(model.phi[m]);
     A[m].imag = model.A[m]*sin(model.phi[m]);
     sig += model.A[m]*model.A[m];
@@ -156,6 +160,11 @@ float phase_model_first_order(
 
   /* Sample LPC model at harmonics */
 
+  #ifdef NO_LPC_PHASE
+  /* useful for testing with Sn[] an impulse train */
+  for(i=1; i<=PHASE_LPC_ORD; i++)
+     aks[i] = 0;
+  #endif
   G = 1.0;
   aks_to_H(&model,aks,G,H,PHASE_LPC_ORD);
 
@@ -163,16 +172,16 @@ float phase_model_first_order(
 
   Emin = 1E32;
   P = floor(TWO_PI/model.Wo + 0.5);
-  for(i=0; i<P; i++) {
+  for(n=0; n<P; n+=0.25) {
 
     /* determine complex gain */
 
     Am.real = 0.0;
     Am.imag = 0.0;
     den = 0.0;
-    for(m=1; m<=model.L; m++) {
-      Ex[m].real = cos(model.Wo*m*i);
-      Ex[m].imag = sin(-model.Wo*m*i);
+    for(m=1; m<=Lmax; m++) {
+      Ex[m].real = cos(model.Wo*m*n);
+      Ex[m].imag = sin(-model.Wo*m*n);
       A_[m].real = H[m].real*Ex[m].real - H[m].imag*Ex[m].imag;
       A_[m].imag = H[m].imag*Ex[m].real + H[m].real*Ex[m].imag;
       Am.real += A[m].real*A_[m].real + A[m].imag*A_[m].imag;
@@ -186,17 +195,20 @@ float phase_model_first_order(
     /* determine error */
 
     E = 0.0;
-    for(m=1; m<=model.L; m++) {
+    for(m=1; m<=Lmax; m++) {
+      float new_phi;
 
       Em.real = A_[m].real*Am.real - A_[m].imag*Am.imag;
       Em.imag = A_[m].imag*Am.real + A_[m].real*Am.imag;
          
-      E += pow(A[m].real-Em.real,2.0) + pow(A[m].imag-Em.imag,2.0);
+      new_phi = atan2(Em.imag, Em.real+1E-12);
+      E += pow(model.A[m]*(cos(model.phi[m])-cos(new_phi)),2.0);
+      E += pow(model.A[m]*(sin(model.phi[m])-sin(new_phi)),2.0);
     }
 
     if (E < Emin) {
       Emin = E;
-      *i_min = i;
+      *n_min = n;
       minAm->real = Am.real;
       minAm->imag = Am.imag;
     }
@@ -370,7 +382,7 @@ void phase_synth_zero_order(
 void phase_synth_first_order(
   float snr,     /* SNR from first order model */
   COMP  H[],     /* LPC spectra samples        */
-  int   i_min,   /* best pulse position        */
+  float n_min,   /* best pulse position        */
   COMP  minAm    /* best complex gain          */
 )
 {
@@ -397,8 +409,8 @@ void phase_synth_first_order(
     /* generate excitation */
 
     if (m <= Lrand) {
-       Ex[m].real = cos(model.Wo*m*i_min);
-       Ex[m].imag = sin(-model.Wo*m*i_min);
+       Ex[m].real = cos(model.Wo*m*n_min);
+       Ex[m].imag = sin(-model.Wo*m*n_min);
     }
     else {
       float phi = TWO_PI*(float)rand()/RAND_MAX;

diff --git a/codec2/src/phase.h b/codec2/src/phase.h
index e692f93835047643755e7e00a864f0fe32fc8e92..5dd4ad69242136717fe635d73bb940b2d05504bb 100644 (file)
--- a/codec2/src/phase.h
+++ b/codec2/src/phase.h
@@ -31,8 +31,8 @@
 
 #define PHASE_LPC_ORD 10
 
-float phase_model_first_order(float aks[], COMP H[], int *i_min, COMP *min_Am);
+float phase_model_first_order(float aks[], COMP H[], float *n_min, COMP *min_Am);
 void phase_synth_zero_order(float snr, COMP H[], float *prev_Wo, float *ex_phase);
-void phase_synth_first_order(float snr, COMP H[], int i_min, COMP min_Am);
+void phase_synth_first_order(float snr, COMP H[], float n_min, COMP min_Am);
 
 #endif

diff --git a/codec2/src/sinedec.c b/codec2/src/sinedec.c
index 54083330f8149bb8e4f0ba0784427ea2828a5f4c..9d72a274a18759084df6a197052dd5709a3fce34 100644 (file)
--- a/codec2/src/sinedec.c
+++ b/codec2/src/sinedec.c
@@ -200,7 +200,7 @@ int main(int argc, char *argv[])
        float Wn[AW_ENC];               /* windowed speech samples */
        float Rk[PHASE_LPC_ORD+1];      /* autocorrelation coeffs  */
         COMP  H[MAX_AMP];               /* LPC freq domain samples */
-       int   i_min;
+       float n_min;
        COMP  min_Am;
        
        dump_phase(&model.phi[0]);
@@ -217,14 +217,16 @@ int main(int argc, char *argv[])
        else
            assert(order == PHASE_LPC_ORD);
 
-       snr = phase_model_first_order(ak, H, &i_min, &min_Am);
-       snr = 5;
-       if (phase_model == 0)
+       snr = phase_model_first_order(ak, H, &n_min, &min_Am);
+       dump_snr(snr);
+       if (phase_model == 0) {
            phase_synth_zero_order(snr, H, &prev_Wo, &ex_phase);
-       if (phase_model == 1)
-           phase_synth_first_order(snr, H, i_min, min_Am);
-       
-        dump_phase_(&model.phi[0]);
+       }
+
+       if (phase_model == 1) {
+           phase_synth_first_order(snr, H, n_min, min_Am);
+            dump_phase_(&model.phi[0]);
+        }
     }
 
     /* Synthesise speech */