diff --git a/SineTools.py b/SineTools.py index 6379d58bf88918e25cce4849f8099d490b556952..1526ca0ae155955c57be1b7d3435a84c85bef86e 100644 --- a/SineTools.py +++ b/SineTools.py @@ -1,195 +1,1187 @@ -<!DOCTYPE html> -<html class="devise-layout-html"> -<head prefix="og: http://ogp.me/ns#"> -<meta charset="utf-8"> -<link rel="preload" href="/assets/application_utilities-3676200ca543122eb8a1e1722a7139b82fbc787011ec0c4c17ac75145f60120f.css" as="style" type="text/css"> -<link rel="preload" href="/assets/application-4f233d907f30a050ca7e40fbd91742d444d28e50691c51b742714df8181bf4e7.css" as="style" type="text/css"> -<link rel="preload" href="/assets/highlight/themes/white-21f90a158663d6eabb1646d83d9e353d6904978fbb8391ab39ab4d1e4a1042f3.css" as="style" type="text/css"> +# -*- coding: utf-8 -*- +""" +Created on Fri Jun 14 20:36:01 2013 +SineTools.py +auxiliary functions related to sine-approximation -<meta content="IE=edge" http-equiv="X-UA-Compatible"> +@author: Th. Bruns, B.Seeger +""" -<meta content="object" property="og:type"> -<meta content="GitLab" property="og:site_name"> -<meta content="Sign in" property="og:title"> -<meta content="GitLab Community Edition" property="og:description"> -<meta content="https://gitlab1.ptb.de/assets/gitlab_logo-7ae504fe4f68fdebb3c2034e36621930cd36ea87924c11ff65dbcb8ed50dca58.png" property="og:image"> -<meta content="64" property="og:image:width"> -<meta content="64" property="og:image:height"> -<meta content="https://gitlab1.ptb.de/users/sign_in" property="og:url"> -<meta content="summary" property="twitter:card"> -<meta content="Sign in" property="twitter:title"> -<meta content="GitLab Community Edition" property="twitter:description"> -<meta content="https://gitlab1.ptb.de/assets/gitlab_logo-7ae504fe4f68fdebb3c2034e36621930cd36ea87924c11ff65dbcb8ed50dca58.png" property="twitter:image"> +import numpy as np +from scipy import linalg as la +from scipy.sparse import coo_matrix, hstack +from scipy.sparse.linalg import lsqr as sp_lsqr +import matplotlib.pyplot as mp +import warnings -<title>Sign in ยท GitLab</title> -<meta content="GitLab Community Edition" name="description"> +import multiprocessing +from functools import partial +from multiprocessing import shared_memory +from tqdm.contrib.concurrent import process_map +import os -<link rel="shortcut icon" type="image/png" href="/uploads/-/system/appearance/favicon/1/favicon.ico" id="favicon" data-original-href="/uploads/-/system/appearance/favicon/1/favicon.ico" /> -<style> -@keyframes blinking-dot{0%{opacity:1}25%{opacity:0.4}75%{opacity:0.4}100%{opacity:1}}@keyframes blinking-scroll-button{0%{opacity:0.2}50%{opacity:1}100%{opacity:0.2}}@keyframes gl-spinner-rotate{0%{transform:rotate(0)}100%{transform:rotate(360deg)}}body.ui-indigo .navbar-gitlab{background-color:#292961}body.ui-indigo .navbar-gitlab .navbar-collapse{color:#d1d1f0}body.ui-indigo .navbar-gitlab .container-fluid .navbar-toggler{border-left:1px solid #6868b9;color:#d1d1f0}body.ui-indigo .navbar-gitlab .navbar-sub-nav>li>a:hover,body.ui-indigo .navbar-gitlab .navbar-sub-nav>li>a:focus,body.ui-indigo .navbar-gitlab .navbar-sub-nav>li>button:hover,body.ui-indigo .navbar-gitlab .navbar-sub-nav>li>button:focus,body.ui-indigo .navbar-gitlab .navbar-nav>li>a:hover,body.ui-indigo .navbar-gitlab .navbar-nav>li>a:focus,body.ui-indigo .navbar-gitlab .navbar-nav>li>button:hover,body.ui-indigo .navbar-gitlab .navbar-nav>li>button:focus{background-color:rgba(209,209,240,0.2)}body.ui-indigo .navbar-gitlab .navbar-sub-nav>li.active>a,body.ui-indigo .navbar-gitlab .navbar-sub-nav>li.active>button,body.ui-indigo .navbar-gitlab .navbar-sub-nav>li.dropdown.show>a,body.ui-indigo .navbar-gitlab .navbar-sub-nav>li.dropdown.show>button,body.ui-indigo .navbar-gitlab .navbar-nav>li.active>a,body.ui-indigo .navbar-gitlab .navbar-nav>li.active>button,body.ui-indigo .navbar-gitlab .navbar-nav>li.dropdown.show>a,body.ui-indigo .navbar-gitlab .navbar-nav>li.dropdown.show>button{color:#292961;background-color:#fff}body.ui-indigo .navbar-gitlab .navbar-sub-nav>li.line-separator,body.ui-indigo .navbar-gitlab .navbar-nav>li.line-separator{border-left:1px solid rgba(209,209,240,0.2)}body.ui-indigo .navbar-gitlab .navbar-sub-nav{color:#d1d1f0}body.ui-indigo .navbar-gitlab .nav>li{color:#d1d1f0}body.ui-indigo .navbar-gitlab .nav>li>a .notification-dot{border:2px solid #292961}body.ui-indigo .navbar-gitlab .nav>li>a.header-help-dropdown-toggle .notification-dot{background-color:#d1d1f0}body.ui-indigo .navbar-gitlab .nav>li>a.header-user-dropdown-toggle .header-user-avatar{border-color:#d1d1f0}@media (min-width: 576px){body.ui-indigo .navbar-gitlab .nav>li>a:hover,body.ui-indigo .navbar-gitlab .nav>li>a:focus{background-color:rgba(209,209,240,0.2)}}body.ui-indigo .navbar-gitlab .nav>li>a:hover svg,body.ui-indigo .navbar-gitlab .nav>li>a:focus svg{fill:currentColor}body.ui-indigo .navbar-gitlab .nav>li>a:hover .notification-dot,body.ui-indigo .navbar-gitlab .nav>li>a:focus .notification-dot{will-change:border-color, background-color;border-color:#4a4a82}body.ui-indigo .navbar-gitlab .nav>li>a:hover.header-help-dropdown-toggle .notification-dot,body.ui-indigo .navbar-gitlab .nav>li>a:focus.header-help-dropdown-toggle .notification-dot{background-color:#fff}body.ui-indigo .navbar-gitlab .nav>li.active>a,body.ui-indigo .navbar-gitlab .nav>li.dropdown.show>a{color:#292961;background-color:#fff}body.ui-indigo .navbar-gitlab .nav>li.active>a:hover svg,body.ui-indigo .navbar-gitlab .nav>li.dropdown.show>a:hover svg{fill:#292961}body.ui-indigo .navbar-gitlab .nav>li.active>a .notification-dot,body.ui-indigo .navbar-gitlab .nav>li.dropdown.show>a .notification-dot{border-color:#fff}body.ui-indigo .navbar-gitlab .nav>li.active>a.header-help-dropdown-toggle .notification-dot,body.ui-indigo .navbar-gitlab .nav>li.dropdown.show>a.header-help-dropdown-toggle .notification-dot{background-color:#292961}body.ui-indigo .navbar-gitlab .nav>li .impersonated-user svg,body.ui-indigo .navbar-gitlab .nav>li .impersonated-user:hover svg{fill:#292961}body.ui-indigo .navbar .title>a:hover,body.ui-indigo .navbar .title>a:focus{background-color:rgba(209,209,240,0.2)}body.ui-indigo .search form{background-color:rgba(209,209,240,0.2)}body.ui-indigo .search form:hover{background-color:rgba(209,209,240,0.3)}body.ui-indigo .search .search-input::placeholder{color:rgba(209,209,240,0.8)}body.ui-indigo .search .search-input-wrap .search-icon,body.ui-indigo .search .search-input-wrap .clear-icon{fill:rgba(209,209,240,0.8)}body.ui-indigo .search.search-active form{background-color:#fff}body.ui-indigo .search.search-active .search-input-wrap .search-icon{fill:rgba(209,209,240,0.8)}body.ui-indigo .nav-sidebar li.active>a{color:#2f2a6b}body.ui-indigo .nav-sidebar .fly-out-top-item a,body.ui-indigo .nav-sidebar .fly-out-top-item a:hover,body.ui-indigo .nav-sidebar .fly-out-top-item.active a,body.ui-indigo .nav-sidebar .fly-out-top-item .fly-out-top-item-container{background-color:#2f2a6b;color:var(--black, #fff)}body.ui-indigo .nav-links li.active a,body.ui-indigo .nav-links li.md-header-tab.active button,body.ui-indigo .nav-links li a.active{border-bottom:2px solid #6666c4}body.ui-indigo .nav-links li.active a .badge.badge-pill,body.ui-indigo .nav-links li.md-header-tab.active button .badge.badge-pill,body.ui-indigo .nav-links li a.active .badge.badge-pill{font-weight:600}body.ui-indigo .emoji-picker-category-active{border-bottom-color:#6666c4}body.ui-indigo .branch-header-title{color:#4b4ba3}body.ui-indigo .ide-sidebar-link.active{color:#4b4ba3}body.ui-indigo .ide-sidebar-link.active.is-right{box-shadow:inset -3px 0 #4b4ba3} +import os +from numba import njit, prange +canUseMP=False +if os.name == 'nt': + print("SineTools Multiprocessing is not impleneted on Windows ") +if os.name=='posix': + canUseMP=True + manager = multiprocessing.Manager() + mpdata = manager.dict() -*,*::before,*::after{box-sizing:border-box}html{font-family:sans-serif;line-height:1.15}header{display:block}body{margin:0;font-family:-apple-system, BlinkMacSystemFont, "Segoe UI", Roboto, "Noto Sans", Ubuntu, Cantarell, "Helvetica Neue", sans-serif, "Apple Color Emoji", "Segoe UI Emoji", "Segoe UI Symbol", "Noto Color Emoji";font-size:1rem;font-weight:400;line-height:1.5;color:#303030;text-align:left;background-color:#fff}hr{box-sizing:content-box;height:0;overflow:visible}h1,h3{margin-top:0;margin-bottom:0.25rem}p{margin-top:0;margin-bottom:1rem}a{color:#007bff;text-decoration:none;background-color:transparent}a:not([href]):not([class]){color:inherit;text-decoration:none}img{vertical-align:middle;border-style:none}svg{overflow:hidden;vertical-align:middle}label{display:inline-block;margin-bottom:0.5rem}input{margin:0;font-family:inherit;font-size:inherit;line-height:inherit}input{overflow:visible}[type="submit"]:not(:disabled){cursor:pointer}[type="submit"]::-moz-focus-inner{padding:0;border-style:none}fieldset{min-width:0;padding:0;margin:0;border:0}[hidden]{display:none !important}h1,h3{margin-bottom:0.25rem;font-weight:600;line-height:1.2;color:#303030}h1{font-size:2.1875rem}h3{font-size:1.53125rem}hr{margin-top:0.5rem;margin-bottom:0.5rem;border:0;border-top:1px solid rgba(0,0,0,0.1)}.container{width:100%;padding-right:15px;padding-left:15px;margin-right:auto;margin-left:auto}@media (min-width: 576px){.container{max-width:540px}}@media (min-width: 768px){.container{max-width:720px}}@media (min-width: 992px){.container{max-width:960px}}@media (min-width: 1200px){.container{max-width:1140px}}.row{display:flex;flex-wrap:wrap;margin-right:-15px;margin-left:-15px}.col,.col-sm-5,.col-sm-7,.col-sm-12{position:relative;width:100%;padding-right:15px;padding-left:15px}.col{flex-basis:0;flex-grow:1;max-width:100%}.order-1{order:1}.order-12{order:12}@media (min-width: 576px){.col-sm-5{flex:0 0 41.66667%;max-width:41.66667%}.col-sm-7{flex:0 0 58.33333%;max-width:58.33333%}.col-sm-12{flex:0 0 100%;max-width:100%}.order-sm-1{order:1}.order-sm-12{order:12}}.form-control{display:block;width:100%;height:34px;padding:0.375rem 0.75rem;font-size:0.875rem;font-weight:400;line-height:1.5;color:#303030;background-color:#fff;background-clip:padding-box;border:1px solid #dbdbdb;border-radius:0.25rem}.form-control:-moz-focusring{color:transparent;text-shadow:0 0 0 #303030}.form-control::placeholder{color:#5e5e5e;opacity:1}.form-control:disabled{background-color:#fafafa;opacity:1}.form-group{margin-bottom:1rem}.form-row{display:flex;flex-wrap:wrap;margin-right:-5px;margin-left:-5px}.form-row>.col{padding-right:5px;padding-left:5px}.btn{display:inline-block;font-weight:400;color:#303030;text-align:center;vertical-align:middle;-webkit-user-select:none;user-select:none;background-color:transparent;border:1px solid transparent;padding:0.375rem 0.75rem;font-size:1rem;line-height:20px;border-radius:0.25rem}.btn:disabled{opacity:0.65}.btn:not(:disabled):not(.disabled){cursor:pointer}fieldset:disabled a.btn{pointer-events:none}.navbar{position:relative;display:flex;flex-wrap:wrap;align-items:center;justify-content:space-between;padding:0.25rem 0.5rem}.navbar .container{display:flex;flex-wrap:wrap;align-items:center;justify-content:space-between}.d-block{display:block !important}.d-flex{display:flex !important}.flex-wrap{flex-wrap:wrap !important}.justify-content-between{justify-content:space-between !important}.align-items-center{align-items:center !important}.fixed-top{position:fixed;top:0;right:0;left:0;z-index:1030}.ml-2{margin-left:0.5rem !important}.mt-3{margin-top:1rem !important}.mb-3{margin-bottom:1rem !important}.text-center{text-align:center !important}.font-weight-normal{font-weight:400 !important}.gl-form-input,.gl-form-input.form-control{background-color:#fff;font-family:-apple-system, BlinkMacSystemFont, "Segoe UI", Roboto, "Noto Sans", Ubuntu, Cantarell, "Helvetica Neue", sans-serif, "Apple Color Emoji", "Segoe UI Emoji", "Segoe UI Symbol", "Noto Color Emoji";font-size:0.875rem;line-height:1rem;padding-top:0.5rem;padding-bottom:0.5rem;padding-left:0.75rem;padding-right:0.75rem;height:auto;color:#303030;box-shadow:inset 0 0 0 1px #868686;border-style:none;-webkit-appearance:none;appearance:none;-moz-appearance:none}.gl-form-input:disabled,.gl-form-input:not(.form-control-plaintext):read-only,.gl-form-input.form-control:disabled,.gl-form-input.form-control:not(.form-control-plaintext):read-only{background-color:#fafafa;color:#868686;box-shadow:inset 0 0 0 1px #dbdbdb;cursor:not-allowed}.gl-form-input::placeholder,.gl-form-input.form-control::placeholder{color:#868686}.gl-button{display:inline-flex}.gl-button:not(.btn-link):active{text-decoration:none}.gl-button.gl-button{border-width:0;padding-top:0.5rem;padding-bottom:0.5rem;padding-left:0.75rem;padding-right:0.75rem;background-color:transparent;line-height:1rem;color:#303030;fill:currentColor;box-shadow:inset 0 0 0 1px #bfbfbf;justify-content:center;align-items:center;font-size:0.875rem;border-radius:0.25rem}.gl-button.gl-button .gl-button-icon{height:1rem;width:1rem;flex-shrink:0;margin-right:0.25rem;top:auto}.gl-button.gl-button.btn-default{background-color:#fff}.gl-button.gl-button.btn-default:active{box-shadow:inset 0 0 0 2px #5e5e5e,0 0 0 1px rgba(255,255,255,0.4),0 0 0 4px rgba(31,117,203,0.48);outline:none;background-color:#dbdbdb}.gl-button.gl-button.btn-confirm{color:#fff}.gl-button.gl-button.btn-confirm{background-color:#1f75cb;box-shadow:inset 0 0 0 1px #1068bf}.gl-button.gl-button.btn-confirm:active{box-shadow:inset 0 0 0 2px #033464,0 0 0 1px rgba(255,255,255,0.4),0 0 0 4px rgba(31,117,203,0.48);outline:none;background-color:#0b5cad}body,.form-control{font-size:0.875rem}[type="submit"]{cursor:pointer}h1,h3{margin-top:20px;margin-bottom:10px}a{color:#1068bf}hr{overflow:hidden}.hidden{display:none !important;visibility:hidden !important}.hide{display:none}svg{vertical-align:baseline}html{overflow-y:scroll}body{text-decoration-skip:ink}body.navless{background-color:#fff !important}.container{padding-top:0;z-index:5}.container .content{margin:0}@media (max-width: 575.98px){.container .content{margin-top:20px}}.navless-container{margin-top:40px;padding-top:32px}.btn{border-radius:4px;font-size:0.875rem;font-weight:400;padding:6px 10px;background-color:#fff;border-color:#dbdbdb;color:#303030;color:#303030;white-space:nowrap}.btn:active{background-color:#f0f0f0;box-shadow:none}.btn:active{background-color:#eaeaea;border-color:#e3e3e3;color:#303030}.btn svg{height:15px;width:15px}.btn svg:not(:last-child){margin-right:5px}.light{color:#303030}hr{margin:1.5rem 0;border-top:1px solid #eee}.footer-links{margin-bottom:20px}.footer-links a{margin-right:15px}.flash-container{margin:0;margin-bottom:16px;font-size:14px;position:relative;z-index:1}.flash-container.sticky{position:sticky;top:48px;z-index:251}.flash-container.flash-container-page{margin-bottom:0}.flash-container:empty{margin:0}input{border-radius:0.25rem;color:#303030;background-color:#fff}label{font-weight:600}label.label-bold{font-weight:600}.form-control{border-radius:4px;padding:6px 10px}.form-control::placeholder{color:#868686}.gl-show-field-errors .form-control:not(textarea){height:34px}.gl-show-field-errors .gl-field-hint{color:#303030}.navbar-empty{justify-content:center;height:40px;background:#fff;border-bottom:1px solid #dbdbdb}.navbar-empty .tanuki-logo,.navbar-empty .brand-header-logo{max-height:100%}.tanuki-logo .tanuki-left-ear,.tanuki-logo .tanuki-right-ear,.tanuki-logo .tanuki-nose{fill:#e24329}.tanuki-logo .tanuki-left-eye,.tanuki-logo .tanuki-right-eye{fill:#fc6d26}.tanuki-logo .tanuki-left-cheek,.tanuki-logo .tanuki-right-cheek{fill:#fca326}input::-moz-placeholder{color:#868686;opacity:1}input::-ms-input-placeholder{color:#868686}input:-ms-input-placeholder{color:#868686}svg{fill:currentColor}.login-page .container{max-width:960px}.login-page .navbar-gitlab .container{max-width:none}.login-page .flash-container{margin-bottom:16px}.login-page .brand-holder{font-size:18px;line-height:1.5}.login-page .brand-holder p{font-size:16px;color:#888}.login-page .brand-holder h3{font-size:22px}.login-page .brand-holder img{max-width:100%;margin-bottom:30px}.login-page .brand-holder a{font-weight:600}.login-page p{font-size:13px}.login-page .login-box,.login-page .omniauth-container{box-shadow:0 0 0 1px #dbdbdb;border-radius:0.25rem;padding:15px}.login-page .login-box .login-heading h3,.login-page .omniauth-container .login-heading h3{font-weight:400;line-height:1.5;margin:0 0 10px}.login-page .login-box .login-footer,.login-page .omniauth-container .login-footer{margin-top:10px}.login-page .login-box .login-footer p:last-child,.login-page .omniauth-container .login-footer p:last-child{margin-bottom:0}.login-page .login-box a.forgot,.login-page .omniauth-container a.forgot{float:right;padding-top:6px}.login-page .login-box .nav .active a,.login-page .omniauth-container .nav .active a{background:transparent}.login-page .login-box .login-body,.login-page .omniauth-container .login-body{font-size:13px}.login-page .login-box .login-body input+p,.login-page .login-box .login-body input ~ p.field-validation,.login-page .omniauth-container .login-body input+p,.login-page .omniauth-container .login-body input ~ p.field-validation{margin-top:5px}.login-page .login-box .login-body .username .validation-success,.login-page .omniauth-container .login-body .username .validation-success{color:#217645}.login-page .login-box .login-body .username .validation-error,.login-page .omniauth-container .login-body .username .validation-error{color:#dd2b0e}.login-page .omniauth-container{border-radius:0.25rem;font-size:13px}.login-page .omniauth-container p{margin:0}.login-page .omniauth-container form{width:48%;padding:0;border:0;background:none;margin-bottom:16px}@media (max-width: 991.98px){.login-page .omniauth-container form{width:100%}}.login-page .omniauth-container .omniauth-btn{width:100%}.login-page .new-session-tabs{display:flex;box-shadow:0 0 0 1px #dbdbdb;border-top-right-radius:4px;border-top-left-radius:4px}.login-page .new-session-tabs.custom-provider-tabs{flex-wrap:wrap}.login-page .new-session-tabs.custom-provider-tabs li{min-width:85px;flex-basis:auto}.login-page .new-session-tabs.custom-provider-tabs li:nth-child(n+5){border-top:1px solid #dbdbdb}.login-page .new-session-tabs.custom-provider-tabs a{font-size:16px}.login-page .new-session-tabs li{flex:1;text-align:center;border-left:1px solid #dbdbdb}.login-page .new-session-tabs li:first-of-type{border-left:0;border-top-left-radius:4px}.login-page .new-session-tabs li:last-of-type{border-top-right-radius:4px}.login-page .new-session-tabs li:not(.active){background-color:#fafafa}.login-page .new-session-tabs li a{width:100%;font-size:18px}.login-page .new-session-tabs li.active>a{cursor:default}.login-page .form-control:active,.login-page .form-control:focus{background-color:#fff}.login-page .submit-container{margin-top:16px}.login-page input[type="submit"]{margin-bottom:0;display:block;width:100%}.login-page .devise-errors h2{margin-top:0;font-size:14px;color:#ae1800}@media (max-width: 575.98px){.login-page .col-md-5.float-right{float:none !important;margin-bottom:45px}}.devise-layout-html{margin:0;padding:0;height:100%}.devise-layout-html body{height:calc(100% - 51px);margin:0;padding:0}.devise-layout-html body.navless{height:calc(100% - 11px)}.devise-layout-html body .page-wrap{min-height:100%;position:relative}.devise-layout-html body .footer-container,.devise-layout-html body hr.footer-fixed{position:absolute;bottom:0;left:0;right:0;height:40px;background:#fff}.devise-layout-html body .login-page-broadcast{margin-top:40px}.devise-layout-html body .navless-container{padding:65px 15px}@media (max-width: 575.98px){.devise-layout-html body .navless-container{padding:0 15px 65px}}.gl-border-solid{border-style:solid}.gl-border-gray-100{border-color:#dbdbdb}.gl-border-1{border-width:1px}.gl-rounded-base{border-radius:0.25rem}.gl-text-green-600{color:#217645}.gl-text-red-500{color:#dd2b0e}.gl-display-flex{display:flex}.gl-display-block{display:block}.gl-align-items-center{align-items:center}.gl-w-full{width:100%}.gl-p-2{padding:0.25rem}.gl-p-4{padding:0.75rem}.gl-mt-2{margin-top:0.25rem}.gl-mb-2{margin-bottom:0.25rem}.gl-mb-3{margin-bottom:0.5rem}.gl-mb-5{margin-bottom:1rem}@media (min-width: 36rem){.gl-sm-mt-0{margin-top:0}}.gl-text-left{text-align:left}.content-wrapper>.alert-wrapper,#content-body,.modal-dialog{display:none} - -</style> - -<link rel="stylesheet" media="print" href="/assets/application-4f233d907f30a050ca7e40fbd91742d444d28e50691c51b742714df8181bf4e7.css" /> - -<link rel="stylesheet" media="print" href="/assets/application_utilities-3676200ca543122eb8a1e1722a7139b82fbc787011ec0c4c17ac75145f60120f.css" /> - - -<link rel="stylesheet" media="print" href="/assets/highlight/themes/white-21f90a158663d6eabb1646d83d9e353d6904978fbb8391ab39ab4d1e4a1042f3.css" /> -<script> -//<![CDATA[ -document.querySelectorAll('link[media="print"]').forEach(linkTag => { - linkTag.setAttribute('data-startupcss', 'loading'); - const startupLinkLoadedEvent = new CustomEvent('CSSStartupLinkLoaded'); - linkTag.addEventListener('load',function(){this.media='all';this.setAttribute('data-startupcss', 'loaded');document.dispatchEvent(startupLinkLoadedEvent);},{once: true}); -}) - -//]]> -</script> - -<script> -//<![CDATA[ -window.gon={};gon.features={"webauthn":false}; -//]]> -</script> - - - - - -<script src="/assets/webpack/runtime.f3c0b0d6.bundle.js" defer="defer"></script> -<script src="/assets/webpack/main.5565962b.chunk.js" defer="defer"></script> -<script src="/assets/webpack/commons-pages.admin.sessions-pages.ldap.omniauth_callbacks-pages.omniauth_callbacks-pages.profiles.t-f04c18ab.f60a0bce.chunk.js" defer="defer"></script> -<script src="/assets/webpack/commons-pages.admin.sessions-pages.sessions-pages.sessions.new.d49a906e.chunk.js" defer="defer"></script> -<script src="/assets/webpack/pages.sessions.new.099e6dce.chunk.js" defer="defer"></script> - -<meta name="csrf-param" content="authenticity_token" /> -<meta name="csrf-token" content="XkP21DHJ+s6JVfQJRVi97dJqV3VYF4ACRuPAprbMCJ8+O9yUxU6WeoRbZrUKHT5MbIbtwpoyHoeccY19F+Ydng==" /> - -<meta name="action-cable-url" content="/-/cable" /> -<meta content="width=device-width, initial-scale=1, maximum-scale=1" name="viewport"> -<meta content="#292961" name="theme-color"> -<link rel="apple-touch-icon" type="image/x-icon" href="/assets/touch-icon-iphone-5a9cee0e8a51212e70b90c87c12f382c428870c0ff67d1eb034d884b78d2dae7.png" /> -<link rel="apple-touch-icon" type="image/x-icon" href="/assets/touch-icon-ipad-a6eec6aeb9da138e507593b464fdac213047e49d3093fc30e90d9a995df83ba3.png" sizes="76x76" /> -<link rel="apple-touch-icon" type="image/x-icon" href="/assets/touch-icon-iphone-retina-72e2aadf86513a56e050e7f0f2355deaa19cc17ed97bbe5147847f2748e5a3e3.png" sizes="120x120" /> -<link rel="apple-touch-icon" type="image/x-icon" href="/assets/touch-icon-ipad-retina-8ebe416f5313483d9c1bc772b5bbe03ecad52a54eba443e5215a22caed2a16a2.png" sizes="152x152" /> -<link color="rgb(226, 67, 41)" href="/assets/logo-d36b5212042cebc89b96df4bf6ac24e43db316143e89926c0db839ff694d2de4.svg" rel="mask-icon"> -<link href="/search/opensearch.xml" rel="search" title="Search GitLab" type="application/opensearchdescription+xml"> -<meta content="/assets/msapplication-tile-1196ec67452f618d39cdd85e2e3a542f76574c071051ae7effbfde01710eb17d.png" name="msapplication-TileImage"> -<meta content="#30353E" name="msapplication-TileColor"> - - - - -</head> - -<body class="ui-indigo login-page application navless gl-browser-generic gl-platform-other" data-page="sessions:new" data-qa-selector="login_page"> - -<script> -//<![CDATA[ -gl = window.gl || {}; -gl.client = {"isGeneric":true,"isOther":true}; - - -//]]> -</script> -<div class="page-wrap"> -<header class="navbar fixed-top navbar-empty"> -<img class="brand-header-logo lazy" data-src="/uploads/-/system/appearance/header_logo/1/ptb_logo.png" src="data:image/gif;base64,R0lGODlhAQABAAAAACH5BAEKAAEALAAAAAABAAEAAAICTAEAOw==" /> -</header> - -<div class="login-page-broadcast"> - - -</div> -<div class="container navless-container"> -<div class="content"> -<div class="flash-container flash-container-page sticky" data-qa-selector="flash_container"> -<div class="flash-alert mb-2"> -<svg class="s16 align-middle mr-1" data-testid="error-icon"><use xlink:href="/assets/icons-05c4d4d8f3cc1fe0f22064d47d6a57d254ff9686a08abb74993ade21581e46f8.svg#error"></use></svg> -<span>You need to sign in or sign up before continuing.</span> -<div class="close-icon-wrapper js-close-icon"> -<svg class="s16 close-icon gl-vertical-align-baseline!" data-testid="close-icon"><use xlink:href="/assets/icons-05c4d4d8f3cc1fe0f22064d47d6a57d254ff9686a08abb74993ade21581e46f8.svg#close"></use></svg> -</div> -</div> -</div> - -<div class="row mt-3"> -<div class="col-sm-12"> -<h1 class="mb-3 font-weight-normal"> -GitLab -</h1> -</div> -</div> -<div class="row mb-3"> -<div class="col-sm-7 order-12 order-sm-1 brand-holder"> - -<h3 class="gl-sm-mt-0"> -A complete DevOps platform -</h3> -<p> -GitLab is a single application for the entire software development lifecycle. From project planning and source code management to CI/CD, monitoring, and security. -</p> -<p> -This is a self-managed instance of GitLab. -</p> -<p data-sourcepos="1:1-1:66" dir="auto"><strong>New users may self-register with their PTB email address only!</strong></p>
<p data-sourcepos="3:1-4:70" dir="auto"><strong>For external users please send a request to gitlab.admin(at)ptb.de</strong>
Please provide full name, email and the name of the institute/company.</p> - -</div> -<div class="col-sm-5 order-1 new-session-forms-container order-sm-12"> - -<div id="signin-container"> -<div class="tab-content"> -<div class="login-box tab-pane active" id="login-pane" role="tabpanel"> -<div class="login-body"> -<form class="new_user gl-show-field-errors" id="new_user" aria-live="assertive" action="/users/sign_in" accept-charset="UTF-8" method="post"><input type="hidden" name="authenticity_token" value="lzyqLM7u0x5ttq0wUUoAbZB3uy1+I7YrDFaKYeY2XST3RIBsOmm/qmC4P4weD4PMLpsBmrwGKK7WxMe6RxxIJQ==" /><div class="form-group"> -<label for="user_login" class="label-bold">Username or email</label> -<input class="form-control gl-form-input top" autofocus="autofocus" autocapitalize="off" autocorrect="off" required="required" title="This field is required." data-qa-selector="login_field" type="text" name="user[login]" id="user_login" /> -</div> -<div class="form-group"> -<label class="label-bold" for="user_password">Password</label> -<input class="form-control gl-form-input bottom" required="required" title="This field is required." data-qa-selector="password_field" type="password" name="user[password]" id="user_password" /> -</div> -<div class="remember-me"> -<label for="user_remember_me"> -<input name="user[remember_me]" type="hidden" value="0" /><input class="remember-me-checkbox" type="checkbox" value="1" name="user[remember_me]" id="user_remember_me" /> -<span>Remember me</span> -</label> -<div class="float-right"> -<a href="/users/password/new">Forgot your password?</a> -</div> -</div> -<div></div> -<div class="submit-container move-submit-down"> -<input type="submit" name="commit" value="Sign in" class="gl-button btn btn-confirm" data-qa-selector="sign_in_button" data-disable-with="Sign in" /> -</div> -</form> -</div> -</div> - -</div> -<p class="gl-mt-3"> -Don't have an account yet? -<a data-qa-selector="register_link" href="/users/sign_up">Register now</a> -</p> -</div> - -</div> -</div> -</div> -</div> -<hr class="footer-fixed"> -<div class="container footer-container"> -<div class="footer-links"> -<a href="/explore">Explore</a> -<a href="/help">Help</a> -<a href="https://about.gitlab.com/">About GitLab</a> -</div> -</div> - - -</div> -</body> -</html> + +def sampletimes(Fs, T): # + """ + generate a t_i vector with \n + sample rate Fs \n + from 0 to T + """ + num = int(np.ceil(T * Fs)) + return np.linspace(0, T, num, dtype=np.float64) + + +# a = displacement, velocity or acceleration amplitude +def sinewave(f, a, phi, ti, offset=0, noise=0, absnoise=0, drift=0, ampdrift=0): + """ + generate a sampled sine wave s_i = s(t_i) with \n + amplitude a \n + initial phase phi \n + sample times t_i \n + bias offset (default 0)\n + noise as multiple of the amplitude in noise level \n + absnoise as a additive noise component \n + drift as multiples of amplitude per duration in drifting zero \n + ampdrift as a drifting amplitude given as multiple of amplitude \n + """ + Tau = ti[-1] - ti[0] + n = 0 + n = 0 + if noise != 0: + n = a * noise * np.random.randn(len(ti)) + if absnoise != 0: + n = n + absnoise * np.random.randn(len(ti)) + + d = drift * a / Tau + + s = ( + a * (1 + ampdrift / Tau * ti) * np.sin(2 * np.pi * f * ti - phi) + + n + + d * ti + + offset + ) + return s + + +def fm_counter_sine( + fm, f, x, phi, ti, offset=0, noise=0, absnoise=0, drift=0, ampdrift=0, lamb=633.0e-9 +): + """ + calculate counter value of heterodyne signal at \n + carrier freq. fm\n + x = displacement amplitude \n + initial phase phi \n + sample times t_i \n + bias or offset (default 0)\n + noise as multiple of the amplitude in noise level \n + absnoise as a additive noise component \n + drift as multiples of amplitude per duration in drifting zero \n + ampdrift as a drifting amplitude given as multiple of amplitude \n + lamb as wavelength of Laser + """ + Tau = ti[-1] - ti[0] + n = 0 + if noise != 0: + n = x * noise * np.random.randn(len(ti)) + if absnoise != 0: + n = n + absnoise * np.random.randn(len(ti)) + + d = drift * x / Tau + + s = ( + 1.0 + / lamb + * ( + x * (1 + ampdrift / Tau * ti) * np.sin(2 * np.pi * f * ti - phi) + + n + + d * ti + + offset + ) + ) + s = np.floor(s + fm * ti) + + return s + + +# sine fit at known frequency +@njit +def threeparsinefit(y, t, f0): + """ + sine-fit at a known frequency\n + y vector of sample values \n + t vector of sample times\n + f0 known frequency\n + \n + returns a vector of coefficients [a,b,c]\n + for y = a*sin(2*pi*f0*t) + b*cos(2*pi*f0*t) + c + """ + w0 = 2 * np.pi * f0 + lengtht=int(t.size) + a=np.ones((3,lengtht), dtype=np.float64) + a[0,:]=np.cos(w0 * t) + a[1,:]=np.sin(w0 * t) + #a = np.array([np.cos(w0 * t), np.sin(w0 * t), np.ones(t.size)]) + aT=a.transpose() + abc = np.linalg.lstsq(aT, y) + return abc[0][0:3] ## fit vector a*sin+b*cos+c + +def threeparsinefitNOJIT(y, t, f0): + """ + sine-fit at a known frequency\n + y vector of sample values \n + t vector of sample times\n + f0 known frequency\n + \n + returns a vector of coefficients [a,b,c]\n + for y = a*sin(2*pi*f0*t) + b*cos(2*pi*f0*t) + c + """ + w0 = 2 * np.pi * f0 + lengtht=int(t.size) + a=np.ones((3,lengtht), dtype=np.float64) + a[0,:]=np.cos(w0 * t) + a[1,:]=np.sin(w0 * t) + #a = np.array([np.cos(w0 * t), np.sin(w0 * t), np.ones(t.size)]) + aT=a.transpose() + abc = np.linalg.lstsq(aT, y) + return abc[0][0:3] ## fit vector a*sin+b*cos+c + + +# multiprocessing wraper for threeparsinefit +def multiPrcoess3paramFit(i): + """ + mpdata['t_shared_name'] = shmY.name + mpdata['y_shared_name'] = shmT.name + mpdata['t_shared_type'] = y.dtype + mpdata['y_shared_type'] = t.dtype + mpdata['freq'] = f0 + mpdata['samples_per_block'] = N + mpdata['nmumber_of_samples'] = y.size + """ + N=mpdata['samples_per_block'] + + existing_shm_t = shared_memory.SharedMemory(name=mpdata['t_shared_name']) + t = np.ndarray((mpdata['nmumber_of_samples'],), dtype=mpdata['t_shared_type'], buffer=existing_shm_t.buf) + ti = t[i * N: (i + 1) * N] + + existing_shm_y = shared_memory.SharedMemory(name=mpdata['y_shared_name']) + y = np.ndarray((mpdata['nmumber_of_samples'],), dtype=mpdata['y_shared_type'], buffer=existing_shm_y.buf) + yi = y[i * N: (i + 1) * N] + result=threeparsinefit(yi, ti, mpdata['freq']) + # Clean up + del t # Unnecessary; merely emphasizing the array is no longer used + del y # Unnecessary; merely emphasizing the array is no longer used + existing_shm_t.close() + existing_shm_y.close() + return result + +# sine fit at known frequency and detrending +def threeparsinefit_lin(y, t, f0): + """ + sine-fit with detrending at a known frequency\n + y vector of sample values \n + t vector of sample times\n + f0 known frequency\n + \n + returns a vector of coefficients [a,b,c,d]\n + for y = a*sin(2*pi*f0*t) + b*cos(2*pi*f0*t) + c*t + d + """ + w0 = 2 * np.pi * f0 + + a = np.array([np.cos(w0 * t), np.sin(w0 * t), np.ones(t.size), t, np.ones(t.size)]) + + abc = la.lstsq(a.transpose(), y) + return abc[0][0:4] ## fit vector + + +def calc_threeparsine(abc, t, f0): + """ + return y = abc[0]*sin(2*pi*f0*t) + abc[1]*cos(2*pi*f0*t) + abc[2] + """ + w0 = 2 * np.pi * f0 + return abc[0] * np.cos(w0 * t) + abc[1] * np.sin(w0 * t) + abc[2] + + +def amplitude(abc): + """ + return the amplitude given the coefficients of\n + y = a*sin(2*pi*f0*t) + b*cos(2*pi*f0*t) + c + """ + return np.absolute(abc[1] + 1j * abc[0]) + + +def phase(abc, deg=False): + """ + return the (sine-)phase given the coefficients of\n + y = a*sin(2*pi*f0*t) + b*cos(2*pi*f0*t) + c \n + returns angle in rad by default, in degree if deg=True + """ + return np.angle(abc[1] + 1j * abc[0], deg=deg) + + +def magnitude(A1, A2): + """ + return the magnitude of the complex ratio of sines A2/A1\n + given two sets of coefficients \n + A1 = [a1,b1,c1]\n + A2 = [a2,b2,c2] + """ + return amplitude(A2) / amplitude(A1) + + +def phase_delay(A1, A2, deg=False): + """ + return the phase difference of the complex ratio of sines A2/A1\n + given two sets of coefficients \n + A1 = [a1,b1,c1]\n + A2 = [a2,b2,c2]\n + returns angle in rad by default, in degree if deg=True + """ + return phase(A2, deg=deg) - phase(A1, deg=deg) + + +# periodical sinefit at known frequency +def seq_threeparsinefit(y, t, f0, periods=1,multiTasiking=False,returnSectionStartTimes=False): + """ + period-wise sine-fit at a known frequency\n + y vector of sample values \n + t vector of sample times\n + f0 known frequency\n + periods number of full periods for every fit (default=1)\n + \n + returns a (n,3)-matrix of coefficient-triplets [[a,b,c], ...]\n + for y = a*sin(2*pi*f0*t) + b*cos(2*pi*f0*t) + c + """ + if y.size < t.size: + raise ValueError("Dimension mismatch in input data y<t") + if t.size < y.size: + warnings.warn( + "Dimension mismatch in input data y>t. fiting only for t.size values", + RuntimeWarning, + ) + Tau = 1.0 / f0 *periods + dt = np.mean(np.diff(t)) + N = int(Tau // dt) ## samples per section + M = int(t.size // N) ## number of sections or periods + + abc = np.zeros((M, 3)) + + if multiTasiking==False: + for i in range(int(M)): + ti = t[i * N : (i + 1) * N] + yi = y[i * N : (i + 1) * N] + abc[i, :] = threeparsinefitNOJIT(yi, ti, f0) + if returnSectionStartTimes: + sectionStartTimes = t[::N][:-1] + return abc, sectionStartTimes + else: + return abc ## matrix of all fit vectors per period + else: + if canUseMP==False: + raise RuntimeError("SineTools does not support Multi Processing for OS: "+str(os.name)) + print(int(multiprocessing.cpu_count() - 1)) + shmY = shared_memory.SharedMemory(create=True, size=y.nbytes) + shmT = shared_memory.SharedMemory(create=True, size=t.nbytes) + # Now create a NumPy array backed by shared memory + #to acces use this + tmpY = np.ndarray(y.shape, dtype=y.dtype, buffer=shmY.buf) + tmpT = np.ndarray(t.shape, dtype=t.dtype, buffer=shmT.buf) + tmpY[:] = y[:] + tmpT[:] = t[:] + mpdata['t_shared_name'] = shmT.name + mpdata['y_shared_name'] = shmY.name + mpdata['t_shared_type'] = y.dtype + mpdata['y_shared_type'] = t.dtype + mpdata['freq'] = f0 + mpdata['samples_per_block']=N + mpdata['nmumber_of_samples']=y.size + numCoresToBeUsed=int(multiprocessing.cpu_count()-1) + i = np.arange(M) + mpchunksize=10 + #iterMP3paramFit = partial(multiPrcoess3paramFit, mpdata=mpdata) + results = process_map(threeparsinefitNOJIT, i,chunksize=mpchunksize, max_workers=numCoresToBeUsed) + shmT.close() + shmT.unlink() # Free and release the shared memory block at the very end + shmY.close() + shmY.unlink() # Free and release the shared memory block at the very end + for i in np.arange(M): + abc[i,:]=results[i] + if returnSectionStartTimes: + sectionStartTimes = t[::N][:-1] + return abc, sectionStartTimes + else: + return abc + +# periodical sinefit at known frequency +@njit(parallel=True) +def seq_numbathreeparsinefit2(y, t, f0, periods=1,returnSectionStartTimes=False): + """ + period-wise sine-fit at a known frequency\n + y vector of sample values \n + t vector of sample times\n + f0 known frequency\n + periods number of full periods for every fit (default=1)\n + \n + returns a (n,3)-matrix of coefficient-triplets [[a,b,c], ...]\n + for y = a*sin(2*pi*f0*t) + b*cos(2*pi*f0*t) + c + """ + if y.size < t.size: + raise ValueError("Dimension mismatch in input data y<t") + if t.size < y.size: + RuntimeWarning("Dimension mismatch in input data y>t. fiting only for t.size values") + Tau = 1.0 / f0 *periods + dt = np.mean(np.diff(t)) + N = int(Tau // dt) ## samples per section + M = int(t.size // N) ## number of sections or periods + w0 = 2 * np.pi * f0 + abc = np.zeros((M, 3),dtype=np.float64) + for i in prange(int(M)): + ti = t[i * N : (i + 1) * N] + yi = y[i * N : (i + 1) * N] + a = np.ones((3, N), dtype=np.float64) + a[0, :] = np.cos(w0 * ti) + a[1, :] = np.sin(w0 * ti) + # a = np.array([np.cos(w0 * t), np.sin(w0 * t), np.ones(t.size)]) + aT = a.transpose() + abc[i, :] = np.linalg.lstsq(aT, yi)[0][0:3] + return abc ## matrix of all fit vectors per period + +@njit(parallel=True) +def seq_numbathreeparsinefit(y, t, f0, periods=1,returnSectionStartTimes=False): + """ + period-wise sine-fit at a known frequency\n + y vector of sample values \n + t vector of sample times\n + f0 known frequency\n + periods number of full periods for every fit (default=1)\n + \n + returns a (n,3)-matrix of coefficient-triplets [[a,b,c], ...]\n + for y = a*sin(2*pi*f0*t) + b*cos(2*pi*f0*t) + c + """ + if y.size < t.size: + raise ValueError("Dimension mismatch in input data y<t") + if t.size < y.size: + RuntimeWarning("Dimension mismatch in input data y>t. fiting only for t.size values") + Tau = 1.0 / f0 *periods + dt = np.mean(np.diff(t)) + N = int(Tau // dt) ## samples per section + M = int(t.size // N) ## number of sections or periods + + abc = np.zeros((M, 3),dtype=np.float64) + for i in prange(int(M)): + ti = t[i * N : (i + 1) * N] + yi = y[i * N : (i + 1) * N] + abc[i, :] = threeparsinefit(yi, ti, f0) + return abc ## matrix of all fit vectors per period + + + +def getFreqOffSetFromSeqThreeSineFitPhaseSlope(abc, sectionStartTimes): + Complex = abc[:, 1] + 1j * abc[:, 0] + nomalized = Complex / abs(Complex) # all vectors ar normalized + radialCord = np.sum(nomalized) / nomalized.size # all vectors are added to have one vector pointing in the direction of the mean value + deltaAng = np.angle(nomalized) - np.angle(radialCord) # differences of the angles this value can be bigger than -180 -- 180 deg + coef = np.polyfit(sectionStartTimes, deltaAng, 1) # dphi /dt + deltaF=coef[0]/(np.pi*2) + return deltaF + +# four parameter sine-fit (with frequency approximation) +def fourparsinefit(y, t, f0, tol=1.0e-7, nmax=1000): + """ + y sampled data values \n + t sample times of y \n + f0 estimate of sine frequency \n + tol rel. frequency correction where we stop \n + nmax maximum number of iterations taken \n + \n + returns the vector [a, b, c, w] of a*sin(w*t)+b*cos(w*t)+c + """ + abcd = threeparsinefit(y, t, f0) + w = 2 * np.pi * f0 + err = 1 + i = 0 + while (err > tol) and (i < nmax): + D = np.array( + [ + np.cos(w * t), + np.sin(w * t), + np.ones(t.size), + (-1.0) * abcd[0] * t * np.sin(w * t) + abcd[1] * t * np.cos(w * t), + ] + ) + abcd = (la.lstsq(D.transpose(), y))[0] + dw = abcd[3] + w = w + 0.9 * dw + i += 1 + err = np.absolute(dw / w) + + assert i < nmax, "iteration error" + + return np.hstack((abcd[0:3], w / (2 * np.pi))) + +# multiprocessing wraper for threeparsinefit +def multiPrcoess4paramFit(i): + """ + mpdata['t_shared_name'] = shmY.name + mpdata['y_shared_name'] = shmT.name + mpdata['t_shared_type'] = y.dtype + mpdata['y_shared_type'] = t.dtype + mpdata['freq'] = f0 + mpdata['samples_per_block'] = N + mpdata['nmumber_of_samples'] = y.size + """ + N=mpdata['samples_per_block'] + + existing_shm_t = shared_memory.SharedMemory(name=mpdata['t_shared_name']) + t = np.ndarray((mpdata['nmumber_of_samples'],), dtype=mpdata['t_shared_type'], buffer=existing_shm_t.buf) + ti = t[i * N: (i + 1) * N] + + existing_shm_y = shared_memory.SharedMemory(name=mpdata['y_shared_name']) + y = np.ndarray((mpdata['nmumber_of_samples'],), dtype=mpdata['y_shared_type'], buffer=existing_shm_y.buf) + yi = y[i * N: (i + 1) * N] + result=fourparsinefit(yi, ti, mpdata['freq']) + # Clean up + del t # Unnecessary; merely emphasizing the array is no longer used + del y # Unnecessary; merely emphasizing the array is no longer used + existing_shm_t.close() + existing_shm_y.close() + return result + + +def calc_fourparsine(abcf, t): + """ + return y = abc[0]*sin(2*pi*f0*t) + abc[1]*cos(2*pi*f0*t) + abc[2] + """ + w0 = 2 * np.pi * abcf[3] + return abcf[0] * np.cos(w0 * t) + abcf[1] * np.sin(w0 * t) + abcf[2] + + +""" +from octave ... +function abcw = fourParSinefit(data,w0) + abc = threeParSinefit(data,w0); + a=abc(1); + b=abc(2); + c=abc(3); + w = w0; + + do + D = [sin(w.*data(:,1)) , cos(w.*data(:,1)) , ones(rows(data),1) , a.*data(:,1).*cos(w.*data(:,1)) - b.*data(:,1).*sin(w.*data(:,1)) ]; + + s = D \ data(:,2); + dw = s(4); + w = w+0.9*dw; + err = abs(dw/w); + + until (err < 1.0e-8 ); + + abcw = [s(1),s(2),s(3),w]; + +endfunction +""" + +# periodical sinefit at known frequency +def seq_fourparsinefit(y, t, f0, tol=1.0e-7, nmax=1000, periods=1, debug_plot=False, multiTasking=True): + """ + sliced or period-wise sine-fit at a known frequency\n + y vector of sample values \n + t vector of sample times\n + f0 estimate of excitation frequency\n + periods integer number of periods used in each slice for fitting + nmax maximum of iteration to improve f0 \n + debug_plot Flag for plotting the sequential fit for dubugging \n + \n + returns a (n,3)-matrix of coefficient-triplets [[a,b,c], ...]\n + for y = a*sin(2*pi*f0*t) + b*cos(2*pi*f0*t) + c + """ + if y.size < t.size: + raise ValueError("Dimension mismatch in input data y<t") + if t.size < y.size: + warnings.warn( + "Dimension mismatch in input data y>t. fiting only for t.size values", + RuntimeWarning, + ) + Tau = 1.0 / f0 *periods + dt = np.mean(np.diff(t)) + N = int(Tau // dt) ## samples per section + M = int(t.size // N) ## number of sections or periods + + abcd = np.zeros((M, 4)) + if multiTasking == False: + for i in range(M): + ti = t[i * N : (i + 1) * N] + yi = y[i * N : (i + 1) * N] + abcd[i, :] = fourparsinefit(yi, ti, f0, tol=tol, nmax=nmax) + else: + shmY = shared_memory.SharedMemory(create=True, size=y.nbytes) + shmT = shared_memory.SharedMemory(create=True, size=t.nbytes) + # Now create a NumPy array backed by shared memory + #to accesdta use this + tmpY = np.ndarray(y.shape, dtype=y.dtype, buffer=shmY.buf) + tmpT = np.ndarray(t.shape, dtype=t.dtype, buffer=shmT.buf) + tmpY[:] = y[:] + tmpT[:] = t[:] + mpdata['t_shared_name'] = shmT.name + mpdata['y_shared_name'] = shmY.name + mpdata['t_shared_type'] = y.dtype + mpdata['y_shared_type'] = t.dtype + mpdata['freq'] = f0 + mpdata['samples_per_block']=N + mpdata['nmumber_of_samples']=y.size + numCoresToBeUsed=int(multiprocessing.cpu_count()-1) + + i = np.arange(M) + mpchunksize=10 + results = process_map(multiPrcoess4paramFit, i,chunksize=mpchunksize, max_workers=numCoresToBeUsed) + shmT.close() + shmT.unlink() # Free and release the shared memory block at the very end + shmY.close() + shmY.unlink() # Free and release the shared memory block at the very end + for i in np.arange(M): + abcd[i,:]=results[i] + return abcd + + if debug_plot: + mp.ioff() + fig = mp.figure("seq_fourparsinefit") + fig.clear() + p1 = fig.add_subplot(211) + p2 = fig.add_subplot(212, sharex=p1) + + for i in range(M): + p1.plot(t[i * N : (i + 1) * N], y[i * N : (i + 1) * N], ".") + s = calc_fourparsine( + abcd[i, :], t[i * N : (i + 1) * N] + ) # fitted data to plot + p1.plot(t[i * N : (i + 1) * N], s, "-") + r = y[i * N : (i + 1) * N] - s # residuals to plot + p2.plot(t[i * N : (i + 1) * N], r, ".") + yi = y[i * N : (i + 1) * N] + mp.show() + + return abcd ## matrix of all fit vectors per period + +# fitting a pseudo-random multi-sine signal with 2*Nf+1 parameters +def multi_threeparsinefit(y, t, f0): # f0 vector of frequencies + """ + fit a time series of a sum of sine-waveforms with a given set of frequencies\n + y vector of sample values \n + t vector of sample times\n + f0 vector of known frequencies\n + \n + returns a vector of coefficient-triplets [a,b,c] for the frequencies\n + for y = sum_i (a_i*sin(2*pi*f0_i*t) + b_i*cos(2*pi*f0_i*t) + c_i + """ + w0 = 2 * np.pi * f0 + D = np.ones((len(t), 1)) # for the bias + # set up design matrix + for w in w0[::-1]: + D = np.hstack((np.cos(w * t)[:, np.newaxis], np.sin(w * t)[:, np.newaxis], D)) + + abc = np.linalg.lstsq(D, y) + return abc[0] ## fit vector a*cos+b*sin+c + + +def multi_amplitude(abc): # abc = [a1,b1 , a2,b2, ...,bias] + """ + return the amplitudes given the coefficients of a multi-sine\n + abc = [a1,b1 , a2,b2, ...,bias] \n + y = sum_i (a_i*sin(2*pi*f0_i*t) + b_i*cos(2*pi*f0_i*t) + c_i + """ + x = abc[:-1][1::2] + 1j * abc[:-1][0::2] # make complex without Bias (last value) + return np.absolute(x) + + +def multi_phase(abc, deg=False): # abc = [bias, a1,b1 , a2,b2, ...] + """ + return the initial phases given the coefficients of a multi-sine\n + abc = [a1,b1 , a2,b2, ...,bias] \n + y = sum_i (a_i*sin(2*pi*f0_i*t) + b_i*cos(2*pi*f0_i*t) + c_i + """ + x = abc[:-1][1::2] + 1j * abc[:-1][0::2] # make complex without Bias (last value) + return np.angle(x, deg=deg) + + +def multi_waveform_abc(f, abc, t): + """ + generate a sample time series of a multi-sine from coefficients and frequencies\n + f vector of given frequencies \n + abc = [a1,ba, a2,b2, ..., bias]\n + t vector of sample times t_i\n + \n + returns the vector \n + y = sum_i (a_i*sin(2*pi*f0_i*t) + b_i*cos(2*pi*f0_i*t) + bias + """ + ret = 0.0 * t + abc[-1] # bias + for fi, a, b in zip(f, abc[0::2], abc[1::2]): + ret = ret + a * np.cos(2 * np.pi * fi * t) + b * np.sin(2 * np.pi * fi * t) + return ret + + +def multi_waveform_mp(f, m, p, t, bias=0, deg=True): + """ + generate a sample time series of a multi-sine from magnitude/phase and frequencies\n + f vector of given frequencies \n + m vector of magnitudes \n + p vector of phases \n + t vector of sample times t_i \n + bias scalar value for total bias + deg = boolean whether phase is in degree + \n + returns the vector \n + y = sum_i (a_i*sin(2*pi*f0_i*t) + b_i*cos(2*pi*f0_i*t) + bias + """ + ret = 0.0 * t + bias # init and bias + if deg: + p = np.deg2rad(p) + for fi, m_i, p_i in zip(f, m, p): + ret = ret + m_i * np.sin(2 * np.pi * fi * t + p_i) + return ret + + +################################## +# Counter based stuff + +# periodical sinefit to the linearly increasing heterodyne counter +# version based on Blume +def seq_threeparcounterfit(y, t, f0, diff=False): + """ + period-wise (single-)sinefit to the linearly increasing heterodyne counter + version based on "Blume et al. "\n + y vector of sampled counter values + t vector of sample times + f given frequency\n + \n + returns (n,3)-matrix of coefficient-triplets [a,b,c] per period \n + + if diff=True use differentiation to remove carrier (c.f. source) + """ + Tau = 1.0 / f0 + dt = t[1] - t[0] + N = int(np.floor(Tau / dt)) ## samples per section + M = int(np.floor(t.size / N)) ## number of sections or periods + + remove_counter_carrier(y, diff=diff) + + abc = np.zeros((M, 4)) + + for i in range(int(M)): + ti = t[i * N : (i + 1) * N] + yi = y[i * N : (i + 1) * N] + + abc[i, :] = threeparsinefit_lin(yi, ti, f0) + return abc ## matrix of all fit vectors per period + + +def remove_counter_carrier(y, diff=False): + """ + remove the linear increase in the counter signal + generated by the carrier frequency of a heterodyne signal\n + y vector of samples of the signal + """ + if diff: + d = np.diff(y) + d = d - np.mean(d) + y = np.hstack((0, np.cumsum(d))) + else: + slope = y[-1] - y[0] # slope of linear increment + y = y - slope * np.linspace( + 0.0, 1.0, len(y), endpoint=False + ) # removal of linear increment + return y + + +# calculate displacement and acceleration to the same analytical s(t) +# Bsp: fm = 2e7, f=10, s0=0.15, phi0=np.pi/3, ti, drift=0.03, ampdrift=0.03,thd=[0,0.02,0,0.004] +def disp_acc_distorted(fm, f, s0, phi0, ti, drift=0, ampdrift=0, thd=0): + """ + calculate the respective (displacement-) counter and acceleration + for a parmeterized distorted sine-wave motion in order to compare accelerometry with interferometry \n + fm is heterodyne carrier frequency (after mixing)\n + f is mechanical sine frequency (nominal) \n + phi_0 accelerometer phase delay \n + ti vector of sample times \n + drift is displacement zero drift \n + ampdrift is displacement amplitude druft\n + thd is vector of higher harmonic amplitudes (c.f. source) + """ + om = 2 * np.pi * f + om2 = om ** 2 + tau = ti[-1] - ti[0] + disp = np.sin(om * ti + phi0) + if thd != 0: + i = 2 + for h in thd: + disp = disp + h * np.sin(i * om * ti + phi0) + i = i + 1 + if ampdrift != 0: + disp = disp * (1 + ampdrift / tau * ti) + if drift != 0: + disp = disp + s0 * drift / tau * ti + disp = disp * s0 + disp = np.floor((disp * 2 / 633e-9) + fm * ti) + + acc = -s0 * om2 * (1 + ampdrift / tau * ti) * np.sin(om * ti + phi0) + if ampdrift != 0: + acc = acc + (2 * ampdrift * s0 * om * np.cos(om * ti + phi0)) / tau + if thd != 0: + i = 2 + for h in thd: + acc = acc - s0 * h * om2 * (1 + ampdrift / tau * ti) * i ** 2 * np.sin( + i * om * ti + phi0 + ) + if ampdrift != 0: + acc = ( + acc + + (2 * ampdrift * s0 * om * i * h * np.cos(om * ti + phi0)) / tau + ) + i = i + 1 + + return disp, acc + + +################################### +# Generation and adaptation of Parameters of the Multi-Sine considering hardware constraints +def PR_MultiSine_adapt( + f1, + Nperiods, + Nsamples, + Nf=8, + fs_min=0, + fs_max=1e9, + frange=10, + log=True, + phases=None, + sample_inkr=1, +): + """ + Returns an additive normalized Multisine time series. \n + f1 = start frequency (may be adapted by the algorithm) \n + Nperiods = number of periods of f1 (may be increased by algorithm) \n + Nsamples = Minimum Number of samples \n + Nf = number of frequencies in multi frequency mix \n + fs_min = minimum sample rate of used device (default 0) \n + fs_max = maximum sample rate of used device (default 0) \n + frange = range of frequency as a factor relative to f1 (default 10 = decade) \n + log = boolean for logarithmic (True, default) or linear (False) frequency scale \n + phases = float array of given phases for the frequencies (default=None=random) \n + deg= boolean for return phases in deg (True) or rad (False) \n + sample_inkr = minimum block of samples to add to a waveform + \n + returns: freq,phase,fs,ti,multi \n + freq= array of frequencies \n + phase=used phases in deg or rad \n + T1 = (adapted) duration + fs=sample rate \n + """ + if ( + Nsamples // sample_inkr * sample_inkr != Nsamples + ): # check multiplicity of sample_inkr + Nsamples = ( + Nsamples // sample_inkr + 1 + ) * sample_inkr # round to next higher multiple + + T0 = Nperiods / f1 # given duration + fs0 = Nsamples / T0 # (implicitly) given sample rate + + if False: + print("0 Nperiods: " + str(Nperiods)) + print("0 Nsamples: " + str(Nsamples)) + print("0 fs: " + str(fs0)) + print("0 T0: " + str(T0)) + print("0 f1: " + str(f1)) + + fs = fs0 + if fs0 < fs_min: # sample rate too low, then set to minimum + fs = fs_min + print("sample rate increased") + elif fs0 > fs_max: # sample rate too high, set to max-allowed and + fs = fs_max + Nperiods = np.ceil( + Nperiods * fs0 / fs_max + ) # increase number of periods to get at least Nsamples samples + T0 = Nperiods / f1 + print("sample rate reduced, Nperiods=" + str(Nperiods)) + + Nsamples = T0 * fs + if ( + Nsamples // sample_inkr * sample_inkr != Nsamples + ): # check multiplicity of sample_inkr + Nsamples = ( + Nsamples // sample_inkr + 1 + ) * sample_inkr # round to next higher multiple + + T1 = Nsamples / fs # adapt exact duration + f1 = Nperiods / T1 # adapt f1 for complete cycles + if False: + print("Nperiods: " + str(Nperiods)) + print("Nsamples: " + str(Nsamples)) + print("fs: " + str(fs)) + print("T1: " + str(T1)) + print("f1: " + str(f1)) + + f_res = 1 / T1 # frequency resolution + # determine a series of frequencies (freq[]) + if log: + fact = np.power(frange, 1.0 / (Nf - 1)) # factor for logarithmic scale + freq = f1 * np.power(fact, np.arange(Nf)) + else: + step = (frange - 1) * f1 / (Nf - 1) + freq = np.arange(f1, frange * f1 + step, step) + + # auxiliary function to find the nearest available frequency + def find_nearest( + x, possible + ): # match the theoretical freqs to the possible periodic freqs + idx = (np.absolute(possible - x)).argmin() + return possible[idx] + + fi_pos = np.arange(f1, frange * f1 + f_res, f_res) # possible periodic frequencies + f_real = [] + for f in freq: + f_real.append(find_nearest(f, fi_pos)) + freq = np.hstack(f_real) + if True: + print("freq: " + str(freq)) + + if phases is None: # generate random phases + phase = np.random.randn(Nf) * 2 * np.pi # random phase + else: # use given phases + phase = phases + + return freq, phase, T1, fs + + +################################### +# Pseudo-Random-MultiSine for "quick calibration" +def PR_MultiSine( + f1, + Nperiods, + Nsamples, + Nf=8, + fs_min=0, + fs_max=1e9, + frange=10, + log=True, + phases=None, + deg=False, + sample_inkr=1, +): + """ + Returns an additive normalized Multisine time series. \n + f1 = start frequency (may be adapted) \n + Nperiods = number of periods of f1 (may be increased) \n + Nsamples = Minimum Number of samples \n + Nf = number of frequencies in multi frequency mix \n + fs_min = minimum sample rate of used device (default 0) \n + fs_max = maximum sample rate of used device (default 0) \n + frange = range of frequency as a factor relative to f1 (default 10 = decade) \n + log = boolean for logarithmic (True, default) or linear (False) frequency scale \n + phases = float array of given phases for the frequencies (default=None=random) \n + deg= boolean for return phases in deg (True) or rad (False) \n + sample_inkr = minimum block of samples to add to a waveform + \n + returns: freq,phase,fs,ti,multi \n + freq= array of frequencies \n + phase=used phases in deg or rad \n + fs=sample rate \n + ti=timestamps \n + multi=array of time series values \n + """ + + freq, phase, T1, fs = PR_MultiSine_adapt( + f1, + Nperiods, + Nsamples, + Nf=Nf, + fs_min=fs_min, + fs_max=fs_max, + frange=frange, + log=log, + phases=phases, + sample_inkr=sample_inkr, + ) + + if deg: # rad -> deg + phase = phase * np.pi / 180.0 + + ti = np.arange(T1 * fs, dtype=np.float32) / fs + + multi = np.zeros(len(ti), dtype=np.float64) + for f, p in zip(freq, phase): + multi = multi + np.sin(2 * np.pi * f * ti + p) + + multi = multi / np.amax(np.absolute(multi)) # normalize + + if False: + import matplotlib.pyplot as mp + + fig = mp.figure(1) + fig.clear() + pl1 = fig.add_subplot(211) + pl2 = fig.add_subplot(212) + pl1.plot(ti, multi, "-o") + pl2.plot(np.hstack((ti, ti + ti[-1] + ti[1])), np.hstack((multi, multi)), "-o") + mp.show() + + return ( + freq, + phase, + fs, + multi, + ) # frequency series, sample rate, sample timestamps, waveform + + +def seq_multi_threeparam_Dmatrix(f,t,periods=1, progressive=True): + """ + Fit a multi-sinus-signal in slices in one go. + + Parameters + ---------- + f : numpy.array of floats + list of frequencies in the signal + t : numpy.array of floats + timestamps of y (typically seconds) + periods : float, optional + the number of periods of each frequncy used for each fit. The default is 1. + + Returns + ------- + + fr, : 1-d numpy.array of floats + frequencies related to slices + D : 2-d numpy.array of floats + Design matrix for the fit + + """ + T = t[-1]-t[0] + col = 0 + data = np.array([]) + ci = [] + ri = [] + fr = [] + # Designmatrix for sin/cos + for fi, omi in zip(f,2*np.pi*f): + Nri = 0 # counter for the current row index + if progressive: + tau = np.ceil(f[0]/fi*periods)/fi # approximately same abs. slice length for all fi + else: + tau = 1/fi*periods # slice length in seconds for periods of fi + t_sl = np.array_split(t,np.ceil(T/tau)) # array of slices of sample times + fr = fr + [fi]*len(t_sl) # len(t_sl) times frequency fi in Design matrix + for ti in t_sl: # loop over slices fo one frequncy + sin = np.sin(omi * ti) + cos = np.cos(omi * ti) + + data = np.hstack((data,cos)) # data vector + ci = ci+[col]*len(ti) # column index vector + col = col+1 + ri = ri+ [Nri+i for i in range(len(ti))] # row index vector + + data = np.hstack((data,sin)) # data vector + ci = ci+[col]*len(ti) # column index vector + col = col+1 + ri = ri+ [Nri+i for i in range(len(ti))] # row index vector + + Nri = Nri+len(ti) + + # Designmatrix for Bias + data = np.hstack((data,np.ones((len(t))))) + ci = ci = ci+[col]*len(t) + ri = ri + [i for i in range(len(t))] + col = col+1 + + # build sparse matrix, init as coo, map to csr + D = coo_matrix((data,(ri,ci)),shape=(len(t),col)).tocsr() + return np.array(fr), D + +def seq_multi_threeparsinefit(f,y,t,periods=1, D_fr=None, abc0=None, progressive=True): + """ + performs a simultanius, sliced three-parameter fit on a multisine signal y + + Parameters + ---------- + f : 1-d numpy.array of floats + frequencies in the signal + y : 1-d numpy.array + samples of the signal + t : 1-d numpy.array + vector of sample times + periods : float + (fractional) number of periods in one slice. The default is 1. + + Returns + ------- + f_ab_c : 2-d numpy.array of floats + frequencies, a/b-coefficients and bias related to given frequencies [f1,f1,f1,f2,f2, ...fn, fn, f=0=bias] + y0 : 1-d numpy.array of floats + samples of the optimal fit + resid : 1-d numpy.array of floats + Residuals = Difference =(y-y0) + + """ + if D_fr is None: + fr, D = seq_multi_threeparam_Dmatrix(f,t,periods, progressive=progressive) # calculate the design matrix (as sparse matrix) + else: + D = D_fr[0] + fr = D_fr[1] + + abc = sp_lsqr(D,y,x0=abc0, atol=1.0e-9, btol=1.0e-9) + y0 = D.dot(abc[0]) + + #print(abc) + f_ab_c =[] + k=0 + for fi in fr: + f_ab_c = f_ab_c+ [fi, abc[0][k],abc[0][k+1]] + k=k+2 + f_ab_c = f_ab_c + [0.0,abc[0][-1],0.0] + f_ab_c = np.array(f_ab_c).reshape((len(f_ab_c)//3,3)) + + #print(f_ab_c) + return f_ab_c, y0, y-y0 # coefficients, fit signal, residuals + + +def seq_multi_fourparam_sineFit(f,y,t,periods=1, tol=1.0e-6, n_max=100): + """ + Fit a multi-sinus-signal in slices in one go. With a frequency correction + universal for all frequncies (time-stretch) + + Parameters + ---------- + f : numpy.array of floats + list of frequencies in the signal + t : numpy.array of floats + timestamps of y (typically seconds) + periods : float, optional + the number of periods of each frequncy used for each fit. The default is 1. + tol : float, optional + frequency correction factor (1-tol) when convergence is assumed default 1.0e-6 + n_max : uint, optional + maximum number of iterations performed. + + Returns + ------- + + fr, : 1-d numpy.array of floats + frequencies related to slices + D : 2-d numpy.array of floats + Design matrix for the fit + + """ + T = t[-1]-t[0] + dw = 1.0 + lam=1.0 + itr=0 + while (np.abs(dw) > tol) and (itr<n_max): + itr += 1 + # Designmatrix for three paramsin/cos + fr,D = seq_multi_threeparam_Dmatrix(f, t, periods=periods, progressive=False) + # initial ab + f_ab_c = seq_multi_threeparsinefit(f, y, t, D_fr=(D,fr))[0] + last_col = np.zeros(len(t)) + for fi, omi in zip(f,2*np.pi*f): + tau = 1/fi*periods # slice length in seconds + t_sl = np.array_split(t,np.ceil(T/tau)) # array of slices of sample times + sin_cos=np.array([]) + for i,ti in enumerate(t_sl): # loop over slices for one frequncy + si = -f_ab_c[i,1]*omi*ti*np.sin(omi * ti) + co = f_ab_c[i,2]*omi*ti*np.cos(omi * ti) + + sin_cos = np.hstack((sin_cos,si+co)) # data vector + last_col = last_col+sin_cos + last_col = coo_matrix((last_col,([i for i in range(len(last_col))],[0 for i in range(len(last_col))])),shape=(len(last_col),1)).tocsr() + D = hstack((D, last_col)) + abc = sp_lsqr(D,y)[0] + dw = abc[-1] + f = (1-0.9*dw)*f + lam=lam*(1-dw) + if False: + print("abc:") + print(abc) + print("f=") + print(f) + print("cor. = %f" %(lam)) + print("took %d iterations" % (itr)) + f_ab_c, y0, res = seq_multi_threeparsinefit(f, y, t, D_fr=None) + + return f_ab_c, y0, y-y0 # coefficients, fit signal, residuals + +def seq_multi_amplitude(f_ab_c): + """ + return the amplitude(s) of a sequentially fit multi-sine signal, + i.e. amplitudes encoded in the return value of seq_multi_threeparsinefit. + + Parameters + ---------- + f_ab_c : 2-d numpy array of floats (Nx3) + f,a,b in a row for several rows, as returned by seq_multi_threeparsinefit. + + Returns + ------- + 2d-numpy-array of floats (Nx2) + frequency and associated amplitude. + + """ + #print("Test") + N = f_ab_c.shape[0]-1 + amp = np.array([1j*np.abs(f_ab_c[i,1]+f_ab_c[i][2]) for i in range(f_ab_c.shape[0]-1)]).reshape((N,1)) + return np.hstack((f_ab_c[:-1,0].reshape((N,1)), amp)) + +def seq_multi_phase(f_ab_c, deg=True): + """ + Calculates the initial phase of a sequentially fit multi-sine signal, + i.e. initial phases encoded in the return value of seq_multi_threeparsinefit. + result is either in degrees (deg=True) or rad (deg=False). + + Parameters + ---------- + f_ab_c : 2-d numpy array of floats (Nx3) + f,a,b in a row for several rows, as returned by seq_multi_threeparsinefit. + deg : Boolean, optional + Flag whether result is in degrees or rad. The default is True (Degrees). + + Returns + ------- + 2d-numpy-array of floats (Nx2) + frequency and associated initial phase. + + """ + N = f_ab_c.shape[0]-1 + phase = np.array([1j*np.angle(f_ab_c[i,1]+f_ab_c[i][2],deg=deg) for i in range(f_ab_c.shape[0]-1)]).reshape((N,1)) + return np.hstack((f_ab_c[:-1,0].reshape((N,1)), phase)) + +def seq_multi_bias(f_ab_c): + """ + Returns the single bias of a sequentially fit multi-sine signal, + i.e. bias encoded in the return value of seq_multi_threeparsinefit. + + Parameters + ---------- + f_ab_c : 2-d numpy array of floats (Nx3) + f,a,b in a row for several rows, as returned by seq_multi_threeparsinefit. + + Returns + ------- + float + Bias of the signal. + + """ + return f_ab_c[-1][1] + +if __name__ == '__main__': + print("Sine Tools is not Intended to be run as main Programm") + sampleTimes = sampletimes(20e6, 25) + sw = sinewave(1000, 1, 0.1, sampleTimes) + seq_numbathreeparsinefit(sw, sampleTimes, 1000, periods=100) + threeparsinefit.parallel_diagnostics(level=4) + seq_numbathreeparsinefit.parallel_diagnostics(level=4)