text:
// Copyright Twitch Interactive, Inc. or its affiliates. All Rights Reserved.
// SPDX-License-Identifier: MIT

#pragma once

#include <gsl/span>
#include <vector>

namespace Twitch::Utility {
class OutputByteStream {
public:
    void reserve(size_t size)
    {
        _data.reserve(size);
    }

    template<typename ValueType, typename ArgType>
    void write(const ArgType &arg)
    {
        size_t startSize = _data.size();
        size_t valueSize = sizeof(ValueType);
        ValueType value = static_cast<ValueType>(arg);
        _data.resize(startSize + valueSize);
        memcpy_s(_data.data() + startSize, _data.size() - startSize, &value, valueSize);
    }

    void writeBytes(gsl::span<const uint8_t> bytes)
    {
        size_t startSize = _data.size();
        _data.resize(startSize + bytes.size());
        memcpy_s(_data.data() + startSize, _data.size() - startSize, bytes.data(), bytes.size());
    }

    std::vector<uint8_t> consume()
    {
        return std::move(_data);
    }

private:
    std::vector<uint8_t> _data;
};

class InputByteStream {
public:
    InputByteStream(gsl::span<const uint8_t> data)
        : _data(std::move(data))
    {
    }

    template<typename ValueType>
    ValueType read()
    {
        ValueType value;
        size_t valueSize = sizeof(ValueType);
        memcpy_s(&value, sizeof(value), _data.data() + _currentPosition, valueSize);
        _currentPosition += valueSize;
        return value;
    }

    void readBytes(gsl::span<uint8_t> bytes)
    {
        memcpy_s(bytes.data(), bytes.size(), _data.data() + _currentPosition, bytes.size());
        _currentPosition += bytes.size();
    }

    gsl::span<const uint8_t> readBytesView(size_t size)
    {
        gsl::span<const uint8_t> bytesView = _data.subspan(_currentPosition, size);
        _currentPosition += size;
        return bytesView;
    }

    size_t remainingSize() const
    {
        return _data.size() - _currentPosition;
    }

private:
    gsl::span<const uint8_t> _data;
    size_t _currentPosition = 0;
};
} // namespace Twitch::Utility

max_issues_repo_name: twitchtv/libsoundtrackutil


text:
/*
 * Copyright 2016 Maikel Nadolski
 *
 * Licensed under the Apache License, Version 2.0 (the "License");
 * you may not use this file except in compliance with the License.
 * You may obtain a copy of the License at
 *
 *     http://www.apache.org/licenses/LICENSE-2.0

 * Unless required by applicable law or agreed to in writing, software
 * distributed under the License is distributed on an "AS IS" BASIS,
 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
 * See the License for the specific language governing permissions and
 * limitations under the License.
 */

#ifndef HMM_ALGORITHM_H_
#define HMM_ALGORITHM_H_

#include <map>
#include <range/v3/core.hpp>
#include <range/v3/algorithm.hpp>
#include <range/v3/view/map.hpp>
#include <range/v3/view/unique.hpp>
#include <gsl_util.h>

#include "maikel/hmm/algorithm/forward.h"
#include "maikel/hmm/algorithm/backward.h"
#include "maikel/hmm/algorithm/baum_welch.h"

namespace maikel {

  template <class SymbolType, class IndexType>
    // requires UnsignedIntegral<I>
    bool is_bijective_index_map(std::map<SymbolType,IndexType> const& map)
    {
      IndexType max_index = ranges::max( map | ranges::view::values );
      std::vector<std::size_t> histogram(max_index+1);
      for (IndexType index : map|ranges::view::values)
        ++histogram[index];
      return ranges::all_of(histogram, [](std::size_t count){ return count == 1; });
    }

  template <class T>
  using ValueType = typename std::remove_reference<T>::type::value_type;

  template <class Index, class Range, class T = ranges::range_value_t<Range>>
    std::map<T,Index>
    map_from_symbols(Range&& range)
    {
      std::map<T,Index> symbols_to_index;
      Index index { 0 };
      ranges::for_each(range | ranges::view::unique, [&symbols_to_index,&index](T const& s){
          symbols_to_index[s] = index++;
      });
      Ensures(is_bijective_index_map(symbols_to_index));
      return symbols_to_index;
    }

}

#endif /* HMM_ALGORITHM_H_ */

ext: h

max_stars_repo_name: maikel/hidden-markov-model

max_issues_repo_name: maikel/Hidden-Markov-Model

max_forks_repo_name: maikel/Hidden-Markov-Model


text:
///  @file rainbow-gen-userpk.c
///  @brief A command-line tool for generating the user public keys.

#include <rainbow_keypair.h>
#include <stdlib.h>
#include <stdio.h>
#include <api.h>
#include <utils.h>
#include <string.h>
#include <blas.h>

int main(int argc, char **argv) {

    printf("%s\n", CRYPTO_ALGNAME);

    printf("msk size: %lu\n", CRYPTO_MASTER_SECRET_KEY_BYTES);
    printf("mpk size: %lu\n", CRYPTO_MASTER_PUBLIC_KEY_BYTES);
    printf("hash size: %d\n", _HASH_LEN);
    printf("signature size: %d\n\n", CRYPTO_BYTES);

    if (4 != argc) {
        printf("Usage:\n\n\trainbow-gen-userpk mpk_file_name identity upk_file_name\n\n");
        return -1;
    }

    //malloc upk

    uint8_t *_mpk = malloc(CRYPTO_MASTER_PUBLIC_KEY_BYTES);
    FILE *fp;
    unsigned r;

    fp = fopen( argv[1] , "r");
    if( NULL == fp ) {
        printf("fail to open master key file.\n");
        return -1;
    }
    r = byte_fget(fp, _mpk, CRYPTO_MASTER_PUBLIC_KEY_BYTES);
    fclose(fp);
    if (CRYPTO_MASTER_PUBLIC_KEY_BYTES != r) {
        printf("fail to load key file.\n");
        return -1;
    }

    unsigned char * _identity = malloc((_ID+1)/2);
    r = 0;

    generate_identity_hash(_identity, (unsigned char *) argv[2], strlen(argv[2]));
    if( NULL == _identity ) {
        printf("fail to create identity hash.\n");
        return -1;
    }
    printf("identity hash: ");
    for(unsigned i = 0; i<_ID; i++){
        printf("%hhu ",gf16v_get_ele(_identity,i));
    }
    printf("\n\n");

    uint8_t *_upk = malloc(sizeof(upk_t));

    //calculate usk with mpk and ID

    int re = calculate_upk((upk_t *) _upk, (mpk_t *) _mpk, _identity);
    if (0 != re) {
        printf("%s generate user-public-key fails.\n", CRYPTO_ALGNAME);
        return -1;
    }

    //write upk to disk
    fp = fopen(argv[3], "w+");
    if (NULL == fp) {
        printf("fail to open user public key file.\n");
        return -1;
    }
    byte_fdump(fp, CRYPTO_ALGNAME " user-public-key", _upk, sizeof(upk_t)); //pk speichern und beschriften
    fclose(fp);

    //free upk and rest

    free(_mpk);
    free(_identity);
    free(_upk);

    return 0;
}

ext: c

max_stars_repo_name: mnmuser/ibs-rainbow

max_issues_repo_name: mnmuser/ibs-rainbow

max_forks_repo_name: mnmuser/ibs-rainbow


text:
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
//#include <cem.h>
#include <gsl/gsl_linalg.h>
#include <gsl/gsl_matrix.h>
#include <gsl/gsl_vector.h>

void least_square_solver(double *matA, int n_row, int n_col, double *y, double *sol_x) {

    int i, j, k, tmp;
    double sum;
    gsl_matrix * A = gsl_matrix_alloc (n_row,n_col);
    gsl_matrix * V = gsl_matrix_alloc (n_col,n_col);
    gsl_vector * S = gsl_vector_alloc (n_col);
    gsl_vector * work = gsl_vector_alloc (n_col);
//    gsl_matrix * U = gsl_matrix_alloc (n_row-1,n_col_ids);
    gsl_vector * b = gsl_vector_alloc (n_row);
    gsl_vector * x = gsl_vector_alloc (n_col);

/*    //printf A and b
    FILE *fpA, *fpU, *fpb;
    fpA = fopen("A.dat","w");
    fpU = fopen("U.dat","w");
    fpb = fopen("b.dat","w");*/

    // need to copy matA to A, y to b!
    for (i=0;i<n_row;i++) {
        for (j=0; j<n_col; j++) {
            gsl_matrix_set (A, i, j, *(matA+i*n_col+j));
        }
        gsl_vector_set (b, i, y[i]);
    }
/*    for (i=0;i<n_row;i++) {
        for (j=0; j<n_col; j++) {
            fprintf(fpA,"%f ", gsl_matrix_get (A, i, j));
        }
        fprintf(fpA,"\n");
        fprintf(fpb,"%f\n", gsl_vector_get(b,i));
    }*/
    gsl_linalg_SV_decomp(A, V, S, work);    // on output, A is replaced by U!
/*    for (i=0;i<n_row;i++) {
        for (j=0; j<n_col; j++) {
            fprintf(fpU,"%f ", *(matA+i*n_col+j));
        }
        fprintf(fpU,"\n");
    }*/
    gsl_linalg_SV_solve (A, V, S, b, x);    // So, A should be used instead of U, here.

    for (i=0;i<n_col;i++) {
        sol_x[i] = gsl_vector_get(x,i);
    }

    gsl_matrix_free(A);
    gsl_matrix_free(V);
    gsl_vector_free(S);
    gsl_vector_free(work);
    gsl_vector_free(b);
    gsl_vector_free(x);

}


text:
#pragma once

#include <gsl/span>

#include <memory>

#include <qbuffer.h>
#include <qthread.h>

#if QT_VERSION >= 0x060000
    #include <QtMultimedia/QAudioSink>
#else
    #include <QtMultimedia/QAudioOutput>
#endif
namespace contour
{
class Audio: public QObject
{
    Q_OBJECT

  public:
    Audio();
    ~Audio() override;
  signals:
    void play(int volume, int duration, std::vector<int> const& notes);
  private slots:
    void handleStateChanged(QAudio::State state);
    void handlePlayback(int volume, int duration, std::vector<int> const& notes);

  private:
    void fillBuffer(int volume, int duration, gsl::span<const int> notes);
    std::vector<std::int16_t> createMusicalNote(double volume, int duration, int note_) noexcept;

    QByteArray byteArray_;
    QBuffer audioBuffer_;
    QThread soundThread_;
#if QT_VERSION >= 0x060000
    std::unique_ptr<QAudioSink> audio;
#else
    std::unique_ptr<QAudioOutput> audio;
#endif
};
} // namespace contour

max_stars_repo_name: christianparpart/libterminal

max_issues_repo_name: christianparpart/libterminal

max_forks_repo_name: christianparpart/libterminal


text:
/* rng/rng.c
 *
 * Copyright (C) 1996, 1997, 1998, 1999, 2000, 2007 James Theiler, Brian Gough
 *
 * This program is free software; you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation; either version 3 of the License, or (at
 * your option) any later version.
 *
 * This program 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 General Public License
 * along with this program; if not, write to the Free Software
 * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
 */

#include <config.h>
#include <stdlib.h>
#include <stdio.h>
#include <string.h>
#include <gsl/gsl_errno.h>
#include <gsl/gsl_rng.h>

gsl_rng *
gsl_rng_alloc (const gsl_rng_type * T)
{

  gsl_rng *r = (gsl_rng *) malloc (sizeof (gsl_rng));

  if (r == 0)
    {
      GSL_ERROR_VAL ("failed to allocate space for rng struct",
                        GSL_ENOMEM, 0);
    };

  r->state = malloc (T->size);

  if (r->state == 0)
    {
      free (r);         /* exception in constructor, avoid memory leak */

      GSL_ERROR_VAL ("failed to allocate space for rng state",
                        GSL_ENOMEM, 0);
    };

  r->type = T;

  gsl_rng_set (r, gsl_rng_default_seed);        /* seed the generator */

  return r;
}

int
gsl_rng_memcpy (gsl_rng * dest, const gsl_rng * src)
{
  if (dest->type != src->type)
    {
      GSL_ERROR ("generators must be of the same type", GSL_EINVAL);
    }

  memcpy (dest->state, src->state, src->type->size);

  return GSL_SUCCESS;
}

gsl_rng *
gsl_rng_clone (const gsl_rng * q)
{
  gsl_rng *r = (gsl_rng *) malloc (sizeof (gsl_rng));

  if (r == 0)
    {
      GSL_ERROR_VAL ("failed to allocate space for rng struct",
                        GSL_ENOMEM, 0);
    };

  r->state = malloc (q->type->size);

  if (r->state == 0)
    {
      free (r);         /* exception in constructor, avoid memory leak */

      GSL_ERROR_VAL ("failed to allocate space for rng state",
                        GSL_ENOMEM, 0);
    };

  r->type = q->type;

  memcpy (r->state, q->state, q->type->size);

  return r;
}

void
gsl_rng_set (const gsl_rng * r, unsigned long int seed)
{
  (r->type->set) (r->state, seed);
}

unsigned long int
gsl_rng_max (const gsl_rng * r)
{
  return r->type->max;
}

unsigned long int
gsl_rng_min (const gsl_rng * r)
{
  return r->type->min;
}

const char *
gsl_rng_name (const gsl_rng * r)
{
  return r->type->name;
}

size_t
gsl_rng_size (const gsl_rng * r)
{
  return r->type->size;
}

void *
gsl_rng_state (const gsl_rng * r)
{
  return r->state;
}

void
gsl_rng_print_state (const gsl_rng * r)
{
  size_t i;
  unsigned char *p = (unsigned char *) (r->state);
  const size_t n = r->type->size;

  for (i = 0; i < n; i++)
    {
      /* FIXME: we're assuming that a char is 8 bits */
      printf ("%.2x", *(p + i));
    }

}

void
gsl_rng_free (gsl_rng * r)
{
  RETURN_IF_NULL (r);
  free (r->state);
  free (r);
}


text:
//
//  likelihood.h
//  EpiGenMCMC
//
//  Created by Lucy Li on 13/04/2016.
//  Copyright (c) 2016 Lucy Li, Imperial College London. All rights reserved.
//

#ifndef __EpiGenMCMC__likelihood__
#define __EpiGenMCMC__likelihood__

#include <vector>
#include <stdio.h>
#include <limits>
#include <gsl/gsl_cdf.h>
#include "parameter.h"

class Likelihood {
public:
    Likelihood();
    double binomial_lik(double, double, int, int, int, int, bool);
    double binomial_lik(double, double, std::vector<double>::iterator, int, int, int, int, int, bool);
    double coalescent_lik(std::vector<double>::iterator, std::vector<double>::iterator, std::vector<double>::iterator, std::vector<double>::iterator, std::vector<double>::iterator, int, int, int, bool);
};

#endif /* defined(__EpiGenMCMC__likelihood__) */

max_stars_repo_name: lucymli/EpiGenMCMC

max_issues_repo_name: lucymli/EpiGenMCMC

max_forks_repo_name: lucymli/EpiGenMCMC


text:
/* Copyright 2015 The TensorFlow Authors. All Rights Reserved.

Licensed under the Apache License, Version 2.0 (the "License");
you may not use this file except in compliance with the License.
You may obtain a copy of the License at

    http://www.apache.org/licenses/LICENSE-2.0

Unless required by applicable law or agreed to in writing, software
distributed under the License is distributed on an "AS IS" BASIS,
WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
See the License for the specific language governing permissions and
limitations under the License.
==============================================================================*/
// Portions Copyright (c) Microsoft Corporation

#pragma once

#include <functional>
#include <memory>
#include <string>
#include <unordered_map>
#include <vector>
#include <gsl/gsl>

#include "core/common/common.h"
#include "core/framework/callback.h"
#include "core/platform/env_time.h"
#include "core/platform/telemetry.h"
#include "core/session/onnxruntime_c_api.h"  // for ORTCHAR_T

#ifndef _WIN32
#include <sys/types.h>
#include <unistd.h>
#endif

namespace onnxruntime {

#ifdef _WIN32
using PIDType = unsigned long;
using FileOffsetType = int64_t;
#else
using PIDType = pid_t;
using FileOffsetType = off_t;
#endif

/// \brief An interface used by the onnxruntime implementation to
/// access operating system functionality like the filesystem etc.
///
/// Callers may wish to provide a custom Env object to get fine grain
/// control.
///
/// All Env implementations are safe for concurrent access from
/// multiple threads without any external synchronization.
class Env {
 public:
  virtual ~Env() = default;

  /// \brief Returns a default environment suitable for the current operating
  /// system.
  ///
  /// Sophisticated users may wish to provide their own Env
  /// implementation instead of relying on this default environment.
  ///
  /// The result of Default() belongs to this library and must never be deleted.
  static const Env& Default();

  virtual int GetNumCpuCores() const = 0;

  /// \brief Returns the number of micro-seconds since the Unix epoch.
  virtual uint64_t NowMicros() const { return env_time_->NowMicros(); }

  /// \brief Returns the number of seconds since the Unix epoch.
  virtual uint64_t NowSeconds() const { return env_time_->NowSeconds(); }

  /// Sleeps/delays the thread for the prescribed number of micro-seconds.
  /// On Windows, it's the min time to sleep, not the actual one.
  virtual void SleepForMicroseconds(int64_t micros) const = 0;

  /**
   * Gets the length of the specified file.
   */
  virtual common::Status GetFileLength(
      const ORTCHAR_T* file_path, size_t& length) const = 0;

  /**
   * Copies the content of the file into the provided buffer.
   * @param file_path The path to the file.
   * @param offset The file offset from which to start reading.
   * @param length The length in bytes to read.
   * @param buffer The buffer in which to write.
   */
  virtual common::Status ReadFileIntoBuffer(
      const ORTCHAR_T* file_path, FileOffsetType offset, size_t length,
      gsl::span<char> buffer) const = 0;

  using MappedMemoryPtr = std::unique_ptr<char[], OrtCallbackInvoker>;

  /**
   * Maps the content of the file into memory.
   * This is a copy-on-write mapping, so any changes are not written to the
   * actual file.
   * @param file_path The path to the file.
   * @param offset The file offset from which to start the mapping.
   * @param length The length in bytes of the mapping.
   * @param[out] mapped_memory A smart pointer to the mapped memory which
   *             unmaps the memory (unless release()'d) when destroyed.
   */
  virtual common::Status MapFileIntoMemory(
      const ORTCHAR_T* file_path, FileOffsetType offset, size_t length,
      MappedMemoryPtr& mapped_memory) const = 0;

#ifdef _WIN32
  //Mainly for use with protobuf library
  virtual common::Status FileOpenRd(const std::wstring& path, /*out*/ int& fd) const = 0;
  //Mainly for use with protobuf library
  virtual common::Status FileOpenWr(const std::wstring& path, /*out*/ int& fd) const = 0;
#endif
  //Mainly for use with protobuf library
  virtual common::Status FileOpenRd(const std::string& path, /*out*/ int& fd) const = 0;
  //Mainly for use with protobuf library
  virtual common::Status FileOpenWr(const std::string& path, /*out*/ int& fd) const = 0;
  //Mainly for use with protobuf library
  virtual common::Status FileClose(int fd) const = 0;
  //This functions is always successful. It can't fail.
  virtual PIDType GetSelfPid() const = 0;

  // \brief Load a dynamic library.
  //
  // Pass "library_filename" to a platform-specific mechanism for dynamically
  // loading a library.  The rules for determining the exact location of the
  // library are platform-specific and are not documented here.
  //
  // On success, returns a handle to the library in "*handle" and returns
  // OK from the function.
  // Otherwise returns nullptr in "*handle" and an error status from the
  // function.
  virtual common::Status LoadDynamicLibrary(const std::string& library_filename, void** handle) const = 0;

  virtual common::Status UnloadDynamicLibrary(void* handle) const = 0;

  // \brief Get a pointer to a symbol from a dynamic library.
  //
  // "handle" should be a pointer returned from a previous call to LoadDynamicLibrary.
  // On success, store a pointer to the located symbol in "*symbol" and return
  // OK from the function. Otherwise, returns nullptr in "*symbol" and an error
  // status from the function.
  virtual common::Status GetSymbolFromLibrary(void* handle, const std::string& symbol_name, void** symbol) const = 0;

  // \brief build the name of dynamic library.
  //
  // "name" should be name of the library.
  // "version" should be the version of the library or NULL
  // returns the name that LoadDynamicLibrary() can use
  virtual std::string FormatLibraryFileName(const std::string& name, const std::string& version) const = 0;

  // \brief returns a provider that will handle telemetry on the current platform
  virtual const Telemetry& GetTelemetryProvider() const = 0;

 protected:
  Env();

 private:
  ORT_DISALLOW_COPY_ASSIGNMENT_AND_MOVE(Env);
  EnvTime* env_time_ = EnvTime::Default();
};

}  // namespace onnxruntime

max_stars_repo_name: nrandell/onnxruntime

max_issues_repo_name: nrandell/onnxruntime

max_forks_repo_name: nrandell/onnxruntime


text:
/* -*- mode: C; c-basic-offset: 4 -*- */
/* ex: set shiftwidth=4 tabstop=4 expandtab: */
/*
 * Copyright (c) 2018-2019, Colorado School of Mines
 * All rights reserved.
 *
 * Author(s): Neil T. Dantam <ndantam@miens.edu>
 * Georgia Tech Humanoid Robotics Lab
 * Under Direction of Prof. Mike Stilman <mstilman@cc.gatech.edu>
 *
 *
 * This file is provided under the following "BSD-style" License:
 *
 *
 *   Redistribution and use in source and binary forms, with or
 *   without modification, are permitted provided that the following
 *   conditions are met:
 *
 *   * Redistributions of source code must retain the above copyright
 *     notice, this list of conditions and the following disclaimer.
 *
 *   * Redistributions in binary form must reproduce the above
 *     copyright notice, this list of conditions and the following
 *     disclaimer in the documentation and/or other materials provided
 *     with the distribution.
 *
 *   THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND
 *   CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES,
 *   INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF
 *   MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
 *   DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR
 *   CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
 *   SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
 *   LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF
 *   USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED
 *   AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
 *   LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
 *   ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
 *   POSSIBILITY OF SUCH DAMAGE.
 *
 */

#ifndef AMINO_MAT_H
#define AMINO_MAT_H

#include <cblas.h>

/**
 * @file mat.h
 *
 * Block matrix descriptors and linear algebra operations.
 *
 */

typedef size_t aa_la_size;

/**
 * Descriptor for a vector.
 */
struct aa_dvec {
    aa_la_size len;   ///< Number of elements in vector
    double *data; ///< Pointer to data
    aa_la_size inc;   ///< Increment between successive vector elements
};

/**
 * Descriptor for a block matrix.
 */
struct aa_dmat {
    aa_la_size rows;    ///< number of rows in matrix
    aa_la_size cols;    ///< number of columns
    double *data;   ///< Pointer to matrix data
    aa_la_size ld;      ///< Leading dimension of matrix
};

#define AA_DVEC_REF(v,i) ((v)->data[(i) * (v)->inc])

/**
 * Reference a matrix entry.
 *
 * @param M matrix
 * @param i row
 * @param j column
 */
#define AA_DMAT_REF(M,i,j) (((M)->data)[(M)->ld*(j) + (i)])

typedef void
(aa_la_err_fun)( const char *message );

AA_API void
aa_la_err( const char *message );

AA_API void
aa_la_set_err( aa_la_err_fun *fun );

/**
 * BLAS arguments for a vector
 */
#define AA_VEC_ARGS(X) (X->data), ((int)(X->inc))

/**
 * BLAS arguments for a matrix
 */
#define AA_MAT_ARGS(X) (X->data), ((int)(X->ld))

/* Construction */

/**
 * Fill in a vector descriptor.
 *
 * @param len  Number of elements in vector
 * @param data Pointer to vector data
 * @param inc  Increment between sucessive elements
 */
static inline struct aa_dvec
AA_DVEC_INIT( size_t len, double *data, size_t inc )
{
    struct aa_dvec vec;
    vec.len = len;
    vec.data = data;
    vec.inc = inc;
    return vec;
}

/**
 * Fill in a vector descriptor.
 *
 * @param vec  Pointer to descriptor
 * @param len  Number of elements in vector
 * @param data Pointer to vector data
 * @param inc  Increment between sucessive elements
 */
AA_API void
aa_dvec_view( struct aa_dvec *vec, size_t len, double *data, size_t inc );

/**
 * Fill in a matrix descriptor.
 *
 * @param mat   Pointer to descriptor
 * @param rows  Number of rows in matrix
 * @param cols  Number of colums in matrix
 * @param data  Pointer to vector data
 * @param ld    Leading dimension of matrix
 */
AA_API void
aa_dmat_view( struct aa_dmat *mat, size_t rows, size_t cols, double *data, size_t ld );

/**
 * View a block of a matrix.
 */
AA_API void
aa_dmat_view_block( struct aa_dmat *dst,
                    const struct aa_dmat *src,
                    size_t row_start, size_t col_start,
                    size_t rows, size_t cols );

/**
 * View a slice of a vector.
 */
AA_API void
aa_dvec_slice( const struct aa_dvec *src,
               size_t start,
               size_t stop,
               size_t step,
               struct aa_dvec *dst );

/**
 * View a block of a matrix.
 */
AA_API void
aa_dmat_block( const struct aa_dmat *src,
               size_t row_start, size_t col_start,
               size_t row_end, size_t col_end,
               struct aa_dmat *dst );

/**
 * View a row of a matrix as a vector.
 */
AA_API void
aa_dmat_row_vec( const struct aa_dmat *src, size_t row, struct aa_dvec *dst );

/**
 * View a column of a matrix as a vector.
 */
AA_API void
aa_dmat_col_vec( const struct aa_dmat *src, size_t col, struct aa_dvec *dst );

/**
 * View the diagonal of a matrix as a vector.
 */
AA_API void
aa_dmat_diag_vec( const struct aa_dmat *src, struct aa_dvec *dst );

/**
 * Fill in a matrix descriptor.
 *
 * @param rows  Number of rows in matrix
 * @param cols  Number of colums in matrix
 * @param data  Pointer to vector data
 * @param ld    Leading dimension of matrix
 */
static inline struct aa_dmat
AA_DMAT_INIT( size_t rows, size_t cols, double *data, size_t ld )
{
    struct aa_dmat mat;
    mat.rows = rows;
    mat.cols = cols;
    mat.data = data;
    mat.ld = ld;
    return mat;
}

#define AA_MAT_DIAG(VEC,MAT)                                    \
    aa_dvec_view((VEC), (MAT)->cols, (MAT)->data, 1+(MAT)->ld);

/**
 * Region-allocate a vector.
 *
 * When finished, pop the descriptor.
 */
AA_API struct aa_dvec *
aa_dvec_alloc( struct aa_mem_region *reg, size_t len );

/**
 * Duplicate vector out of region
 *
 * When finished, pop the descriptor.
 */
AA_API struct aa_dvec *
aa_dvec_dup( struct aa_mem_region *reg, const struct aa_dvec *src);

/**
 * Duplicate matrix out of region
 *
 * When finished, pop the descriptor.
 */
AA_API struct aa_dmat *
aa_dmat_dup( struct aa_mem_region *reg, const struct aa_dmat *src);

/**
 * Region-allocate a matrix.
 *
 * When finished, pop the descriptor.
 */
AA_API struct aa_dmat *
aa_dmat_alloc( struct aa_mem_region *reg, size_t rows, size_t cols );

/**
 * Heap-allocate a vector.
 *
 * The descriptor and data are contained in a single malloc()'ed block.
 * When finished, call free() on the descriptor.
 */
AA_API struct aa_dvec *
aa_dvec_malloc( size_t len );

/**
 * Heap-allocate a matrix.
 *
 * The descriptor and data are contained in a single malloc()'ed block.
 * When finished, call free() on the descriptor.
 */
AA_API struct aa_dmat *
aa_dmat_malloc( size_t rows, size_t cols );

/**
 * Zero a vector.
 */
AA_API void
aa_dvec_zero( struct aa_dvec *vec );

/**
 * Fill a vector.
 */
AA_API void
aa_dvec_set( struct aa_dvec *vec, double alpha );

/**
 * Fill a matrix diagonal and off-diagonal elements.
 *
 * @param[in] A        The matrix to fill
 * @param[in] alpha    off-diagonal entries of A
 * @param[in] beta     diagonal entries of A
 */
AA_API void
aa_dmat_set( struct aa_dmat *A, double alpha, double beta );

/**
 * Zero a matrix.
 */
AA_API void
aa_dmat_zero( struct aa_dmat *mat );

/* Level 1 BLAS */

/**
 * Swap x and y
 *
 * \f[ \mathbf{x} \leftrightarrow \mathbf{y} \f]
 */
AA_API void
aa_dvec_swap( struct aa_dvec *x, struct aa_dvec *y );

/**
 * Scale x by alpha.
 *
 * \f[ \mathbf{x} \leftarrow \alpha \mathbf{x} \f]
 */
AA_API void
aa_dvec_scal( double alpha, struct aa_dvec *x );

/**
 * Increment x by alpha.
 *
 * \f[ \mathbf{x} \leftarrow \alpha + \mathbf{x} \f]
 */
AA_API void
aa_dvec_inc( double alpha, struct aa_dvec *x );

/**
 * Copy x to y.
 *
 * \f[ \mathbf{y} \leftarrow \mathbf{x} \f]
 */
AA_API void
aa_dvec_copy( const struct aa_dvec *x, struct aa_dvec *y );

/**
 * Alpha x plus y.
 *
 * \f[ \mathbf{y} \leftarrow \alpha \mathbf{x} + \mathbf{y} \f]
 */
AA_API void
aa_dvec_axpy( double a, const struct aa_dvec *x, struct aa_dvec *y );

/**
 * Dot product
 *
 * \f[ \mathbf{x}^T \mathbf{y} \f]
 */
AA_API double
aa_dvec_dot( const struct aa_dvec *x, struct aa_dvec *y );

/**
 * Euclidean Norm
 *
 * \f[ \left\Vert \mathbf{x} \right\Vert_2  \f]
 */
AA_API double
aa_dvec_nrm2( const struct aa_dvec *x );

/* Level 2 BLAS */

/**
 * General Matrix-Vector multiply
 *
 * \f[ \mathbf{y} \leftarrow \alpha \mathbf{A}^{\rm op} \mathbf{x} + \beta \mathbf{y}  \f]
 */
AA_API void
aa_dmat_gemv( CBLAS_TRANSPOSE trans,
              double alpha, const struct aa_dmat *A,
              const struct aa_dvec *x,
              double beta, struct aa_dvec *y );

/* Level 3 BLAS */

/**
 * General Matrix-Matrix multiply
 *
 * \f[ \mathbf{y} \leftarrow \alpha \mathbf{A}^{\rm opA} \mathbf{B}^\rm{opB} + \beta \mathbf{C}  \f]
 */
AA_API void
aa_dmat_gemm( CBLAS_TRANSPOSE transA, CBLAS_TRANSPOSE transB,
              double alpha, const struct aa_dmat *A,
              const struct aa_dmat *B,
              double beta, struct aa_dmat *C );

/* LAPACK */

/**
 * Copies all or part of a two-dimensional matrix A to another
 * matrix B.
 *
 *  @param[in] UPLO
 *          Specifies the part of the matrix A to be copied to B.
 *          - = 'U':      Upper triangular part
 *          - = 'L':      Lower triangular part
 *          - Otherwise:  All of the matrix A
 *
 *  @param[in] A
 *          dimension (LDA,N)
 *          The m by n matrix A.  If UPLO = 'U', only the upper triangle
 *          or trapezoid is accessed; if UPLO = 'L', only the lower
 *          triangle or trapezoid is accessed.
 *
 *  @param[out] B
 *          dimension (LDB,N)
 *          On exit, B = A in the locations specified by UPLO.
 *
 */
AA_API void
aa_dmat_lacpy( const char uplo[1],
               const struct aa_dmat *A,
               struct aa_dmat *B );

/* Matrix/Vector Functions */

/**
 * sum-square-differences of two vectors
 */
AA_API double
aa_dvec_ssd( const struct aa_dvec *x, const struct aa_dvec *y);

/**
 * Y += alpha * X
 */
AA_API void
aa_dmat_axpy( double alpha, const struct aa_dmat *X, struct aa_dmat *Y);

/**
 * sum-square-differences of two matrices
 */
AA_API double
aa_dmat_ssd( const struct aa_dmat *x, const struct aa_dmat *y);

/**
 * Scale the matrix A by alpha
 */
AA_API void
aa_dmat_scal( struct aa_dmat *A, double alpha );

/**
 * Increment the matrix A by alpha
 */
AA_API void
aa_dmat_inc( struct aa_dmat *A, double alpha );

/**
 * Euclidean Norm
 */
AA_API double
aa_dmat_nrm2( const struct aa_dmat *x );

/**
 * Matrix transpose.
 */
AA_API void
aa_dmat_trans( const struct aa_dmat *A, struct aa_dmat *At);

/**
 * Matrix inverse, in-place.
 */
AA_API int
aa_dmat_inv( struct aa_dmat *A);

/**
 * Pseudo-inverse.
 *
 * Singular values less than tol are ignored.  If tol < 0, then a sane
 * default is used.
 */
AA_API int
aa_dmat_pinv( const struct aa_dmat *A, double tol, struct aa_dmat *As);

/**
 * Damped pseudo-inverse.
 */
AA_API int
aa_dmat_dpinv( const struct aa_dmat *A, double k, struct aa_dmat *As);

/**
 * Dead-zone damped pseudo-inverse.
 */
AA_API int
aa_dmat_dzdpinv(  const struct aa_dmat *A, double s_min, struct aa_dmat *As);

/**
 * Copy a matrix
 */
AA_API void
aa_dmat_copy(  const struct aa_dmat *A, struct aa_dmat *B);

#endif /* AMINO_MAT_H */

max_stars_repo_name: dyalab/amino

max_issues_repo_name: dyalab/amino

max_forks_repo_name: dyalab/amino


text:
/* vector/gsl_vector_uchar.h
 *
 * Copyright (C) 1996, 1997, 1998, 1999, 2000, 2007 Gerard Jungman, Brian Gough
 *
 * This program is free software; you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation; either version 3 of the License, or (at
 * your option) any later version.
 *
 * This program 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 General Public License
 * along with this program; if not, write to the Free Software
 * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
 */

#ifndef __GSL_VECTOR_UCHAR_H__
#define __GSL_VECTOR_UCHAR_H__

#if !defined( GSL_FUN )
#  if !defined( GSL_DLL )
#    define GSL_FUN extern
#  elif defined( BUILD_GSL_DLL )
#    define GSL_FUN extern __declspec(dllexport)
#  else
#    define GSL_FUN extern __declspec(dllimport)
#  endif
#endif

#include <stdlib.h>
#include <gsl/gsl_types.h>
#include <gsl/gsl_errno.h>
#include <gsl/gsl_inline.h>
#include <gsl/gsl_check_range.h>
#include <gsl/gsl_block_uchar.h>

#undef __BEGIN_DECLS
#undef __END_DECLS
#ifdef __cplusplus
# define __BEGIN_DECLS extern "C" {
# define __END_DECLS }
#else
# define __BEGIN_DECLS /* empty */
# define __END_DECLS /* empty */
#endif

__BEGIN_DECLS

typedef struct
{
  size_t size;
  size_t stride;
  unsigned char *data;
  gsl_block_uchar *block;
  int owner;
}
gsl_vector_uchar;

typedef struct
{
  gsl_vector_uchar vector;
} _gsl_vector_uchar_view;

typedef _gsl_vector_uchar_view gsl_vector_uchar_view;

typedef struct
{
  gsl_vector_uchar vector;
} _gsl_vector_uchar_const_view;

typedef const _gsl_vector_uchar_const_view gsl_vector_uchar_const_view;

/* Allocation */

GSL_FUN gsl_vector_uchar *gsl_vector_uchar_alloc (const size_t n);
GSL_FUN gsl_vector_uchar *gsl_vector_uchar_calloc (const size_t n);

GSL_FUN gsl_vector_uchar *gsl_vector_uchar_alloc_from_block (gsl_block_uchar * b,
                                                     const size_t offset,
                                                     const size_t n,
                                                     const size_t stride);

GSL_FUN gsl_vector_uchar *gsl_vector_uchar_alloc_from_vector (gsl_vector_uchar * v,
                                                      const size_t offset,
                                                      const size_t n,
                                                      const size_t stride);

GSL_FUN void gsl_vector_uchar_free (gsl_vector_uchar * v);

/* Views */

GSL_FUN _gsl_vector_uchar_view
gsl_vector_uchar_view_array (unsigned char *v, size_t n);

GSL_FUN _gsl_vector_uchar_view
gsl_vector_uchar_view_array_with_stride (unsigned char *base,
                                         size_t stride,
                                         size_t n);

GSL_FUN _gsl_vector_uchar_const_view
gsl_vector_uchar_const_view_array (const unsigned char *v, size_t n);

GSL_FUN _gsl_vector_uchar_const_view
gsl_vector_uchar_const_view_array_with_stride (const unsigned char *base,
                                               size_t stride,
                                               size_t n);

GSL_FUN _gsl_vector_uchar_view
gsl_vector_uchar_subvector (gsl_vector_uchar *v,
                            size_t i,
                            size_t n);

GSL_FUN _gsl_vector_uchar_view
gsl_vector_uchar_subvector_with_stride (gsl_vector_uchar *v,
                                        size_t i,
                                        size_t stride,
                                        size_t n);

GSL_FUN _gsl_vector_uchar_const_view
gsl_vector_uchar_const_subvector (const gsl_vector_uchar *v,
                                  size_t i,
                                  size_t n);

GSL_FUN _gsl_vector_uchar_const_view
gsl_vector_uchar_const_subvector_with_stride (const gsl_vector_uchar *v,
                                              size_t i,
                                              size_t stride,
                                              size_t n);

/* Operations */

GSL_FUN void gsl_vector_uchar_set_zero (gsl_vector_uchar * v);
GSL_FUN void gsl_vector_uchar_set_all (gsl_vector_uchar * v, unsigned char x);
GSL_FUN int gsl_vector_uchar_set_basis (gsl_vector_uchar * v, size_t i);

GSL_FUN int gsl_vector_uchar_fread (FILE * stream, gsl_vector_uchar * v);
GSL_FUN int gsl_vector_uchar_fwrite (FILE * stream, const gsl_vector_uchar * v);
GSL_FUN int gsl_vector_uchar_fscanf (FILE * stream, gsl_vector_uchar * v);
GSL_FUN int gsl_vector_uchar_fprintf (FILE * stream, const gsl_vector_uchar * v,
                              const char *format);

GSL_FUN int gsl_vector_uchar_memcpy (gsl_vector_uchar * dest, const gsl_vector_uchar * src);

GSL_FUN int gsl_vector_uchar_reverse (gsl_vector_uchar * v);

GSL_FUN int gsl_vector_uchar_swap (gsl_vector_uchar * v, gsl_vector_uchar * w);
GSL_FUN int gsl_vector_uchar_swap_elements (gsl_vector_uchar * v, const size_t i, const size_t j);

GSL_FUN unsigned char gsl_vector_uchar_max (const gsl_vector_uchar * v);
GSL_FUN unsigned char gsl_vector_uchar_min (const gsl_vector_uchar * v);
GSL_FUN void gsl_vector_uchar_minmax (const gsl_vector_uchar * v, unsigned char * min_out, unsigned char * max_out);

GSL_FUN size_t gsl_vector_uchar_max_index (const gsl_vector_uchar * v);
GSL_FUN size_t gsl_vector_uchar_min_index (const gsl_vector_uchar * v);
GSL_FUN void gsl_vector_uchar_minmax_index (const gsl_vector_uchar * v, size_t * imin, size_t * imax);

GSL_FUN int gsl_vector_uchar_add (gsl_vector_uchar * a, const gsl_vector_uchar * b);
GSL_FUN int gsl_vector_uchar_sub (gsl_vector_uchar * a, const gsl_vector_uchar * b);
GSL_FUN int gsl_vector_uchar_mul (gsl_vector_uchar * a, const gsl_vector_uchar * b);
GSL_FUN int gsl_vector_uchar_div (gsl_vector_uchar * a, const gsl_vector_uchar * b);
GSL_FUN int gsl_vector_uchar_scale (gsl_vector_uchar * a, const unsigned char x);
GSL_FUN int gsl_vector_uchar_add_constant (gsl_vector_uchar * a, const double x);
GSL_FUN int gsl_vector_uchar_axpby (const unsigned char alpha, const gsl_vector_uchar * x, const unsigned char beta, gsl_vector_uchar * y);
GSL_FUN unsigned char gsl_vector_uchar_sum (const gsl_vector_uchar * a);

GSL_FUN int gsl_vector_uchar_equal (const gsl_vector_uchar * u,
                            const gsl_vector_uchar * v);

GSL_FUN int gsl_vector_uchar_isnull (const gsl_vector_uchar * v);
GSL_FUN int gsl_vector_uchar_ispos (const gsl_vector_uchar * v);
GSL_FUN int gsl_vector_uchar_isneg (const gsl_vector_uchar * v);
GSL_FUN int gsl_vector_uchar_isnonneg (const gsl_vector_uchar * v);

GSL_FUN INLINE_DECL unsigned char gsl_vector_uchar_get (const gsl_vector_uchar * v, const size_t i);
GSL_FUN INLINE_DECL void gsl_vector_uchar_set (gsl_vector_uchar * v, const size_t i, unsigned char x);
GSL_FUN INLINE_DECL unsigned char * gsl_vector_uchar_ptr (gsl_vector_uchar * v, const size_t i);
GSL_FUN INLINE_DECL const unsigned char * gsl_vector_uchar_const_ptr (const gsl_vector_uchar * v, const size_t i);

#ifdef HAVE_INLINE

INLINE_FUN
unsigned char
gsl_vector_uchar_get (const gsl_vector_uchar * v, const size_t i)
{
#if GSL_RANGE_CHECK
  if (GSL_RANGE_COND(i >= v->size))
    {
      GSL_ERROR_VAL ("index out of range", GSL_EINVAL, 0);
    }
#endif
  return v->data[i * v->stride];
}

INLINE_FUN
void
gsl_vector_uchar_set (gsl_vector_uchar * v, const size_t i, unsigned char x)
{
#if GSL_RANGE_CHECK
  if (GSL_RANGE_COND(i >= v->size))
    {
      GSL_ERROR_VOID ("index out of range", GSL_EINVAL);
    }
#endif
  v->data[i * v->stride] = x;
}

INLINE_FUN
unsigned char *
gsl_vector_uchar_ptr (gsl_vector_uchar * v, const size_t i)
{
#if GSL_RANGE_CHECK
  if (GSL_RANGE_COND(i >= v->size))
    {
      GSL_ERROR_NULL ("index out of range", GSL_EINVAL);
    }
#endif
  return (unsigned char *) (v->data + i * v->stride);
}

INLINE_FUN
const unsigned char *
gsl_vector_uchar_const_ptr (const gsl_vector_uchar * v, const size_t i)
{
#if GSL_RANGE_CHECK
  if (GSL_RANGE_COND(i >= v->size))
    {
      GSL_ERROR_NULL ("index out of range", GSL_EINVAL);
    }
#endif
  return (const unsigned char *) (v->data + i * v->stride);
}
#endif /* HAVE_INLINE */

__END_DECLS

#endif /* __GSL_VECTOR_UCHAR_H__ */

max_stars_repo_name: zzpwahaha/Chimera-Control-Trim

max_issues_repo_name: zzpwahaha/Chimera-Control-Trim

max_forks_repo_name: zzpwahaha/Chimera-Control-Trim


text:
//
// Created by chris on 07.05.21.
//
#pragma once

#include <gsl/gsl_errno.h>
#include <gsl/gsl_fft_real.h>
#include <gsl/gsl_fft_halfcomplex.h>
#include <numeric>
#include <vector>

class test {
public:
    int n;
};

/**
 * Calculates autocorrelation with FFT
 * @param vec
 * @return
 *//*
std::vector<double> ac(const std::vector<double> &vec){
    std::vector<double> fftReturn=vec;

    gsl_fft_real_wavetable *real;
    gsl_fft_halfcomplex_wavetable *hc;
    gsl_fft_real_workspace *work;

    //FFT
    work = gsl_fft_real_workspace_alloc(vec.size());
    real = gsl_fft_real_wavetable_alloc(vec.size());

    gsl_fft_real_transform(&fftReturn[0], 1, vec.size(), real, work);
    gsl_fft_real_wavetable_free(real);

    hc = gsl_fft_halfcomplex_wavetable_alloc(vec.size());

    gsl_fft_halfcomplex_inverse(&fftReturn[0], 1, vec.size(), hc, work);

    gsl_fft_halfcomplex_wavetable_free(hc);
    return fftReturn;
}*/

/**
 * @tparam T
 * @param vec
 * @return maximum element of given vector
 */
template<typename T>
T max(const std::vector<T> &vec) {
    T maxVal = vec[0];
    for (int i = 1; i < vec.size(); ++i) {
        maxVal = std::max(maxVal, vec[i]);
    }
    return maxVal;
}

/**
 * @tparam T
 * @param vec
 * @return minimum element of given vector
 */
template<typename T>
T min(const std::vector<T> &vec) {
    T minVal = vec[0];
    for (int i = 1; i < vec.size(); ++i) {
        minVal = std::min(minVal, vec[i]);
    }
    return minVal;
}

/**
 * Calculates the arithmetic mean of given vector-elements
 * @tparam T
 * @param vec
 * @return the arithmetic mean
 */
template<typename T>
T mean(const std::vector<T> &vec) {
    return std::accumulate(vec.begin(), vec.end(), 0.0)/vec.size();
}

/**
 * Calculates the arithmetic mean of given vector-of-vector-of-elements
 * @tparam T
 * @param vec vector of vector of elements (tensor of rank 2)
 * @return the arithmetic mean of all values
 */
template<typename T>
T mean(const std::vector<std::vector<T>> &vec) {
    std::vector<T> firstMeans;
    firstMeans.reserve(vec.size());
    for(auto &i:vec){
        firstMeans.push_back(mean(i));
    }
    return mean(firstMeans);
}

/**
 * Calculates the standard deviation of given vector-elements
 * @tparam T
 * @param vec
 * @return standard deviation
 */
template<typename T>
T stdev(const std::vector<T> &vec) {
    T summation=0;
    const T vec_mean=mean(vec);
    for (size_t i = 0; i < vec.size(); ++i) {
        summation+=pow(vec[i]-vec_mean,2);
    }
    return sqrt(summation*1.0/(vec.size()-1));
}

/**
 * Calculates the autocorrelation-function of given data-sample
 * @tparam T
 * @param vec of data-sample with length n
 * @return autocorrelation-function in vector with length n/2
 */
template<typename T>
std::vector<T> autoCorr(const std::vector<T> &vec) {
    float mean = std::accumulate(vec.begin(), vec.end(), 0.0f) / vec.size();

    std::vector<T> autocorrelation(vec.size() / 2);
    for (size_t t = 0; t < autocorrelation.size(); t++) {
        float n = 0; // Numerator
        float d = 0; // Denominator
        for (size_t i = 0; i < vec.size() - t; i++) {
            float xim = vec[i] - mean;
            n += xim * (vec[(i + t) % vec.size()] - mean);
            d += xim * xim;
        }
        autocorrelation[t] = n / d;
    }
    return autocorrelation;
}

/**
 * normalizes given vector call by reference
 * @tparam T
 * @param vec vector which should be normalized
 */
template<typename T>
void normalize(std::vector<T> &vec) {
    const T vec_max = max(vec);
    for (T &i : vec) {
        i /= 1.0 * vec_max;
    }
}