1 #ifndef SIGMATRANSFORM_UTIL_H
2 #define SIGMATRANSFORM_UTIL_H
26 namespace SigmaTransform {
28 const double PI = 3.141592654;
36 std::chrono::steady_clock::time_point
_tic;
75 point<N>(
double const& x) {
for(
auto& _x :
_data ) _x = x; }
76 point<N>( std::array<double,N> dat ) :
_data( dat ) {}
79 void operator=( std::array<int,N> dat ) {
for(
int k=0;k<N;++k)
_data[k]=(
double)dat[k]; }
82 point<N> operator+( point<N>
const& r )
const { point<N> p(
_data);
auto R=r.begin();
for(
auto& x : p ) x += *R++;
return std::move( p ); }
83 point<N> operator-( point<N>
const& r )
const { point<N> p(
_data);
auto R=r.begin();
for(
auto& x : p ) x -= *R++;
return std::move( p ); }
84 point<N> operator*( point<N>
const& r )
const { point<N> p(
_data);
auto R=r.begin();
for(
auto& x : p ) x *= *R++;
return std::move( p ); }
85 point<N> operator/( point<N>
const& r )
const { point<N> p(
_data);
auto R=r.begin();
for(
auto& x : p ) x /= *R++;
return std::move( p ); }
86 point<N> operator/(
double const& r )
const { point<N> p(
_data);
for(
auto& x : p ) x /= r;
return std::move( p ); }
87 point<N> operator*(
double const& r )
const { point<N> p(
_data);
for(
auto& x : p ) x *= r;
return std::move( p ); }
88 point<N> operator+(
double const& r )
const { point<N> p(
_data);
for(
auto& x : p ) x += r;
return std::move( p ); }
89 point<N> operator-(
double const& r )
const { point<N> p(
_data);
for(
auto& x : p ) x -= r;
return std::move( p ); }
92 void operator+=( point<N>
const& r ) {
auto R=r.begin();
for(
auto& x :
_data ) x += *R++; }
93 void operator-=( point<N>
const& r ) {
auto R=r.begin();
for(
auto& x :
_data ) x -= *R++; }
94 void operator*=( point<N>
const& r ) {
auto R=r.begin();
for(
auto& x :
_data ) x *= *R++; }
95 void operator/=( point<N>
const& r ) {
auto R=r.begin();
for(
auto& x :
_data ) x /= *R++; }
98 double sum()
const {
double s = 0;
for(
auto const& x :
_data ) s+=x;
return s; }
99 double prod()
const {
double p = 1;
for(
auto const& x :
_data ) p*=x;
return p; }
100 point<N> sq()
const { point<N> p(
_data );
for(
auto& x : p ) x *= x;
return std::move( p ); }
101 point<N> abs()
const { point<N> p(
_data );
for(
auto& x : p ) x = (x>=0)?x:-x;
return std::move( p ); }
104 point<N> apply(
double(*fun)(
double const&) )
const { point<N> p;
for(
int k=0;k<N;++k) p[k] = fun(
_data[k]);
return std::move( p ); }
105 point<N> apply(
double(*fun)(
double) )
const { point<N> p;
for(
int k=0;k<N;++k) p[k] = fun(
_data[k]);
return std::move( p ); }
108 point<N> operator>(
double const& d )
const { point<N> p(
_data);
for(
auto& x : p ) x = x>d;
return std::move( p ); }
109 point<N> operator<(
double const& d )
const { point<N> p(
_data);
for(
auto& x : p ) x = x<d;
return std::move( p ); }
110 point<N> operator<=(
double const& d )
const { point<N> p(
_data);
for(
auto& x : p ) x = x<=d;
return std::move( p ); }
111 point<N> operator>=(
double const& d )
const { point<N> p(
_data);
for(
auto& x : p ) x = x>=d;
return std::move( p ); }
112 point<N> operator==(
double const& d )
const { point<N> p(
_data);
for(
auto& x : p ) x = x==d;
return std::move( p ); }
113 point<N> operator>( point<N>
const& r )
const { point<N> p(
_data);
auto R=r.begin();
for(
auto&x : p) x = x>*R++;
return std::move( p ); }
114 point<N> operator<( point<N>
const& r )
const { point<N> p(
_data);
auto R=r.begin();
for(
auto&x : p) x = x<*R++;
return std::move( p ); }
115 point<N> operator>=( point<N>
const& r )
const { point<N> p(
_data);
auto R=r.begin();
for(
auto&x : p) x = x>=*R++;
return std::move( p ); }
116 point<N> operator<=( point<N>
const& r )
const { point<N> p(
_data);
auto R=r.begin();
for(
auto&x : p) x = x<=*R++;
return std::move( p ); }
117 point<N> operator==( point<N>
const& r )
const { point<N> p(
_data);
auto R=r.begin();
for(
auto&x : p) x = x==*R++;
return std::move( p ); }
118 operator bool() {
bool b=
true;
for(
auto&x:
_data)
if(x!=1) b=
false;
return b; }
121 double& operator[](
size_t const& ind) {
return _data[ind]; }
122 double const& operator[](
size_t const& ind)
const {
return _data[ind]; }
123 double*
get() {
return _data; }
126 typedef typename std::array<double,N>::iterator
iterator;
147 friend std::ostream& operator<<( std::ostream& os , point<N>
const& p ) { os<<
"(";
for(
auto& x : p ) os << x <<
" ";
return os <<
"\b)"; }
148 friend std::istream&
operator>>( std::istream& is ,
point<N>& p ) {
for(
auto& x : p ) is >> x;
return is; }
153 using cmpx = std::complex<double>;
154 using cxVec = std::vector<cmpx>;
156 template<
size_t N>
using diffFunc = std::function<point<N>(
const point<N>&)>;
157 template<
size_t N>
using winFunc = std::function<cmpx(const point<N>&)>;
158 template<
size_t N>
using actFunc = std::function<point<N>(
const point<N>&,
const point<N>&)>;
180 cxVec
loadAscii2D( std::string
const& filename,
int &x ,
int &y );
190 std::vector<double>
linspace(
const double &L ,
const double &R ,
const int &N );
199 std::vector<double>
FourierAxis(
const double &fs ,
const unsigned &len );
208 void StartRecursiveLoop(
const std::vector<int>& max , std::function<
void(
const std::vector<int>&)> toDo );
220 void recursiveLoop( std::vector<int> &index,
const std::vector<int>& max ,
const int& curr_ind ,
221 std::function<
void(
const std::vector<int>&)> toDo );
230 std::vector<point<N>>
meshgridN( std::array<std::vector<double>,N>
const& dom ) {
232 std::vector<int> sizeVec(0);
234 for(
auto& d : dom ) {
236 sizeVec.push_back( d.size() );
239 std::vector<point<N>> out; out.reserve( p );
244 for(
int k = 0 ; k < index.size() ; ++k )
245 p[k] = dom[k][index[k]];
247 out.emplace_back( p );
250 return std::move( out );
260 std::vector<point<N>>
meshgridN( std::vector<double>
const& dom ) {
262 std::vector<int> sizeVec(0);
264 for(
int k = 0 ; k < N;++k) {
266 sizeVec.push_back( dom.size() );
269 std::vector<point<N>> out; out.reserve( p );
274 for(
int k = 0 ; k < index.size() ; ++k )
275 p[k] = dom[index[k]];
277 out.emplace_back( p );
280 return std::move( out );
294 std::stringstream ss;
298 for(
const auto &c : out )
299 ss << c.real() <<
"," << c.imag() << std::endl;
301 std::ofstream os(filename);
303 if( dim.prod() == out.size() )
304 os <<
"dim=" << dim <<
"\n";
306 os <<
"dim=(" <<out.size()<<
")\n";
318 void save2file_bin( std::string
const& filename ,
const cxVec &out );
392 #endif //SIGMATRANSFORM_H