mirror
/
libdspl-2.0
kopia lustrzana https://github.com/Dsplib/libdspl-2.0
1
0
Forkuj 0
Kod Zgłoszenia Projekty Wydania Wiki Aktywność

added matrix_t

pull/2/head
Dsplib 2018-09-20 23:22:29 +03:00
rodzic 4ec07a1b4f
commit 997fe0d8e7
9 zmienionych plików z 351 dodań i 20 usunięć
Pokaż statystyki Ściągnij plik aktualizacji Ściągnij plik porównania Expand all files Collapse all files

58
dspl/dox/ru/math.dox 100644
Wyświetl plik

@ -0,0 +1,58 @@
/*! ****************************************************************************
\ingroup SPEC_MATH_COMMON_GROUP
\fn int sinc(double* x, int n, double a, double* y)
\brief Функция \f$ \textrm{sinc}(x,a) = \frac{\sin(ax)}{ax}\f$
Функция рассчитывает занчения значения функции для вещественного вектора `x`.
<BR>
\param[in] x Указатель на вектор переменной \f$ x \f$.<BR>
Размер вектора `[n x 1]`.<BR>
Память должна быть выделена.<BR><BR>
\param[in] n Размер входного вектора `x`. <BR><BR>
\param[in] a Параметр функции
\f$ \textrm{sinc}(x,a) = \frac{\sin(ax)}{ax}\f$
\param[out] y Указатель на вектор значений функции.<BR>
Размер вектора `[n x 1]`.<BR>
Память должна быть выделена.<BR><BR>
\return
`RES_OK` Расчет произведен успешно.<BR>
В противном случае
\ref ERROR_CODE_GROUP "код ошибки".<BR>
Пример использования функции `sinc`:
\include sinc_test.c
<BR><BR>
В каталоге `dat` будут созданы три файла:<BR>
<pre>
sinc_test_1.0.txt параметр a = 1.0
sinc_test_pi.txt параметр a = pi
sinc_test_2pi.txt параметр a = 2*pi
</pre>
Кроме того программа GNUPLOT произведет построение графика
по сохраненным в файлах данным:
График: `img/sinc_test.png`
\image html sinc_test.png
Скрипт GNUPLOT для построения графиков из текстовых файлов:
\include sinc_test.plt
\author
Бахурин Сергей
www.dsplib.org
***************************************************************************** */

4
dspl/src/math.c
Wyświetl plik

@ -41,7 +41,7 @@ double DSPL_API dmod (double x, double y)
/*******************************************************************************
sinc(x) = sin(pi*x)/(pi*x)
*******************************************************************************/
int DSPL_API sinc(double* x, int n, double* y)
int DSPL_API sinc(double* x, int n, double a, double* y)
{
int k;
@ -51,7 +51,7 @@ int DSPL_API sinc(double* x, int n, double* y)
return ERROR_SIZE;
for(k = 0; k < n; k++)
y[k] = (x[k]==0.0) ? 1.0 : sin(M_PI*x[k])/(M_PI*x[k]);
y[k] = (x[k]==0.0) ? 1.0 : sin(a*x[k])/(a*x[k]);
return RES_OK;

194
dspl/src/matrix.c 100644
Wyświetl plik

@ -0,0 +1,194 @@
/*
* Copyright (c) 2015-2018 Sergey Bakhurin
* Digital Signal Processing Library [http://dsplib.org]
*
* This file is part of libdspl-2.0.
*
* is free software: you can redistribute it and/or modify
* it under the terms of the GNU Lesser General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* DSPL is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU Lesser General Public License
* along with Foobar. If not, see <http://www.gnu.org/licenses/>.
*/
#include <stdlib.h>
#include <string.h>
#include "dspl.h"
void transpose(double* a, int n, int m, double* b);
void transpose_cmplx(complex_t* a, int n, int m, complex_t* b);
void transpose_hermite(complex_t* a, int n, int m, complex_t* b);
/*******************************************************************************
Real matrx transpose
*******************************************************************************/
int DSPL_API matrix_create(matrix_t* a, int n, int m, int type)
{
if(!a)
return ERROR_PTR;
if(n < 1 || m < 1)
return ERROR_MATRIX_SIZE;
if(a->dat)
{
a->dat = (type & DAT_MASK) ?
(void*) realloc(a->dat, n*m*sizeof(complex_t)):
(void*) realloc(a->dat, n*m*sizeof(double));
}
else
{
a->dat = (type & DAT_MASK) ?
(void*) malloc(n*m*sizeof(complex_t)):
(void*) malloc(n*m*sizeof(double));
}
a->n = n;
a->m = m;
a->type = type;
return RES_OK;
}
/*******************************************************************************
Real matrx transpose
*******************************************************************************/
void transpose(double* a, int n, int m, double* b)
{
int p, q, i, j, aind, bind;
for(p = 0; p < n - DSPL_MATRIX_BLOCK; p+=DSPL_MATRIX_BLOCK)
{
for(q = 0; q < m - DSPL_MATRIX_BLOCK; q+=DSPL_MATRIX_BLOCK)
{
for(i = 0; i < DSPL_MATRIX_BLOCK; i++)
{
for(j = 0; j < DSPL_MATRIX_BLOCK; j++)
{
aind = (q+j) * n + p + i;
bind = (p+i) * m + q + j;
b[bind] = a[aind];
}
}
}
}
for(i = p; i < n; i++)
for(j = 0; j < m; j++)
b[i*m + j] = a[j*n+i];
for(i = 0; i < p; i++)
for(j = q; j < m; j++)
b[i*m + j] = a[j*n+i];
}
/*******************************************************************************
Complex matrx transpose
*******************************************************************************/
void transpose_cmplx(complex_t* a, int n, int m, complex_t* b)
{
int p, q, i, j, aind, bind;
for(p = 0; p < n - DSPL_MATRIX_BLOCK; p+=DSPL_MATRIX_BLOCK)
{
for(q = 0; q < m - DSPL_MATRIX_BLOCK; q+=DSPL_MATRIX_BLOCK)
{
for(i = 0; i < DSPL_MATRIX_BLOCK; i++)
{
for(j = 0; j < DSPL_MATRIX_BLOCK; j++)
{
aind = (q+j) * n + p + i;
bind = (p+i) * m + q + j;
RE(b[bind]) = RE(a[aind]);
IM(b[bind]) = IM(a[aind]);
}
}
}
}
for(i = p; i < n; i++)
{
for(j = 0; j < m; j++)
{
RE(b[i*m + j]) = RE(a[j*n+i]);
IM(b[i*m + j]) = IM(a[j*n+i]);
}
}
for(i = 0; i < p; i++)
{
for(j = q; j < m; j++)
{
RE(b[i*m + j]) = RE(a[j*n+i]);
IM(b[i*m + j]) = IM(a[j*n+i]);
}
}
}
/*******************************************************************************
Hermite matrx transpose
*******************************************************************************/
void transpose_hermite(complex_t* a, int n, int m, complex_t* b)
{
int p, q, i, j, aind, bind;
for(p = 0; p < n - DSPL_MATRIX_BLOCK; p+=DSPL_MATRIX_BLOCK)
{
for(q = 0; q < m - DSPL_MATRIX_BLOCK; q+=DSPL_MATRIX_BLOCK)
{
for(i = 0; i < DSPL_MATRIX_BLOCK; i++)
{
for(j = 0; j < DSPL_MATRIX_BLOCK; j++)
{
aind = (q+j) * n + p + i;
bind = (p+i) * m + q + j;
RE(b[bind]) = RE(a[aind]);
IM(b[bind]) = -IM(a[aind]);
}
}
}
}
for(i = p; i < n; i++)
{
for(j = 0; j < m; j++)
{
RE(b[i*m + j]) = RE(a[j*n+i]);
IM(b[i*m + j]) = IM(a[j*n+i]);
}
}
for(i = 0; i < p; i++)
{
for(j = q; j < m; j++)
{
RE(b[i*m + j]) = RE(a[j*n+i]);
IM(b[i*m + j]) = -IM(a[j*n+i]);
}
}
}

2
include/dspl.c
Wyświetl plik

@ -97,6 +97,7 @@ p_logspace logspace ;
p_low2bp low2bp ;
p_low2high low2high ;
p_low2low low2low ;
p_matrix_create matrix_create ;
p_poly_z2a_cmplx poly_z2a_cmplx ;
p_polyval polyval ;
p_polyval_cmplx polyval_cmplx ;
@ -229,6 +230,7 @@ void* dspl_load()
LOAD_FUNC(low2bp);
LOAD_FUNC(low2high);
LOAD_FUNC(low2low);
LOAD_FUNC(matrix_create);
LOAD_FUNC(poly_z2a_cmplx);
LOAD_FUNC(polyval);
LOAD_FUNC(polyval_cmplx);

33
include/dspl.h
Wyświetl plik

@ -43,15 +43,23 @@ typedef double complex_t[2];
typedef struct
{
complex_t* w;
complex_t* t0;
complex_t* t1;
int n;
int p2;
complex_t* w;
complex_t* t0;
complex_t* t1;
int n;
int p2;
} fft_t;
typedef struct
{
void* dat;
int n;
int m;
int type;
} matrix_t;
#define RE(x) (x[0])
#define IM(x) (x[1])
@ -99,6 +107,7 @@ typedef struct
/* K 0x11xxxxxx*/
/* L 0x12xxxxxx*/
/* M 0x13xxxxxx*/
#define ERROR_MATRIX_SIZE 0x13011909
/* N 0x14xxxxxx*/
#define ERROR_NEGATIVE 0x14050701
/* O 0x15xxxxxx*/
@ -124,11 +133,11 @@ typedef struct
/* Y 0x25xxxxxx*/
/* Z 0x26xxxxxx*/
#define DAT_MASK 0x00000001
#define DAT_DOUBLE 0x00000000
#define DAT_COMPLEX 0x00000001
#define DAT_DOUBLE 0
#define DAT_COMPLEX 1
#define DSPL_MATRIX_BLOCK 32
#define DSPL_SYMMETRIC 0x00000000
@ -547,6 +556,11 @@ DECLARE_FUNC(int, low2low, double* b
COMMA double* beta
COMMA double* alpha);
//------------------------------------------------------------------------------
DECLARE_FUNC(int, matrix_create, matrix_t* a
COMMA int n
COMMA int m
COMMA int type);
//------------------------------------------------------------------------------
DECLARE_FUNC(int, poly_z2a_cmplx, complex_t*
COMMA int
COMMA int
@ -606,6 +620,7 @@ DECLARE_FUNC(int, sin_cmplx, complex_t*
//------------------------------------------------------------------------------
DECLARE_FUNC(int, sinc, double* x
COMMA int n
COMMA double a
COMMA double* y);
//------------------------------------------------------------------------------
DECLARE_FUNC(int, sqrt_cmplx, complex_t*

2
release/include/dspl.c
Wyświetl plik

@ -97,6 +97,7 @@ p_logspace logspace ;
p_low2bp low2bp ;
p_low2high low2high ;
p_low2low low2low ;
p_matrix_create matrix_create ;
p_poly_z2a_cmplx poly_z2a_cmplx ;
p_polyval polyval ;
p_polyval_cmplx polyval_cmplx ;
@ -229,6 +230,7 @@ void* dspl_load()
LOAD_FUNC(low2bp);
LOAD_FUNC(low2high);
LOAD_FUNC(low2low);
LOAD_FUNC(matrix_create);
LOAD_FUNC(poly_z2a_cmplx);
LOAD_FUNC(polyval);
LOAD_FUNC(polyval_cmplx);

33
release/include/dspl.h
Wyświetl plik

@ -43,15 +43,23 @@ typedef double complex_t[2];
typedef struct
{
complex_t* w;
complex_t* t0;
complex_t* t1;
int n;
int p2;
complex_t* w;
complex_t* t0;
complex_t* t1;
int n;
int p2;
} fft_t;
typedef struct
{
void* dat;
int n;
int m;
int type;
} matrix_t;
#define RE(x) (x[0])
#define IM(x) (x[1])
@ -99,6 +107,7 @@ typedef struct
/* K 0x11xxxxxx*/
/* L 0x12xxxxxx*/
/* M 0x13xxxxxx*/
#define ERROR_MATRIX_SIZE 0x13011909
/* N 0x14xxxxxx*/
#define ERROR_NEGATIVE 0x14050701
/* O 0x15xxxxxx*/
@ -124,11 +133,11 @@ typedef struct
/* Y 0x25xxxxxx*/
/* Z 0x26xxxxxx*/
#define DAT_MASK 0x00000001
#define DAT_DOUBLE 0x00000000
#define DAT_COMPLEX 0x00000001
#define DAT_DOUBLE 0
#define DAT_COMPLEX 1
#define DSPL_MATRIX_BLOCK 32
#define DSPL_SYMMETRIC 0x00000000
@ -547,6 +556,11 @@ DECLARE_FUNC(int, low2low, double* b
COMMA double* beta
COMMA double* alpha);
//------------------------------------------------------------------------------
DECLARE_FUNC(int, matrix_create, matrix_t* a
COMMA int n
COMMA int m
COMMA int type);
//------------------------------------------------------------------------------
DECLARE_FUNC(int, poly_z2a_cmplx, complex_t*
COMMA int
COMMA int
@ -606,6 +620,7 @@ DECLARE_FUNC(int, sin_cmplx, complex_t*
//------------------------------------------------------------------------------
DECLARE_FUNC(int, sinc, double* x
COMMA int n
COMMA double a
COMMA double* y);
//------------------------------------------------------------------------------
DECLARE_FUNC(int, sqrt_cmplx, complex_t*

13
test/bin/gnuplot/sinc_test.plt 100644
Wyświetl plik

@ -0,0 +1,13 @@
set grid
set xlabel "x"
set lmargin at screen 0.10
set terminal pngcairo size 560,280 enhanced font 'Verdana,8'
set output 'img/sinc_test.png'
set ylabel "sinc(x,a)"
set yrange [-0.25:1.1]
plot 'dat/sinc_test_1.0.txt' with lines title "a = 1.0", \
'dat/sinc_test_pi.txt' with lines title "a = pi", \
'dat/sinc_test_2pi.txt' with lines title "a = 2pi"

32
test/src/sinc_test.c 100644
Wyświetl plik

@ -0,0 +1,32 @@
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include "dspl.h"
#define N 1000
int main()
{
void* handle; // DSPL handle
handle = dspl_load(); // Load DSPL function
double x[N], y[N];
linspace(-10.0, 10.0, N , DSPL_SYMMETRIC, x);
sinc(x, N, 1.0, y);
writetxt(x, y, N, "dat/sinc_test_1.0.txt");
sinc(x, N, M_PI, y);
writetxt(x, y, N, "dat/sinc_test_pi.txt");
sinc(x, N, 2.0*M_PI, y);
writetxt(x, y, N, "dat/sinc_test_2pi.txt");
dspl_free(handle); // free dspl handle
// выполнить скрипт GNUPLOT для построения графиков
// по рассчитанным данным
return system("gnuplot gnuplot/sinc_test.plt");;
}