summaryrefslogtreecommitdiff
path: root/cesar/cp/spoc
diff options
context:
space:
mode:
Diffstat (limited to 'cesar/cp/spoc')
-rw-r--r--cesar/cp/spoc/Module1
-rw-r--r--cesar/cp/spoc/spoc.h50
-rw-r--r--cesar/cp/spoc/src/spoc.c175
-rw-r--r--cesar/cp/spoc/test/ApproxOrder1ParamMCoeff.sci36
-rw-r--r--cesar/cp/spoc/test/ApproxOrder2ParamMCoeff.sci34
-rwxr-xr-xcesar/cp/spoc/test/Makefile11
-rw-r--r--cesar/cp/spoc/test/coeff_generation.scilab22
-rw-r--r--cesar/cp/spoc/test/define_spoc_c.h20
-rw-r--r--cesar/cp/spoc/test/get_Approx1.sci21
-rw-r--r--cesar/cp/spoc/test/src/test_coeff_check.c194
10 files changed, 564 insertions, 0 deletions
diff --git a/cesar/cp/spoc/Module b/cesar/cp/spoc/Module
new file mode 100644
index 0000000000..3eadb9ee3f
--- /dev/null
+++ b/cesar/cp/spoc/Module
@@ -0,0 +1 @@
+SOURCES := spoc.c
diff --git a/cesar/cp/spoc/spoc.h b/cesar/cp/spoc/spoc.h
new file mode 100644
index 0000000000..36cfca523e
--- /dev/null
+++ b/cesar/cp/spoc/spoc.h
@@ -0,0 +1,50 @@
+#ifndef spoc_h
+#define spoc_h
+/* Cesar project {{{
+ *
+ * Copyright (C) 2008 Spidcom
+ *
+ * <<<Licence>>>
+ *
+ * }}} */
+/**
+ * \file spoc.h
+ * \brief « brief description »
+ * \ingroup « module »
+ *
+ * « long description »
+ */
+
+//typedef float spoc_prec_t;
+typedef double spoc_prec_t; // en flottant sinx/x avec x=0.1ppm retourne 1
+
+
+#ifdef SPOC_COMPARE_SCILAB
+spoc_prec_t args_haut[41];
+spoc_prec_t args_bas[41];
+spoc_prec_t sinc_haut[41];
+spoc_prec_t sinc_bas[41];
+spoc_prec_t divisor[41];
+spoc_prec_t args_haut_rx[41];
+spoc_prec_t args_bas_rx[41];
+spoc_prec_t sinc_haut_rx[41];
+spoc_prec_t sinc_bas_rx[41];
+#endif
+spoc_prec_t reg_pente_tx[41];
+spoc_prec_t reg_ordo_tx[41];
+spoc_prec_t reg_pente_rx[41];
+spoc_prec_t reg_ordo_rx[41];
+
+spoc_prec_t reg_CG_N[21];
+spoc_prec_t reg_CG_FC10[21];
+spoc_prec_t reg_CG_PR1[21];
+spoc_prec_t reg_CG_PR2[21];
+
+
+void
+spoc_MCoeff (spoc_prec_t rho);
+
+void
+spoc_CG (spoc_prec_t rho, int n, spoc_prec_t *reg);
+
+#endif /* spoc_h */
diff --git a/cesar/cp/spoc/src/spoc.c b/cesar/cp/spoc/src/spoc.c
new file mode 100644
index 0000000000..ef258570a0
--- /dev/null
+++ b/cesar/cp/spoc/src/spoc.c
@@ -0,0 +1,175 @@
+/* Cesar project {{{
+ *
+ * Copyright (C) 2008 Spidcom
+ *
+ * <<<Licence>>>
+ *
+ * }}} */
+/**
+ * \file src/spoc.c
+ * \brief « brief description »
+ * \ingroup « module »
+ *
+ * « long description »
+ */
+#include "common/std.h"
+#include "hal/phy/defs.h"
+#include "cp/spoc/spoc.h"
+#include "math.h"
+
+#define SPOC_PMIN PHY_CARRIER_OFFSET
+#define SPOC_PMAX (PHY_CARRIER_NB + SPOC_PMIN - 1)
+#define SPOC_DIAGONAL_LENGTH PHY_CARRIER_NB
+
+spoc_prec_t kaiser[10] = { 0.9905, 0.9622, 0.9164, 0.8546, 0.7792, 0.6929, 0.5989, 0.5005, 0.4011, 0.304 };
+
+void
+spoc_CG (spoc_prec_t rho, int N, spoc_prec_t *reg)
+{
+ spoc_prec_t arg_p = N*M_PI*rho;
+ spoc_prec_t arg_n = arg_p;
+ spoc_prec_t sinus = sin (arg_p);
+ if (arg_p == 0) reg[10] = 1;
+ else reg[10] = sinus/arg_p;
+ int k;
+ spoc_prec_t kk;
+ for(k=1;k<=10; k++)
+ {
+ sinus = -sinus;
+ kk = kaiser[k-1]*sinus;
+ arg_p+=M_PI;
+ arg_n-=M_PI;
+ /* Remarque sur la division par arg_p = pi(k-N*rho).
+ * si |rho|<300ppm, N<3072, |k|<=1 alors arg_p != 0
+ * et lalors arg_p=0 que si rho=0 et k=0. voir au dessus. */
+ reg[k+10] = kk/arg_p;
+ reg[10-k] = kk/arg_n;
+ }
+}
+
+void
+spoc_MCoeff (spoc_prec_t rho)
+{
+ int k;
+ spoc_prec_t xi,xe;
+ spoc_prec_t sin_xi,sin_xe;
+ spoc_prec_t cos_xi,cos_xe;
+ spoc_prec_t pirho = M_PI * rho;
+ spoc_prec_t SIN_PIRHO = sin(pirho);
+ spoc_prec_t COS_PIRHO = cos(pirho);
+ spoc_prec_t delta;
+ uint len;
+ spoc_prec_t sin_prev;
+
+ spoc_prec_t rho_abs = rho;
+ if (rho < 0) rho_abs = -rho;
+ if (rho_abs >= 1e-8)
+ {
+ //printf ("normal computation\n");
+ xi = pirho * SPOC_PMIN;
+ xe = pirho * SPOC_PMAX;
+ sin_xi = sin(xi);
+ sin_xe = sin(xe);
+ cos_xi = cos(xi);
+ len = SPOC_DIAGONAL_LENGTH - 1; // interval nb in the diagonal.
+
+#ifdef SPOC_COMPARE_SCILAB
+ args_haut[20] = xi;
+ args_bas[20] = xe;
+ sinc_haut[20] = sin_xi/xi;
+ sinc_bas[20] = sin_xe/xe;
+ divisor[20] = len*len;
+ args_haut_rx[20] = -xi;
+ args_bas_rx[20] = -xe;
+ sinc_haut_rx[20] = sin_xi/xi;
+ sinc_bas_rx[20] = sin_xe/xe;
+#endif
+ if (xi == 0 )
+ {
+ reg_ordo_tx[20] = reg_ordo_rx[20] = 1;
+ reg_pente_tx[20] = reg_pente_rx[20] = 0;
+ }
+ else
+ {
+ reg_ordo_tx[20] = reg_ordo_rx[20] = sin_xi/xi;
+ reg_pente_tx[20] = reg_pente_rx[20] = (sin_xe/xe - reg_ordo_tx[20])/(len*len);
+ }
+
+ delta = 0;
+ /* Diagonales superieures */
+ for (k=1; k<=20; k++)
+ {
+ xi += (M_PI + pirho);
+ xe += M_PI;
+ delta += 2*M_PI;
+ sin_prev = sin_xi;
+ sin_xi = -sin_prev * COS_PIRHO - cos_xi * SIN_PIRHO;
+ cos_xi = sin_prev * SIN_PIRHO - cos_xi * COS_PIRHO;
+ sin_xe = -sin_xe;
+ sin_xi = sin(xi);
+ len -= 1;
+ spoc_prec_t mux = 1.0/len;
+ reg_ordo_tx[k+20] = sin_xi/xi;
+ reg_ordo_rx[k+20] = -sin_xi/(-xi+delta);
+ reg_pente_tx[k+20] = (sin_xe/xe - reg_ordo_tx[k+20])*mux;
+ reg_pente_rx[k+20] = (-sin_xe/(-xe+delta) - reg_ordo_rx[k+20])*mux;
+#ifdef SPOC_COMPARE_SCILAB
+ divisor[k+20] = len;
+ args_haut[k+20] = xi;
+ args_bas[k+20] = xe;
+ sinc_haut[k+20] = sin_xi/xi;
+ sinc_bas[k+20] = sin_xe/xe;
+ args_haut_rx[k+20] = -xi+delta;
+ args_bas_rx[k+20] = -xe+delta;
+ sinc_haut_rx[k+20] = -sin_xi/(-xi+delta);
+ sinc_bas_rx[k+20] = -sin_xe/(-xe+delta);
+#endif
+ }
+ xi = pirho * SPOC_PMIN;
+ xe = pirho * SPOC_PMAX;
+ sin_xi = reg_ordo_tx[20] * xi;
+ sin_xe = sin(xe);
+ cos_xe = cos(xe);
+ len = SPOC_DIAGONAL_LENGTH - 1;
+ delta=0;
+ /* Diagonales inferieures */
+ for (k=1; k<=20; k++)
+ {
+ xe -= (M_PI + pirho);
+ xi -= M_PI;
+ delta -= 2*M_PI;
+ sin_prev = sin_xe;
+ sin_xe = -sin_prev * COS_PIRHO + cos_xe * SIN_PIRHO;
+ cos_xe = -sin_prev * SIN_PIRHO - cos_xe * COS_PIRHO;
+ sin_xi = -sin_xi;
+ len -= 1;
+ spoc_prec_t mux = 1.0/len;
+ reg_ordo_tx[20-k] = sin_xi/xi;
+ reg_ordo_rx[20-k] = -sin_xi/(-xi+delta);
+ reg_pente_tx[20-k] = (sin_xe/xe - reg_ordo_tx[20-k])*mux;
+ reg_pente_rx[20-k] = (-sin_xe/(-xe+delta) - reg_ordo_rx[20-k])*mux;
+#ifdef SPOC_COMPARE_SCILAB
+ divisor[20-k] = len;
+ args_haut[20-k] = xi;
+ args_bas[20-k] = xe;
+ sinc_haut[20-k] = sin_xi/xi;
+ sinc_bas[20-k] = sin_xe/xe;
+ args_haut_rx[20-k] = -xi+delta;
+ args_bas_rx[20-k] = -xe+delta;
+ sinc_haut_rx[20-k] = -sin_xi/(-xi+delta);
+ sinc_bas_rx[20-k] = -sin_xe/(-xe+delta);
+#endif
+ }
+ }
+ else
+ {
+ //printf ("fast computation\n");
+ for (k=1; k<=20; k++)
+ {
+ reg_ordo_tx[k+20] = reg_ordo_tx[20-k] = reg_ordo_rx[k+20] = reg_ordo_rx [20-k] = 0;
+ reg_pente_tx[k+20] = reg_pente_tx[20-k]= reg_pente_rx[k+20] = reg_pente_rx [20-k] = 0;
+ }
+ reg_ordo_tx[20] = reg_ordo_rx[20] = 1;
+ reg_pente_tx[20] = reg_pente_rx[20] = 0;
+ }
+}
diff --git a/cesar/cp/spoc/test/ApproxOrder1ParamMCoeff.sci b/cesar/cp/spoc/test/ApproxOrder1ParamMCoeff.sci
new file mode 100644
index 0000000000..b76f650c69
--- /dev/null
+++ b/cesar/cp/spoc/test/ApproxOrder1ParamMCoeff.sci
@@ -0,0 +1,36 @@
+function [arg_i, arg_e, sinc_i, sinc_e, div, ordo, pente] = ApproxOrder1ParamMCoeff(rho,X0,N,Coeff_number)
+
+if Coeff_number<0
+ xi=X0 - Coeff_number;
+ xe=N-1;
+ qi= X0-Coeff_number;
+ qe=N-1;
+ ki=X0+ 0;
+ ke=N+Coeff_number-1;
+elseif Coeff_number>0
+ xi=X0;
+ xe=(N-1)-abs(Coeff_number);
+ qi=X0+ 0;
+ qe=N-Coeff_number-1;
+ ki=X0+ Coeff_number;
+ ke=N-1;
+else
+ xi=X0;
+ xe=(N-1)-abs(Coeff_number);
+ qi=X0;
+ qe=N-1;
+ ki=X0;
+ ke=N-1;
+end
+xi;
+xe;
+arg_i = %pi*(ki*(1+rho)-qi);
+arg_e = %pi*(ke*(1+rho) - qe);
+sinc_i=sinc( arg_i);
+sinc_e=sinc( arg_e);
+div = xe-xi;
+pente=(sinc_e-sinc_i)/(xe-xi);
+//pente=(sinc_e-sinc_i);
+ordo=sinc_i;
+
+endfunction
diff --git a/cesar/cp/spoc/test/ApproxOrder2ParamMCoeff.sci b/cesar/cp/spoc/test/ApproxOrder2ParamMCoeff.sci
new file mode 100644
index 0000000000..fca72b620f
--- /dev/null
+++ b/cesar/cp/spoc/test/ApproxOrder2ParamMCoeff.sci
@@ -0,0 +1,34 @@
+function [args_i,args_e, sinc_i, sinc_e, div, ordo, pente] = ApproxOrder2ParamMCoeff(rho,X0,N,Coeff_number)
+
+
+// calcul des points pour approximer
+if Coeff_number<21
+ x1= X0 + 21-Coeff_number;
+ x3=N-1;
+elseif Coeff_number>21
+ x1= X0;
+ x3=N-1 -(Coeff_number-21);
+else
+ x1= X0;
+ x3=N-1;
+
+end
+
+args_i = %pi*((x1 + (Coeff_number-21))*(1+rho) - x1 );
+args_e = %pi*((x3 + (Coeff_number-21))*(1+rho) - x3 );
+sinc_i=sinc(args_i);
+sinc_e=sinc(args_e);
+
+
+// calcul des coefficients when y=ax²+c
+div = (x3 - x1)^2;
+pente= (sinc_e-sinc_i)/(div);
+ordo = sinc_i;
+
+endfunction
+
+
+
+
+
+
diff --git a/cesar/cp/spoc/test/Makefile b/cesar/cp/spoc/test/Makefile
new file mode 100755
index 0000000000..439b2e2572
--- /dev/null
+++ b/cesar/cp/spoc/test/Makefile
@@ -0,0 +1,11 @@
+BASE = ../../..
+
+INCLUDES = cp/spoc/test
+
+EXTRA_HOST_CFLAGS=-O3
+EXTRA_HOST_LDFLAGS=-lm
+HOST_PROGRAMS = test_coeff_check
+test_coeff_check_SOURCES = test_coeff_check.c
+test_coeff_check_MODULES = lib
+
+include $(BASE)/common/make/top.mk
diff --git a/cesar/cp/spoc/test/coeff_generation.scilab b/cesar/cp/spoc/test/coeff_generation.scilab
new file mode 100644
index 0000000000..f26946f208
--- /dev/null
+++ b/cesar/cp/spoc/test/coeff_generation.scilab
@@ -0,0 +1,22 @@
+pwd
+args = sciargs();
+rho = msscanf(args(5),"%lg")
+getf ("ApproxOrder1ParamMCoeff.sci");
+getf ("ApproxOrder2ParamMCoeff.sci");
+getf ("get_Approx1.sci");
+[args_tx_i,args_tx_e,sinc_tx_i, sinc_tx_e, div_tx, ordo_tx, pente_tx] = get_Approx1( rho,74,1155+74);
+[args_rx_i,args_rx_e,sinc_rx_i, sinc_rx_e, div_rx, ordo_rx, pente_rx] = get_Approx1(-rho,74,1155+74);
+fp = file("open","./obj/args_of_sinc.txt","unknown");
+fprintf (fp, "%5.15g %5.15g %5.15g %5.15g\n", args_tx_i, args_tx_e, args_rx_i, args_rx_e);
+file ("close", fp);
+fp = file("open","./obj/sinc.txt","unknown");
+fprintf (fp, "%5.15g %5.15g %5.15g %5.15g\n", sinc_tx_i, sinc_tx_e, sinc_rx_i, sinc_rx_e);
+file ("close", fp);
+fp = file("open","./obj/divisor.txt","unknown");
+fprintf (fp, "%5.15g %5.15g\n", div_tx, div_rx);
+file ("close", fp);
+fp = file("open","./obj/coeff.txt","unknown");
+fprintf (fp, "%5.15g %5.15g %5.15g %5.15g\n", ordo_tx, pente_tx, ordo_rx, pente_rx);
+file ("close", fp);
+exit;
+
diff --git a/cesar/cp/spoc/test/define_spoc_c.h b/cesar/cp/spoc/test/define_spoc_c.h
new file mode 100644
index 0000000000..325e6aa214
--- /dev/null
+++ b/cesar/cp/spoc/test/define_spoc_c.h
@@ -0,0 +1,20 @@
+#ifndef define_spoc_c_h
+#define define_spoc_c_h
+/* Cesar project {{{
+ *
+ * Copyright (C) 2008 Spidcom
+ *
+ * <<<Licence>>>
+ *
+ * }}} */
+/**
+ * \file define_spoc_c.h
+ * \brief « brief description »
+ * \ingroup « module »
+ *
+ * « long description »
+ */
+
+#define SPOC_COMPARE_SCILAB
+#include "cp/spoc/src/spoc.c"
+#endif /* define_spoc_c_h */
diff --git a/cesar/cp/spoc/test/get_Approx1.sci b/cesar/cp/spoc/test/get_Approx1.sci
new file mode 100644
index 0000000000..4ec86afc2f
--- /dev/null
+++ b/cesar/cp/spoc/test/get_Approx1.sci
@@ -0,0 +1,21 @@
+function [args_i, args_e, sinc_i, sinc_e, divs, ordos, pentes] = get_Approx1(rho,X0,N)
+for i=-20:20
+ if (i==0)
+ [arg_i, arg_e, si, se, div, ordo,pente] = ApproxOrder2ParamMCoeff(rho,X0,N,21);
+ else
+ [arg_i, arg_e, si, se, div, ordo,pente] = ApproxOrder1ParamMCoeff(rho,X0,N,i);
+ end
+ args_i(i+21)=arg_i;
+ args_e(i+21)=arg_e;
+ sinc_i(i+21)=si;
+ sinc_e(i+21)=se;
+ divs(i+21) = div;
+ ordos(i+21)=ordo;
+ pentes(i+21)=pente;
+end
+
+
+
+
+endfunction
+
diff --git a/cesar/cp/spoc/test/src/test_coeff_check.c b/cesar/cp/spoc/test/src/test_coeff_check.c
new file mode 100644
index 0000000000..5a50d41d0b
--- /dev/null
+++ b/cesar/cp/spoc/test/src/test_coeff_check.c
@@ -0,0 +1,194 @@
+/* Cesar project {{{
+ *
+ * Copyright (C) 2008 Spidcom
+ *
+ * <<<Licence>>>
+ *
+ * }}} */
+/**
+ * \file src/test_spoc.c
+ * \brief « brief description »
+ * \ingroup « module »
+ *
+ * « long description »
+ */
+#include "define_spoc_c.h"
+#include <stdio.h>
+#include <stdlib.h>
+#include "math.h"
+#include "lib/test.h"
+
+#define FAIL_LIMIT 1e-6
+#define SPOC_TEST_DBG 0
+
+double kaiser_d[10] = { 0.9905, 0.9622, 0.9164, 0.8546, 0.7792, 0.6929, 0.5989, 0.5005, 0.4011, 0.304 };
+
+double
+err (double a, double ref)
+{
+ double ret = 0;
+ //printf ("DBG : a=%5.15g ref=%5.15g\n", a, ref);
+ if (ref == 0)
+ {
+ ret = a;
+ }
+ else
+ {
+ if (a == 0) ret = ref;
+ else ret = (a-ref)/ref;
+ }
+ if (ret<0) ret = -ret;
+ return(ret);
+}
+
+double
+data_Mcheck_dbg_print (char* filestring, char* haut_str, char *bas_str,
+ spoc_prec_t *tx_coeff_haut, spoc_prec_t *tx_coeff_bas,
+ spoc_prec_t *rx_coeff_haut, spoc_prec_t *rx_coeff_bas )
+{
+ double tx_haut, tx_bas;
+ double rx_haut, rx_bas;
+ double tx_haut_err, tx_bas_err;
+ double rx_haut_err, rx_bas_err;
+ double ret=0;
+ FILE *fref = fopen (filestring, "r");
+ int i;
+ for (i=0; i<41; i++)
+ {
+ fscanf (fref, "%lg %lg %lg %lg", &tx_haut, &tx_bas, &rx_haut, &rx_bas);
+ tx_haut_err = err (tx_coeff_haut[i], tx_haut);
+ tx_bas_err = err(tx_coeff_bas[i], tx_bas);
+ rx_haut_err = err (rx_coeff_haut[i], rx_haut);
+ rx_bas_err = err(rx_coeff_bas[i], rx_bas);
+ if (ret < rx_haut_err) ret = rx_haut_err;
+ if (ret < rx_bas_err) ret = rx_bas_err;
+ if (ret < tx_haut_err) ret = tx_haut_err;
+ if (ret < tx_bas_err) ret = tx_bas_err;
+#if defined (SPOC_TEST_DBG) && (SPOC_TEST_DBG==2)
+ printf ("TX_%s[%d] %5.15g<-->%5.15g\t %5.15g\n", haut_str, i-20, tx_coeff_haut[i], tx_haut, tx_haut_err);
+ printf ("TX_%s[%d] %5.15g<-->%5.15g\t %5.15g ", bas_str, i-20, tx_coeff_bas[i], tx_bas, tx_bas_err);
+ if ((tx_bas_err > FAIL_LIMIT) || (tx_haut_err > FAIL_LIMIT)) printf (" * \n");
+ else printf ("\n");
+ printf ("RX_%s[%d] %5.15g<-->%5.15g\t %5.15g\n", haut_str, i-20, rx_coeff_haut[i], rx_haut, rx_haut_err);
+ printf ("RX_%s[%d] %5.15g<-->%5.15g\t %5.15g ", bas_str, i-20, rx_coeff_bas[i], rx_bas, rx_bas_err);
+ if ((rx_bas_err > FAIL_LIMIT) || (rx_haut_err > FAIL_LIMIT)) printf (" * \n");
+ else printf ("\n");
+#endif
+ }
+ fclose (fref);
+ return ret;
+}
+
+void
+test_Mcoeff_check (test_t test, double rho)
+{
+ char cmd_scilab[128];
+ char test_string[64];
+ sprintf (test_string, "rho=%lg\n", rho);
+ test_begin (test, test_string)
+ {
+ spoc_MCoeff (rho);
+ sprintf (cmd_scilab,"scilab -f ./coeff_generation.scilab -nogui -args %lg > ./obj/toto.log 2>./obj/titi.log </dev/null", rho);
+ system (cmd_scilab);
+#if defined (SPOC_TEST_DBG) && (SPOC_TEST_DBG==2)
+ printf ("\nargs\n");
+ data_Mcheck_dbg_print ("./obj/args_of_sinc.txt", "ARGS_HAUT", "ARGS_BAS ",
+ args_haut, args_bas, args_haut_rx, args_bas_rx);
+ printf ("\nsinc\n");
+ data_Mcheck_dbg_print ("./obj/sinc.txt", "SINC_HAUT", "SINC_BAS ",
+ sinc_haut, sinc_bas, sinc_haut_rx, sinc_bas_rx);
+ printf ("\ncoeff\n");
+#endif
+ double err = data_Mcheck_dbg_print ("./obj/coeff.txt", "ORDON", "PENTE",
+ reg_ordo_tx, reg_pente_tx, reg_ordo_rx, reg_pente_rx);
+#if defined (SPOC_TEST_DBG) && (SPOC_TEST_DBG>=1)
+ printf ("Matrix M ... err = %5.15g\n", err);
+#endif
+ test_fail_if (err > FAIL_LIMIT);
+ } test_end;
+}
+
+double
+sinc (double arg)
+{
+ if (arg == 0) return 1;
+ else return (sin(arg)/arg);
+}
+
+void
+test_CG_check (test_t test, double rho)
+{
+ char test_string[64];
+ sprintf (test_string, "rho=%lg\n", rho);
+ test_begin (test, test_string)
+ {
+ double CG_N[21], CG_FC10[21], CG_PR1[21], CG_PR2[21];
+ int k;
+ for (k=-10; k<=10; k++)
+ {
+ double kk, arg;
+ if (k==0) kk=1.0; else kk = kaiser_d[abs(k)-1];
+ arg = M_PI*(k-3072*rho); CG_N[k+10] = kk * sinc(arg);
+ arg = M_PI*(k+3072*rho); CG_FC10[k+10] = kk * sinc(arg);
+ arg = M_PI*(k+372*rho); CG_PR1[k+10] = kk * sinc(arg);
+ arg = M_PI*(k+1140*rho); CG_PR2[k+10] = kk * sinc(arg);
+ }
+ spoc_CG (rho, -3072, reg_CG_N);
+ spoc_CG (rho, 3072, reg_CG_FC10);
+ spoc_CG (rho, 372, reg_CG_PR1);
+ spoc_CG (rho, 1140, reg_CG_PR2);
+ double err_N, err_FC10, err_PR1, err_PR2, error_global;
+ error_global = 0;
+ for (k=-10; k<=10; k++)
+ {
+ double error = 0;
+ int i = 10+k;
+ err_N = err (CG_N[i], reg_CG_N[i]);
+ if (error < err_N) error = err_N;
+ err_FC10 = err (CG_FC10[i], reg_CG_FC10[i]);
+ if (error < err_FC10) error = err_FC10;
+ err_PR1 = err (CG_PR1[i], reg_CG_PR1[i]);
+ if (error < err_PR1) error = err_PR1;
+ err_PR2 = err (CG_PR2[i], reg_CG_PR2[i]);
+ if (error < err_PR2) error = err_PR2;
+ if (error_global< error) error_global=error;
+#if defined (SPOC_TEST_DBG) && (SPOC_TEST_DBG==2)
+ printf ("N[%d] %5.15g<-->%5.15g %5.15g\n", k, CG_N[i], reg_CG_N[i], err_N);
+ printf ("FC10[%d] %5.15g<-->%5.15g %5.15g\n", k, CG_FC10[i], reg_CG_FC10[i], err_FC10);
+ printf ("PR1[%d] %5.15g<-->%5.15g %5.15g\n", k, CG_PR1[i], reg_CG_PR1[i], err_PR1);
+ printf ("PR2[%d] %5.15g<-->%5.15g %5.15g", k, CG_PR2[i], reg_CG_PR2[i], err_PR2);
+ if (error > FAIL_LIMIT) printf (" *\n"); else printf ("\n");
+#endif
+ }
+#if defined (SPOC_TEST_DBG) && (SPOC_TEST_DBG>=1)
+ printf ("Vectors [N,FC10,PR1,PR2] ...err = %5.15g\n", error_global);
+#endif
+ test_fail_if (error_global > FAIL_LIMIT);
+ } test_end;
+}
+
+int
+main (int argc, char **argv)
+{
+ test_t test;
+ test_init (test, argc, argv);
+ int step, rho_eps, max_eps;
+ for (step = 1; step<=100000; step = step*10)
+ {
+ rho_eps = -10 * step;
+ max_eps = 10 * step;
+ for (; rho_eps<=max_eps; rho_eps+=step)
+ {
+ double rho = rho_eps * 1e-9;
+#if defined (SPOC_TEST_DBG) && (SPOC_TEST_DBG>=1)
+ printf ("\nSPOC rho=%g \n", rho);
+#endif
+ test_case_begin(test, "check Mcoeff");
+ test_Mcoeff_check(test, rho);
+ test_case_begin(test, "check CGcoeff");
+ test_CG_check(test, rho);
+ }
+ }
+ test_result (test);
+ return (test_nb_failed (test) == 0 ? 0 : 1);
+}