This commit is contained in:
rysertio 2023-10-30 06:52:43 +06:00
commit e010ebd84e
34 changed files with 97931 additions and 0 deletions

1
.gitignore vendored Normal file
View File

@ -0,0 +1 @@
build/

674
LICENSE.txt Normal file
View File

@ -0,0 +1,674 @@
GNU GENERAL PUBLIC LICENSE
Version 3, 29 June 2007
Copyright (C) 2007 Free Software Foundation, Inc. <https://fsf.org/>
Everyone is permitted to copy and distribute verbatim copies
of this license document, but changing it is not allowed.
Preamble
The GNU General Public License is a free, copyleft license for
software and other kinds of works.
The licenses for most software and other practical works are designed
to take away your freedom to share and change the works. By contrast,
the GNU General Public License is intended to guarantee your freedom to
share and change all versions of a program--to make sure it remains free
software for all its users. We, the Free Software Foundation, use the
GNU General Public License for most of our software; it applies also to
any other work released this way by its authors. You can apply it to
your programs, too.
When we speak of free software, we are referring to freedom, not
price. Our General Public Licenses are designed to make sure that you
have the freedom to distribute copies of free software (and charge for
them if you wish), that you receive source code or can get it if you
want it, that you can change the software or use pieces of it in new
free programs, and that you know you can do these things.
To protect your rights, we need to prevent others from denying you
these rights or asking you to surrender the rights. Therefore, you have
certain responsibilities if you distribute copies of the software, or if
you modify it: responsibilities to respect the freedom of others.
For example, if you distribute copies of such a program, whether
gratis or for a fee, you must pass on to the recipients the same
freedoms that you received. You must make sure that they, too, receive
or can get the source code. And you must show them these terms so they
know their rights.
Developers that use the GNU GPL protect your rights with two steps:
(1) assert copyright on the software, and (2) offer you this License
giving you legal permission to copy, distribute and/or modify it.
For the developers' and authors' protection, the GPL clearly explains
that there is no warranty for this free software. For both users' and
authors' sake, the GPL requires that modified versions be marked as
changed, so that their problems will not be attributed erroneously to
authors of previous versions.
Some devices are designed to deny users access to install or run
modified versions of the software inside them, although the manufacturer
can do so. This is fundamentally incompatible with the aim of
protecting users' freedom to change the software. The systematic
pattern of such abuse occurs in the area of products for individuals to
use, which is precisely where it is most unacceptable. Therefore, we
have designed this version of the GPL to prohibit the practice for those
products. If such problems arise substantially in other domains, we
stand ready to extend this provision to those domains in future versions
of the GPL, as needed to protect the freedom of users.
Finally, every program is threatened constantly by software patents.
States should not allow patents to restrict development and use of
software on general-purpose computers, but in those that do, we wish to
avoid the special danger that patents applied to a free program could
make it effectively proprietary. To prevent this, the GPL assures that
patents cannot be used to render the program non-free.
The precise terms and conditions for copying, distribution and
modification follow.
TERMS AND CONDITIONS
0. Definitions.
"This License" refers to version 3 of the GNU General Public License.
"Copyright" also means copyright-like laws that apply to other kinds of
works, such as semiconductor masks.
"The Program" refers to any copyrightable work licensed under this
License. Each licensee is addressed as "you". "Licensees" and
"recipients" may be individuals or organizations.
To "modify" a work means to copy from or adapt all or part of the work
in a fashion requiring copyright permission, other than the making of an
exact copy. The resulting work is called a "modified version" of the
earlier work or a work "based on" the earlier work.
A "covered work" means either the unmodified Program or a work based
on the Program.
To "propagate" a work means to do anything with it that, without
permission, would make you directly or secondarily liable for
infringement under applicable copyright law, except executing it on a
computer or modifying a private copy. Propagation includes copying,
distribution (with or without modification), making available to the
public, and in some countries other activities as well.
To "convey" a work means any kind of propagation that enables other
parties to make or receive copies. Mere interaction with a user through
a computer network, with no transfer of a copy, is not conveying.
An interactive user interface displays "Appropriate Legal Notices"
to the extent that it includes a convenient and prominently visible
feature that (1) displays an appropriate copyright notice, and (2)
tells the user that there is no warranty for the work (except to the
extent that warranties are provided), that licensees may convey the
work under this License, and how to view a copy of this License. If
the interface presents a list of user commands or options, such as a
menu, a prominent item in the list meets this criterion.
1. Source Code.
The "source code" for a work means the preferred form of the work
for making modifications to it. "Object code" means any non-source
form of a work.
A "Standard Interface" means an interface that either is an official
standard defined by a recognized standards body, or, in the case of
interfaces specified for a particular programming language, one that
is widely used among developers working in that language.
The "System Libraries" of an executable work include anything, other
than the work as a whole, that (a) is included in the normal form of
packaging a Major Component, but which is not part of that Major
Component, and (b) serves only to enable use of the work with that
Major Component, or to implement a Standard Interface for which an
implementation is available to the public in source code form. A
"Major Component", in this context, means a major essential component
(kernel, window system, and so on) of the specific operating system
(if any) on which the executable work runs, or a compiler used to
produce the work, or an object code interpreter used to run it.
The "Corresponding Source" for a work in object code form means all
the source code needed to generate, install, and (for an executable
work) run the object code and to modify the work, including scripts to
control those activities. However, it does not include the work's
System Libraries, or general-purpose tools or generally available free
programs which are used unmodified in performing those activities but
which are not part of the work. For example, Corresponding Source
includes interface definition files associated with source files for
the work, and the source code for shared libraries and dynamically
linked subprograms that the work is specifically designed to require,
such as by intimate data communication or control flow between those
subprograms and other parts of the work.
The Corresponding Source need not include anything that users
can regenerate automatically from other parts of the Corresponding
Source.
The Corresponding Source for a work in source code form is that
same work.
2. Basic Permissions.
All rights granted under this License are granted for the term of
copyright on the Program, and are irrevocable provided the stated
conditions are met. This License explicitly affirms your unlimited
permission to run the unmodified Program. The output from running a
covered work is covered by this License only if the output, given its
content, constitutes a covered work. This License acknowledges your
rights of fair use or other equivalent, as provided by copyright law.
You may make, run and propagate covered works that you do not
convey, without conditions so long as your license otherwise remains
in force. You may convey covered works to others for the sole purpose
of having them make modifications exclusively for you, or provide you
with facilities for running those works, provided that you comply with
the terms of this License in conveying all material for which you do
not control copyright. Those thus making or running the covered works
for you must do so exclusively on your behalf, under your direction
and control, on terms that prohibit them from making any copies of
your copyrighted material outside their relationship with you.
Conveying under any other circumstances is permitted solely under
the conditions stated below. Sublicensing is not allowed; section 10
makes it unnecessary.
3. Protecting Users' Legal Rights From Anti-Circumvention Law.
No covered work shall be deemed part of an effective technological
measure under any applicable law fulfilling obligations under article
11 of the WIPO copyright treaty adopted on 20 December 1996, or
similar laws prohibiting or restricting circumvention of such
measures.
When you convey a covered work, you waive any legal power to forbid
circumvention of technological measures to the extent such circumvention
is effected by exercising rights under this License with respect to
the covered work, and you disclaim any intention to limit operation or
modification of the work as a means of enforcing, against the work's
users, your or third parties' legal rights to forbid circumvention of
technological measures.
4. Conveying Verbatim Copies.
You may convey verbatim copies of the Program's source code as you
receive it, in any medium, provided that you conspicuously and
appropriately publish on each copy an appropriate copyright notice;
keep intact all notices stating that this License and any
non-permissive terms added in accord with section 7 apply to the code;
keep intact all notices of the absence of any warranty; and give all
recipients a copy of this License along with the Program.
You may charge any price or no price for each copy that you convey,
and you may offer support or warranty protection for a fee.
5. Conveying Modified Source Versions.
You may convey a work based on the Program, or the modifications to
produce it from the Program, in the form of source code under the
terms of section 4, provided that you also meet all of these conditions:
a) The work must carry prominent notices stating that you modified
it, and giving a relevant date.
b) The work must carry prominent notices stating that it is
released under this License and any conditions added under section
7. This requirement modifies the requirement in section 4 to
"keep intact all notices".
c) You must license the entire work, as a whole, under this
License to anyone who comes into possession of a copy. This
License will therefore apply, along with any applicable section 7
additional terms, to the whole of the work, and all its parts,
regardless of how they are packaged. This License gives no
permission to license the work in any other way, but it does not
invalidate such permission if you have separately received it.
d) If the work has interactive user interfaces, each must display
Appropriate Legal Notices; however, if the Program has interactive
interfaces that do not display Appropriate Legal Notices, your
work need not make them do so.
A compilation of a covered work with other separate and independent
works, which are not by their nature extensions of the covered work,
and which are not combined with it such as to form a larger program,
in or on a volume of a storage or distribution medium, is called an
"aggregate" if the compilation and its resulting copyright are not
used to limit the access or legal rights of the compilation's users
beyond what the individual works permit. Inclusion of a covered work
in an aggregate does not cause this License to apply to the other
parts of the aggregate.
6. Conveying Non-Source Forms.
You may convey a covered work in object code form under the terms
of sections 4 and 5, provided that you also convey the
machine-readable Corresponding Source under the terms of this License,
in one of these ways:
a) Convey the object code in, or embodied in, a physical product
(including a physical distribution medium), accompanied by the
Corresponding Source fixed on a durable physical medium
customarily used for software interchange.
b) Convey the object code in, or embodied in, a physical product
(including a physical distribution medium), accompanied by a
written offer, valid for at least three years and valid for as
long as you offer spare parts or customer support for that product
model, to give anyone who possesses the object code either (1) a
copy of the Corresponding Source for all the software in the
product that is covered by this License, on a durable physical
medium customarily used for software interchange, for a price no
more than your reasonable cost of physically performing this
conveying of source, or (2) access to copy the
Corresponding Source from a network server at no charge.
c) Convey individual copies of the object code with a copy of the
written offer to provide the Corresponding Source. This
alternative is allowed only occasionally and noncommercially, and
only if you received the object code with such an offer, in accord
with subsection 6b.
d) Convey the object code by offering access from a designated
place (gratis or for a charge), and offer equivalent access to the
Corresponding Source in the same way through the same place at no
further charge. You need not require recipients to copy the
Corresponding Source along with the object code. If the place to
copy the object code is a network server, the Corresponding Source
may be on a different server (operated by you or a third party)
that supports equivalent copying facilities, provided you maintain
clear directions next to the object code saying where to find the
Corresponding Source. Regardless of what server hosts the
Corresponding Source, you remain obligated to ensure that it is
available for as long as needed to satisfy these requirements.
e) Convey the object code using peer-to-peer transmission, provided
you inform other peers where the object code and Corresponding
Source of the work are being offered to the general public at no
charge under subsection 6d.
A separable portion of the object code, whose source code is excluded
from the Corresponding Source as a System Library, need not be
included in conveying the object code work.
A "User Product" is either (1) a "consumer product", which means any
tangible personal property which is normally used for personal, family,
or household purposes, or (2) anything designed or sold for incorporation
into a dwelling. In determining whether a product is a consumer product,
doubtful cases shall be resolved in favor of coverage. For a particular
product received by a particular user, "normally used" refers to a
typical or common use of that class of product, regardless of the status
of the particular user or of the way in which the particular user
actually uses, or expects or is expected to use, the product. A product
is a consumer product regardless of whether the product has substantial
commercial, industrial or non-consumer uses, unless such uses represent
the only significant mode of use of the product.
"Installation Information" for a User Product means any methods,
procedures, authorization keys, or other information required to install
and execute modified versions of a covered work in that User Product from
a modified version of its Corresponding Source. The information must
suffice to ensure that the continued functioning of the modified object
code is in no case prevented or interfered with solely because
modification has been made.
If you convey an object code work under this section in, or with, or
specifically for use in, a User Product, and the conveying occurs as
part of a transaction in which the right of possession and use of the
User Product is transferred to the recipient in perpetuity or for a
fixed term (regardless of how the transaction is characterized), the
Corresponding Source conveyed under this section must be accompanied
by the Installation Information. But this requirement does not apply
if neither you nor any third party retains the ability to install
modified object code on the User Product (for example, the work has
been installed in ROM).
The requirement to provide Installation Information does not include a
requirement to continue to provide support service, warranty, or updates
for a work that has been modified or installed by the recipient, or for
the User Product in which it has been modified or installed. Access to a
network may be denied when the modification itself materially and
adversely affects the operation of the network or violates the rules and
protocols for communication across the network.
Corresponding Source conveyed, and Installation Information provided,
in accord with this section must be in a format that is publicly
documented (and with an implementation available to the public in
source code form), and must require no special password or key for
unpacking, reading or copying.
7. Additional Terms.
"Additional permissions" are terms that supplement the terms of this
License by making exceptions from one or more of its conditions.
Additional permissions that are applicable to the entire Program shall
be treated as though they were included in this License, to the extent
that they are valid under applicable law. If additional permissions
apply only to part of the Program, that part may be used separately
under those permissions, but the entire Program remains governed by
this License without regard to the additional permissions.
When you convey a copy of a covered work, you may at your option
remove any additional permissions from that copy, or from any part of
it. (Additional permissions may be written to require their own
removal in certain cases when you modify the work.) You may place
additional permissions on material, added by you to a covered work,
for which you have or can give appropriate copyright permission.
Notwithstanding any other provision of this License, for material you
add to a covered work, you may (if authorized by the copyright holders of
that material) supplement the terms of this License with terms:
a) Disclaiming warranty or limiting liability differently from the
terms of sections 15 and 16 of this License; or
b) Requiring preservation of specified reasonable legal notices or
author attributions in that material or in the Appropriate Legal
Notices displayed by works containing it; or
c) Prohibiting misrepresentation of the origin of that material, or
requiring that modified versions of such material be marked in
reasonable ways as different from the original version; or
d) Limiting the use for publicity purposes of names of licensors or
authors of the material; or
e) Declining to grant rights under trademark law for use of some
trade names, trademarks, or service marks; or
f) Requiring indemnification of licensors and authors of that
material by anyone who conveys the material (or modified versions of
it) with contractual assumptions of liability to the recipient, for
any liability that these contractual assumptions directly impose on
those licensors and authors.
All other non-permissive additional terms are considered "further
restrictions" within the meaning of section 10. If the Program as you
received it, or any part of it, contains a notice stating that it is
governed by this License along with a term that is a further
restriction, you may remove that term. If a license document contains
a further restriction but permits relicensing or conveying under this
License, you may add to a covered work material governed by the terms
of that license document, provided that the further restriction does
not survive such relicensing or conveying.
If you add terms to a covered work in accord with this section, you
must place, in the relevant source files, a statement of the
additional terms that apply to those files, or a notice indicating
where to find the applicable terms.
Additional terms, permissive or non-permissive, may be stated in the
form of a separately written license, or stated as exceptions;
the above requirements apply either way.
8. Termination.
You may not propagate or modify a covered work except as expressly
provided under this License. Any attempt otherwise to propagate or
modify it is void, and will automatically terminate your rights under
this License (including any patent licenses granted under the third
paragraph of section 11).
However, if you cease all violation of this License, then your
license from a particular copyright holder is reinstated (a)
provisionally, unless and until the copyright holder explicitly and
finally terminates your license, and (b) permanently, if the copyright
holder fails to notify you of the violation by some reasonable means
prior to 60 days after the cessation.
Moreover, your license from a particular copyright holder is
reinstated permanently if the copyright holder notifies you of the
violation by some reasonable means, this is the first time you have
received notice of violation of this License (for any work) from that
copyright holder, and you cure the violation prior to 30 days after
your receipt of the notice.
Termination of your rights under this section does not terminate the
licenses of parties who have received copies or rights from you under
this License. If your rights have been terminated and not permanently
reinstated, you do not qualify to receive new licenses for the same
material under section 10.
9. Acceptance Not Required for Having Copies.
You are not required to accept this License in order to receive or
run a copy of the Program. Ancillary propagation of a covered work
occurring solely as a consequence of using peer-to-peer transmission
to receive a copy likewise does not require acceptance. However,
nothing other than this License grants you permission to propagate or
modify any covered work. These actions infringe copyright if you do
not accept this License. Therefore, by modifying or propagating a
covered work, you indicate your acceptance of this License to do so.
10. Automatic Licensing of Downstream Recipients.
Each time you convey a covered work, the recipient automatically
receives a license from the original licensors, to run, modify and
propagate that work, subject to this License. You are not responsible
for enforcing compliance by third parties with this License.
An "entity transaction" is a transaction transferring control of an
organization, or substantially all assets of one, or subdividing an
organization, or merging organizations. If propagation of a covered
work results from an entity transaction, each party to that
transaction who receives a copy of the work also receives whatever
licenses to the work the party's predecessor in interest had or could
give under the previous paragraph, plus a right to possession of the
Corresponding Source of the work from the predecessor in interest, if
the predecessor has it or can get it with reasonable efforts.
You may not impose any further restrictions on the exercise of the
rights granted or affirmed under this License. For example, you may
not impose a license fee, royalty, or other charge for exercise of
rights granted under this License, and you may not initiate litigation
(including a cross-claim or counterclaim in a lawsuit) alleging that
any patent claim is infringed by making, using, selling, offering for
sale, or importing the Program or any portion of it.
11. Patents.
A "contributor" is a copyright holder who authorizes use under this
License of the Program or a work on which the Program is based. The
work thus licensed is called the contributor's "contributor version".
A contributor's "essential patent claims" are all patent claims
owned or controlled by the contributor, whether already acquired or
hereafter acquired, that would be infringed by some manner, permitted
by this License, of making, using, or selling its contributor version,
but do not include claims that would be infringed only as a
consequence of further modification of the contributor version. For
purposes of this definition, "control" includes the right to grant
patent sublicenses in a manner consistent with the requirements of
this License.
Each contributor grants you a non-exclusive, worldwide, royalty-free
patent license under the contributor's essential patent claims, to
make, use, sell, offer for sale, import and otherwise run, modify and
propagate the contents of its contributor version.
In the following three paragraphs, a "patent license" is any express
agreement or commitment, however denominated, not to enforce a patent
(such as an express permission to practice a patent or covenant not to
sue for patent infringement). To "grant" such a patent license to a
party means to make such an agreement or commitment not to enforce a
patent against the party.
If you convey a covered work, knowingly relying on a patent license,
and the Corresponding Source of the work is not available for anyone
to copy, free of charge and under the terms of this License, through a
publicly available network server or other readily accessible means,
then you must either (1) cause the Corresponding Source to be so
available, or (2) arrange to deprive yourself of the benefit of the
patent license for this particular work, or (3) arrange, in a manner
consistent with the requirements of this License, to extend the patent
license to downstream recipients. "Knowingly relying" means you have
actual knowledge that, but for the patent license, your conveying the
covered work in a country, or your recipient's use of the covered work
in a country, would infringe one or more identifiable patents in that
country that you have reason to believe are valid.
If, pursuant to or in connection with a single transaction or
arrangement, you convey, or propagate by procuring conveyance of, a
covered work, and grant a patent license to some of the parties
receiving the covered work authorizing them to use, propagate, modify
or convey a specific copy of the covered work, then the patent license
you grant is automatically extended to all recipients of the covered
work and works based on it.
A patent license is "discriminatory" if it does not include within
the scope of its coverage, prohibits the exercise of, or is
conditioned on the non-exercise of one or more of the rights that are
specifically granted under this License. You may not convey a covered
work if you are a party to an arrangement with a third party that is
in the business of distributing software, under which you make payment
to the third party based on the extent of your activity of conveying
the work, and under which the third party grants, to any of the
parties who would receive the covered work from you, a discriminatory
patent license (a) in connection with copies of the covered work
conveyed by you (or copies made from those copies), or (b) primarily
for and in connection with specific products or compilations that
contain the covered work, unless you entered into that arrangement,
or that patent license was granted, prior to 28 March 2007.
Nothing in this License shall be construed as excluding or limiting
any implied license or other defenses to infringement that may
otherwise be available to you under applicable patent law.
12. No Surrender of Others' Freedom.
If conditions are imposed on you (whether by court order, agreement or
otherwise) that contradict the conditions of this License, they do not
excuse you from the conditions of this License. If you cannot convey a
covered work so as to satisfy simultaneously your obligations under this
License and any other pertinent obligations, then as a consequence you may
not convey it at all. For example, if you agree to terms that obligate you
to collect a royalty for further conveying from those to whom you convey
the Program, the only way you could satisfy both those terms and this
License would be to refrain entirely from conveying the Program.
13. Use with the GNU Affero General Public License.
Notwithstanding any other provision of this License, you have
permission to link or combine any covered work with a work licensed
under version 3 of the GNU Affero General Public License into a single
combined work, and to convey the resulting work. The terms of this
License will continue to apply to the part which is the covered work,
but the special requirements of the GNU Affero General Public License,
section 13, concerning interaction through a network will apply to the
combination as such.
14. Revised Versions of this License.
The Free Software Foundation may publish revised and/or new versions of
the GNU General Public License from time to time. Such new versions will
be similar in spirit to the present version, but may differ in detail to
address new problems or concerns.
Each version is given a distinguishing version number. If the
Program specifies that a certain numbered version of the GNU General
Public License "or any later version" applies to it, you have the
option of following the terms and conditions either of that numbered
version or of any later version published by the Free Software
Foundation. If the Program does not specify a version number of the
GNU General Public License, you may choose any version ever published
by the Free Software Foundation.
If the Program specifies that a proxy can decide which future
versions of the GNU General Public License can be used, that proxy's
public statement of acceptance of a version permanently authorizes you
to choose that version for the Program.
Later license versions may give you additional or different
permissions. However, no additional obligations are imposed on any
author or copyright holder as a result of your choosing to follow a
later version.
15. Disclaimer of Warranty.
THERE IS NO WARRANTY FOR THE PROGRAM, TO THE EXTENT PERMITTED BY
APPLICABLE LAW. EXCEPT WHEN OTHERWISE STATED IN WRITING THE COPYRIGHT
HOLDERS AND/OR OTHER PARTIES PROVIDE THE PROGRAM "AS IS" WITHOUT WARRANTY
OF ANY KIND, EITHER EXPRESSED OR IMPLIED, INCLUDING, BUT NOT LIMITED TO,
THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
PURPOSE. THE ENTIRE RISK AS TO THE QUALITY AND PERFORMANCE OF THE PROGRAM
IS WITH YOU. SHOULD THE PROGRAM PROVE DEFECTIVE, YOU ASSUME THE COST OF
ALL NECESSARY SERVICING, REPAIR OR CORRECTION.
16. Limitation of Liability.
IN NO EVENT UNLESS REQUIRED BY APPLICABLE LAW OR AGREED TO IN WRITING
WILL ANY COPYRIGHT HOLDER, OR ANY OTHER PARTY WHO MODIFIES AND/OR CONVEYS
THE PROGRAM AS PERMITTED ABOVE, BE LIABLE TO YOU FOR DAMAGES, INCLUDING ANY
GENERAL, SPECIAL, INCIDENTAL OR CONSEQUENTIAL DAMAGES ARISING OUT OF THE
USE OR INABILITY TO USE THE PROGRAM (INCLUDING BUT NOT LIMITED TO LOSS OF
DATA OR DATA BEING RENDERED INACCURATE OR LOSSES SUSTAINED BY YOU OR THIRD
PARTIES OR A FAILURE OF THE PROGRAM TO OPERATE WITH ANY OTHER PROGRAMS),
EVEN IF SUCH HOLDER OR OTHER PARTY HAS BEEN ADVISED OF THE POSSIBILITY OF
SUCH DAMAGES.
17. Interpretation of Sections 15 and 16.
If the disclaimer of warranty and limitation of liability provided
above cannot be given local legal effect according to their terms,
reviewing courts shall apply local law that most closely approximates
an absolute waiver of all civil liability in connection with the
Program, unless a warranty or assumption of liability accompanies a
copy of the Program in return for a fee.
END OF TERMS AND CONDITIONS
How to Apply These Terms to Your New Programs
If you develop a new program, and you want it to be of the greatest
possible use to the public, the best way to achieve this is to make it
free software which everyone can redistribute and change under these terms.
To do so, attach the following notices to the program. It is safest
to attach them to the start of each source file to most effectively
state the exclusion of warranty; and each file should have at least
the "copyright" line and a pointer to where the full notice is found.
<one line to give the program's name and a brief idea of what it does.>
Copyright (C) <year> <name of author>
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, see <https://www.gnu.org/licenses/>.
Also add information on how to contact you by electronic and paper mail.
If the program does terminal interaction, make it output a short
notice like this when it starts in an interactive mode:
<program> Copyright (C) <year> <name of author>
This program comes with ABSOLUTELY NO WARRANTY; for details type `show w'.
This is free software, and you are welcome to redistribute it
under certain conditions; type `show c' for details.
The hypothetical commands `show w' and `show c' should show the appropriate
parts of the General Public License. Of course, your program's commands
might be different; for a GUI interface, you would use an "about box".
You should also get your employer (if you work as a programmer) or school,
if any, to sign a "copyright disclaimer" for the program, if necessary.
For more information on this, and how to apply and follow the GNU GPL, see
<https://www.gnu.org/licenses/>.
The GNU General Public License does not permit incorporating your program
into proprietary programs. If your program is a subroutine library, you
may consider it more useful to permit linking proprietary applications with
the library. If this is what you want to do, use the GNU Lesser General
Public License instead of this License. But first, please read
<https://www.gnu.org/licenses/why-not-lgpl.html>.

49
Makefile Normal file
View File

@ -0,0 +1,49 @@
TARGET_EXEC := final_program.so
BUILD_DIR := ./build
SRC_DIRS := .
# Find all the C and C++ files we want to compile
# Note the single quotes around the * expressions. The shell will incorrectly expand these otherwise, but we want to send the * directly to the find command.
SRCS := $(shell find $(SRC_DIRS) -name '*.cpp' -or -name '*.c' -or -name '*.s')
# Prepends BUILD_DIR and appends .o to every src file
# As an example, ./your_dir/hello.cpp turns into ./build/./your_dir/hello.cpp.o
OBJS := $(SRCS:%=$(BUILD_DIR)/%.o)
# String substitution (suffix version without %).
# As an example, ./build/hello.cpp.o turns into ./build/hello.cpp.d
DEPS := $(OBJS:.o=.d)
# Every folder in ./src will need to be passed to GCC so that it can find header files
INC_DIRS := $(shell find $(SRC_DIRS) -path .git -prune -type d)
# Add a prefix to INC_DIRS. So moduleA would become -ImoduleA. GCC understands this -I flag
INC_FLAGS := $(addprefix -I,$(INC_DIRS))
# The -MMD and -MP flags together generate Makefiles for us!
# These files will have .d instead of .o as the output.
CPPFLAGS := $(INC_FLAGS) -MMD -MP -fPIC
# The final build step.
$(BUILD_DIR)/$(TARGET_EXEC): $(OBJS)
$(CXX) -shared $(OBJS) -o $@ $(LDFLAGS)
# Build step for C source
$(BUILD_DIR)/%.c.o: %.c
mkdir -p $(dir $@)
$(CC) $(CPPFLAGS) $(CFLAGS) -c $< -o $@
# Build step for C++ source
$(BUILD_DIR)/%.cpp.o: %.cpp
mkdir -p $(dir $@)
$(CXX) $(CPPFLAGS) $(CXXFLAGS) -c $< -o $@
.PHONY: clean
clean:
rm -r $(BUILD_DIR)
# Include the .d makefiles. The - at the front suppresses the errors of missing
# Makefiles. Initially, all the .d files will be missing, and we don't want those
# errors to show up.
-include $(DEPS)

2
README.md Normal file
View File

@ -0,0 +1,2 @@
# Mengubah
this repo simlifies the mengubah olugins build system.

23
common.h Normal file
View File

@ -0,0 +1,23 @@
#ifndef MENGA_DSP_COMMON
#define MENGA_DSP_COMMON
#include <complex>
typedef std::complex<float> Complex;
namespace Mengu {
// returns a Complex number with the same norm as x but with the angle phase
template <typename T>
std::complex<T> with_phase(const std::complex<T> &x, T phase);
// A complex number with the same phase as x, but with distance from the origin amp
template <typename T>
std::complex<T> with_amp(const std::complex<T> &x, T amp) {
return x / std::sqrt(std::norm(x)) * amp;
}
}
#endif

92
correlation.cpp Normal file
View File

@ -0,0 +1,92 @@
#include "correlation.h"
#include "common.h"
#include <iostream>
#include "mengumath.h"
#include <cstdint>
float Mengu::dsp::correlation(const Complex *s1, const Complex *s2, const int length, const int n) {
float total = 0;
for (uint32_t i = 0; i < length; i++) {
total += s1[i].real() * s2[n + i].real();
}
return total;
}
float Mengu::dsp::autocorrelation(const Complex *s, const int length, const int n) {
return correlation(s, s, length, n);
}
int Mengu::dsp::find_max_correlation(const Complex *s1, const Complex *s2, const int length, const int search_window_size) {
float max_corr = correlation(s1, s2, length, 0);
int max_lag = 0;
for (int i = 1; i < search_window_size; i++) {
float corr = correlation(s1, s2, length, i);
if (corr > max_corr) {
max_corr = corr;
max_lag = i;
}
}
return max_lag;
}
int Mengu::dsp::find_max_correlation_quad(const Complex *s1, const Complex *s2, const int length, const int search_window_size) {
float max_corr = -1e10;
int max_lag = 0;
float *scaled_s1 = new float[length];
for (int i = 0; i < length; i++) {
scaled_s1[i] = s1[i].real() * i * (length - i);
}
// do max iteration
for (int i = 0; i < search_window_size; i++) {
float corr = 0.0f;
for (int j = 0; j < length; j++) {
corr += scaled_s1[j] * s2[i + j].real();
}
if (corr > max_corr) {
max_corr = corr;
max_lag = i;
}
}
delete[] scaled_s1;
return max_lag;
}
std::vector<float> Mengu::dsp::calc_srhs(const float *envelope,
const int &size,
const int &min_freq_ind,
const int &max_freq_ind,
const int &n_harm,
const int &step) {
std::vector<float> output;
for (int freq_ind = min_freq_ind; freq_ind < max_freq_ind; freq_ind += step) {
output.push_back(calc_srh(envelope, size, freq_ind, n_harm));
}
return output;
}
float Mengu::dsp::calc_srh(const float *envelope, const int &size, const int &freq_ind, const int &n_harm) {
const int pos_n_harm = MIN(n_harm, size / freq_ind);
float pos_interference = 0;
for (int k = 1; k < pos_n_harm; k++) {
pos_interference += envelope[freq_ind * k];
}
const int neg_n_harm = MIN(n_harm, (int) ((float) size / freq_ind + 0.5));
float neg_interference = 0;
for (int k = 2; k < neg_n_harm; k++) {
neg_interference += envelope[(int) (freq_ind * k - 0.5)];
}
return pos_interference - neg_interference;
}

177
correlation.h Normal file
View File

@ -0,0 +1,177 @@
/**
* @file correlation.h
* @author 9exa
* @brief Correlation/Convolution functions between signals
* @version 0.1
* @date 2023-04-30
*
* @copyright Copyright (c) 2023
*/
#ifndef MENGA_CORRELATION
#define MENGA_CORRELATION
#include "common.h"
#include "fft.h"
#include "linalg.h"
#include "mengumath.h"
#include <algorithm>
#include <array>
#include <cmath>
#include <complex>
#include <cstdint>
#include <vector>
namespace Mengu {
namespace dsp {
// The (real finite non-circular) cross-correlation of two signals of equal length, on the offset n
// basically just a dot product
float correlation(const Complex *s1, const Complex *s2, const int length, const int n);
// The (real finite non-circular) cross-correlation of of a signal on itself, on the offset n
float autocorrelation(const Complex *s, const int length, const int n);
// Find the offset/lag that corresponds to the max cross-correlation between s1 and s2 explored up to length
// s2 is assumed to be at least length + search_window_size long
int find_max_correlation(const Complex *s1, const Complex *s2, const int length, const int search_window_size);
// Max correlation where portions toward the center are weighted more
int find_max_correlation_quad(const Complex *s1, const Complex *s2, const int length, const int search_window_size);
// find the sr harmonics of (the positive half of) a frequency amplitude spectrum
std::vector<float> calc_srhs(const float *envelope,
const int &size,
const int &min_freq_ind,
const int &max_freq_ind,
const int &n_harm = 8,
const int &step = 1);
// the srh of only one frequency
float calc_srh(const float *envelope, const int &size, const int &freq_ind, const int &n_harm);
// performs and stores results of LinearPredictiveCoding. expects a fixed process size so it can be put on the stack
template<uint32_t SampleSize, uint32_t NParams>
class LPC {
public:
LPC():
_fft(SampleSize),
// _b(NParams + 1),
// _autocovariance_slice(NParams + 1) {
_b{0},
_autocovariance_slice{0} {
_b[0] = 1;
}
// perform LPC on a sample and set up the intermediate variables
void load_sample(const Complex *sample) {
_fft.transform(sample, _freq_spectrum.data());
// assumes that sample are all real numbers; In general-purpose dsp this will cause bugs
// multiplication in the frequency domain is convolution (reversed correlation) in the real domain
std::array<Complex, SampleSize> freq_squared;
std::array<Complex, SampleSize> autocovariance_comp;
std::transform(
_freq_spectrum.cbegin(),
_freq_spectrum.cend(),
freq_squared.begin(),
[] (Complex f) { return std::norm(f); }
);
_fft.inverse_transform(freq_squared.data(), autocovariance_comp.data());
std::transform(
autocovariance_comp.cbegin(),
autocovariance_comp.cend(),
_autocovariance.begin(),
[] (Complex c) { return c.real(); }
);
std::copy(
_autocovariance.cbegin(),
_autocovariance.cbegin() + NParams + 1,
_autocovariance_slice.begin()
);
std::array<float, NParams + 1> a = solve_sym_toeplitz(_autocovariance_slice, _b);
std::array<Complex, SampleSize> a_comp{0};
const float a0 = a[0];
std::transform(a.cbegin(), a.cend(), a_comp.begin(),
[a0] (float f) { return Complex(f / a0); }
);
std::array<Complex, SampleSize> A;
_fft.transform(a_comp.data(), A.data());
// calc envelope
std::transform(A.cbegin(), A.cend(), _envelope.begin(),
// try to prevent infs
[] (Complex c) { return 1.0f / (sqrt(std::norm(c))); }
);
// calce residuals
for (uint32_t i = 0; i < SampleSize; i++) {
_residuals[i] = std::sqrt(std::norm(_freq_spectrum[i] * A[i]));
}
// std::transform(_freq_spectrum.cbegin(), _freq_spectrum.cend(), A.cbegin(), _residuals.begin(),
// [] (Complex x_val, Complex a_val) { return std::sqrt(std::norm(x_val * a_val)); }
// );
}
// The dft of the loaded samples
const std::array<Complex, SampleSize> &get_freq_spectrum() const {
return _freq_spectrum;
}
// The correlation of the signal with itself
const std::array<float, SampleSize> &get_autocovariance() const {
return _autocovariance;
}
// Normalized copy of autocovariance
std::array<float, SampleSize> get_autocorrelation() const {
std::array<float, SampleSize> autocorrelation;
float max_cov = *std::max_element(_autocovariance);
std::transform(
_autocovariance.cbegin(),
_autocovariance.cend(),
autocorrelation.begin(),
[max_cov] (float cov) {return cov / max_cov;}
);
}
// Envelope of the frequency spectrum
const std::array<float, SampleSize> &get_envelope() const {
return _envelope;
}
// LCP residuals of the frequency
const std::array<float, SampleSize> &get_residuals() const {
return _residuals;
}
// useful for inversion
const FFT &get_fft() const {
return _fft;
}
private:
// results to be getted
std::array<Complex, SampleSize> _freq_spectrum;
std::array<float, SampleSize> _autocovariance;
std::array<float, SampleSize> _envelope;
std::array<float, SampleSize> _residuals;
// intermediates
FFT _fft;
std::array<float, NParams + 1> _b;
std::array<float, NParams + 1> _autocovariance_slice;
// std::vector<float> _b;
// std::vector<float> _autocovariance_slice;
};
}
}
#endif

263
cyclequeue.h Normal file
View File

@ -0,0 +1,263 @@
/**
* @file cyclevector.h
* @author 9exa
* @brief An user-determined-size array where pushing an item removes one from the other end.
* Useful for not having to reallocate memory for rapidly appended contiguious data (i.e. sampling)
* @version 0.1
* @date 2023-02-21
*
* @copyright Copyright (c) 2023
*
*/
#ifndef MENGA_CYCLE_QUEUE
#define MENGA_CYCLE_QUEUE
#include <cstdint>
#include <vector>
#include <cstring>
#include <iostream>
#include "mengumath.h"
namespace Mengu {
template <class T>
class CycleQueue {
private:
T *_data = nullptr;
uint32_t _size = 0;
// what part of the data array is the front of the queue and first to be replaced on a push_back()
uint32_t _front = 0;
uint32_t _capacity = 0;
inline uint32_t start_inc_down1() {
// assert(size != 0)
// equiv to posmod (_front + _size - 1) % _size but faster?????
return _front == 0 ? _size - 1 : _front - 1;
}
public:
CycleQueue() {}
CycleQueue(uint32_t size) {
if (size > 0) {
resize(size);
}
}
// implement "copy" so they don't share the same data array
CycleQueue(const CycleQueue &from) {
_capacity = from._capacity;
_data = new T[_capacity];
std::cout << "from something\n";
memcpy(_data, from._data, _capacity * sizeof(T));
_size = from._size;
_front = from._front;
}
~CycleQueue() {
delete[] _data;
}
inline uint32_t size() const {
return _size;
}
inline void push_back(const T &x) {
if (_size == 0) return;
_data[_front] = x;
_front = (_front + 1) % _size;
}
inline void push_front(const T &x) {
if (_size == 0) return;
_front = (_front == 0) ? _size - 1 : _front - 1;
_data[_front] = x;
}
void resize(const uint32_t &new_size) {
if (_capacity < new_size) {
uint32_t new_cap = MAX(_capacity, 1);
while (new_cap < new_size) {
new_cap = new_cap << 1; // if you leave out the new_cap = the optimizer just skips this loop.
// Which makes this infinite loop bug hard to spot
}
reserve(new_cap);
} else if (new_size < _size) {
//move element of an array such that those at the front are removed
uint32_t shift;
uint32_t i;
if (_front <= new_size) {
// shift elements after _front down
shift = _size - new_size;
i = _front;
}
else {
// shift elements before _front up
shift = _front - new_size;
i = 0;
_front = 0;
}
for (; i < new_size; i++) {
_data[i] = _data[i + shift];
}
// set the deleted slots to be ready for future resizes
for (; i < _size; i++) {
_data[i] = T();
}
}
else { // (new_size > size) but reserve() has yet to initialize values
uint32_t i = _size;
for (; i < _front; i++) {
_data[(i + _size) % new_size] = _data[i];
_data[i] = T();
}
}
_size = new_size;
}
inline uint32_t capacity() const {
return _capacity;
}
const T *data() const {
return _data;
}
void reserve(const uint32_t &new_cap) {
if (_capacity < new_cap) {
T *new_data = new T[new_cap];
// copy and initialise new array
uint32_t i = 0;
if (_data != nullptr) {
for (; i < _capacity; i++) {
new_data[i] = std::move(_data[i]);
}
delete[] _data;
}
for (; i< new_cap; i++) {
new_data[i] = T();
}
_data = new_data;
_capacity = new_cap;
if (_size > _capacity) {
resize(_capacity);
}
}
}
// rotates the data array so that ir begins with _friont
void make_contiguous() {
T *new_data = new T[_capacity];
uint32_t i = 0;
if (_data != nullptr) {
for (; i < _size; i++) {
new_data[i] = std::move(_data[(i + _front) % _size]);
}
delete[] _data;
}
for (; i< _capacity; i++) {
new_data[i] = T();
}
_data = new_data;
}
void set(int i, const T &item) {
if (i > _size) {
std::cout << "tried to set item " << i << "on CycleQueue of size i. Out of range"<< std::endl;
}
_data[posmod(i +_front, _size)] = item;
}
const T get(int i) const {
if (i > _size) {
std::cout << "tried to get item " << i << "on CycleQueue of size i. Out of range"<< std::endl;
}
return _data[posmod(i + _front, _size)];
}
// just moves the front of the Ccle queue by an amount
void rotate(int by) {
_front = posmod(_front + by, _size);
}
//// Operators
inline T &operator[](int i) {
return _data[posmod(i + _front, _size)];
}
inline const T &operator[](int i) const {
// std::cout << (i + _front) % _size << std::endl;
return _data[posmod(i + _front, _size)];
}
CycleQueue &operator=(const CycleQueue &from) {
if (from._size != _size) {
resize(from._size);
}
_front = from._front;
memcpy(_data, from._data, _size);
return *this;
};
//// Conversions
// converts the first 'size' items into a contiguous array. -1 does the whole queue
T* to_array(T *out, int size = -1) const {
if (size == -1) {
size = _size;
}
for (uint32_t i = 0; i < size; i++) {
out[i] = _data[(i + _front) % _size];
}
return out;
}
// converts the first 'size' items into a vector. -1 does the whole queue
std::vector<T> to_vector(int size = -1) const {
if (size == -1) {
size = _size;
}
std::vector<T> out;
out.resize(size);
for (uint32_t i = 0; i < size; i++) {
out[i] = _data[(i + _front) % _size];
}
return out;
}
};
template<typename T>
std::string to_string(const CycleQueue<T> &cq) {
using std::to_string;
std::string out_string;
if (cq.size() == 0) {
return out_string;
}
out_string.reserve(to_string(cq[0]).size() * cq.size());
for (uint32_t i = 0; i < cq.size(); i++) {
out_string += to_string(cq[i]);
}
return out_string;
}
};
#endif

31
effect.cpp Normal file
View File

@ -0,0 +1,31 @@
#include "effect.h"
#include "mengumath.h"
Mengu::dsp::EffectChain::EffectChain(uint32_t buffer_size) {
_buffer_size = Mengu::last_pow_2(buffer_size);
_input_buffer.resize(_buffer_size);
_transformed_buffer.resize(_buffer_size);
}
Mengu::dsp::EffectChain::~EffectChain() {
for (auto effect: _effects) {
}
}
void Mengu::dsp::EffectChain::push_signal(const Complex *input, const uint32_t &size) {
for (uint32_t i = 0; i < size; i++) {
_input_buffer.push_back(input[i]);
}
}
void Mengu::dsp::EffectChain::pop_transformed_signal(Complex *output, const uint32_t &size) {
for (uint32_t i = _buffer_size - size; i < _buffer_size; i++) {
output[i] = _transformed_buffer[i];
}
}
void Mengu::dsp::EffectChain::append_effect(Effect *effect) {
_effects.push_back(effect);
}

113
effect.h Normal file
View File

@ -0,0 +1,113 @@
#ifndef MENGA_EFFECT
#define MENGA_EFFECT
#include <cstdint>
#include "common.h"
#include "cyclequeue.h"
#include <vector>
namespace Mengu {
namespace dsp {
// Types of Properties that an Effect has and how they can be edited by gui
// If this was Rust (a better language) this would all be one enum
enum EffectPropType {
Toggle, // Edited With a ToggleButton
Slider, // Edited with A Slider. A min and max value must be declared. Stores real number
Knob, // Edited with A Knob. A min and max value must be declared. Stores real number
Counter, // Edited with a counter. A step size must be declared. Stores real numbers, but often casted into an int
};
enum EffectPropContScale {
Linear,
Exp,
};
// Description for the editable properties of an Effect.
struct EffectPropDesc {
EffectPropType type;
const char *name;
const char *desc;
union {
struct {
float min_value;
float max_value;
float step_size;
EffectPropContScale scale;
} slider_data; // use by slider, knob and counter
};
};
// Data used to get and set EffectProperty data
struct EffectPropPayload {
EffectPropType type;
union {
bool on; // Used by Toggle
float value; // used by slider, knob, and counter
};
};
//
// object that takes in the next value signal (time or frequency domain)
// and can be queried for the next value in the process signal.
class Effect {
public:
virtual ~Effect() = default;
// tells an EffectChain what type of input the effect expects
enum InputDomain {
Time = 0,
Spectral = 1,
Frequency = 1
};
virtual InputDomain get_input_domain() = 0;
// push new value of signal
virtual void push_signal(const Complex *input, const uint32_t &size) = 0;
// Last value of transformed signal
virtual uint32_t pop_transformed_signal(Complex *output, const uint32_t &size) = 0;
// number of samples that can be output given the current pushed signals of the Effect
virtual uint32_t n_transformed_ready() const = 0;
// resets state of effect to make it reading to take in a new sample
virtual void reset() = 0;
// The properties that this Effect exposes to be changed by GUI.
// The index that they are put in is considered the props id
virtual std::vector<EffectPropDesc> get_property_descs() const = 0;
// Sets a property with the specified id the value declared in the payload
virtual void set_property(uint32_t id, EffectPropPayload data) = 0;
// Gets the value of a property with the specified id
virtual EffectPropPayload get_property(uint32_t id) const = 0;
};
// Represents a series of effects chained consequtivly. Processed on demand
class EffectChain {
public:
EffectChain(uint32_t buffer_size);
~EffectChain();
// push a new signal. Unlike an Effect the pushed signal can be an arbitrary size as it is stored in a ringbuffer
void push_signal(const Complex *input, const uint32_t &size);
// Last values of transformed signal
void pop_transformed_signal(Complex *output, const uint32_t &size);
// add an Effect
void append_effect(Effect *effect);
// apply all effects
void process();
private:
uint32_t _buffer_size;
CycleQueue<Complex> _input_buffer;
std::vector<Complex> _transformed_buffer;
std::vector<Effect *> _effects;
};
}
}
#endif

41
fastmath.h Normal file
View File

@ -0,0 +1,41 @@
#ifndef MENGA_FAST_MATH
#define MENGA_FAST_MATH
#include "mengumath.h"
#include <cmath>
namespace Mengu {
namespace dsp {
// PADE approximations. Only accurate on certain ranges
// use between -pi to +pi
template<typename FloatType>
FloatType sin(FloatType x) {
x = std::fmod(x + MATH_PI, MATH_TAU) - MATH_PI;
FloatType x2 = x * x;
FloatType numerator = x * ((11 * x2) * x2 + 2520);
FloatType denominator = 60 * x2 + 2520;
return numerator / denominator;
}
// use between -pi to +pi
template<typename FloatType>
FloatType cos(FloatType x) {
x = std::fmod(x + MATH_PI, MATH_TAU) - MATH_PI;
FloatType x2 = x * x;
FloatType numerator = 131040 + x2 * (62160 + x2 * (3814 - x2 * 59));
FloatType denominator = 131040 + x2 * (3360 + x2 * 34);
return numerator / denominator;
}
template<typename FloatType>
FloatType tan(FloatType x) {
x = std::fmod(x + MATH_PI, MATH_TAU) - MATH_PI;
FloatType x2 = x * x;
FloatType numerator = x * (-135135 + x2 * (17325 + x2 * (-378 + x2)));
FloatType denominator = -135135 + x2 * (62370 + x2 * (-3150 + 28 * x2));
return numerator / denominator;
}
} // namespace dsp
} // namespace Mengu
#endif

155
fft.cpp Normal file
View File

@ -0,0 +1,155 @@
#include "fft.h"
#include "common.h"
#include "fastmath.h"
#include "mengumath.h"
#include <cstdint>
#include <vector>
using namespace Mengu;
// helpers
static void conjugate_arr(Complex *arr, const uint32_t a_size) {
for (uint32_t i = 0; i < a_size; i++) {
arr[i].imag(-arr[i].imag());
}
}
std::vector<Complex> dsp::SFT::perform(const std::vector<Complex> &input) {
const size_t outsize = input.size();
const Complex mitau(0, -MATH_TAU);
Complex *es = new Complex[outsize];
for (size_t j = 0; j < outsize; j++) {
float n = (float) j / outsize;
es[j] = std::exp(mitau * n);
// std::cout << es[j] << std::endl;
}
std::vector<Complex> outvec;
outvec.resize(outsize);
for (size_t k = 0; k < outsize; k++) {
for (size_t n = 0; n < outsize; n++) {
outvec[k] += input[n] * es[n * k % outsize];
}
outvec[k] /= outsize;
}
delete[] es;
return outvec;
}
dsp::FFT::FFT(uint32_t size) {
_fft_size = is_pow_2(size) ? size : next_pow_2(size);
_size = size;
_es = new complex<float>[_fft_size];
for (uint32_t i = 0; i < _fft_size; i++) {
_es[i] = std::polar(1.0f, -(float) MATH_TAU * i / _fft_size);
}
_inp_vec = new Complex[_fft_size];
_out_vec = new Complex[_fft_size];
}
dsp::FFT::~FFT() {
delete[] _es;
delete[] _inp_vec;
delete[] _out_vec;
}
void dsp::FFT::transform(const Complex *input, Complex *output) const {
for (uint32_t i = 0; i < _size; i++) {
_inp_vec[i] = input[i];
}
for (uint32_t i = _size; i < _fft_size; i++) {
_inp_vec[i] = Complex(0.0f);
}
_transform_rec(_inp_vec, _out_vec, _fft_size, 1);
for (uint32_t i = 0; i < _size; i++) {
output[i] = _out_vec[i] / sqrtf((float) _fft_size);
}
}
void dsp::FFT::transform(const CycleQueue<Complex>&input, Complex *output) const {
// _transform_rec(input, 0, output, _size, 1);
input.to_array(_inp_vec);
// zero pad the input
for (uint32_t i = input.size(); i < _fft_size; i++) {
_inp_vec[i] = 0.0f;
}
_transform_rec(_inp_vec, _out_vec, _fft_size, 1);
for (uint32_t i = 0; i < _size / 2; i++) {
output[i] = _out_vec[i] / sqrtf((float) _fft_size);
}
}
void dsp::FFT::inverse_transform(const Complex *input, Complex *output) const {
for (uint32_t i = 0; i < _size; i++) {
_inp_vec[i] = input[i];
}
// zero pad input buffer
for (uint32_t i = _size; i < _fft_size; i++) {
_inp_vec[i] = 0.0f;
}
conjugate_arr(_inp_vec, _size);
_transform_rec(_inp_vec, _out_vec, _fft_size, 1);
conjugate_arr(_out_vec, _size);
for (uint32_t i = 0; i < _size; i++) {
// only inverse-transforming first half of the freq spectrum. so multiply by 2
output[i] = _out_vec[i] / sqrtf((float) _fft_size);
}
}
void dsp::FFT::_transform_rec(const Complex *input, Complex *output, const uint32_t N, const uint32_t stride) const {
// recursive implementaion of radix-2 fft
// Thanks wikipedia
if (N == 1) {
output[0] = input[0];
return;
}
_transform_rec(input, output, N / 2, 2 * stride);
_transform_rec(input + stride, output + N / 2, N / 2, 2 * stride);
for (uint32_t k = 0; k < N / 2; k++) {
const Complex e = _es[k * stride];
const Complex p = output[k]; // kth even
const Complex q = e * output[k + N / 2]; //kth odd
output[k] = p + q;
output[k + N / 2] = p - q;
}
}
void dsp::FFT::_transform_rec(const CycleQueue<Complex> &input, uint32_t inp_ind, Complex *output, const uint32_t N, const uint32_t stride) const{
// recursive implementaion of radix-2 fft
// Thanks wikipedia
if (N == 1) {
output[0] = input[inp_ind];
return;
}
_transform_rec(input, inp_ind, output, N / 2, 2 * stride);
_transform_rec(input, inp_ind + stride, output + N / 2, N / 2, 2 * stride);
for (uint32_t k = 0; k < N / 2; k++) {
const Complex e = _es[k * stride];
const Complex p = output[k]; // kth even
const Complex q = e * output[k + N / 2]; //kth odd
output[k] = p + q;
output[k + N / 2] = p - q;
}
}
const Complex *dsp::FFT::get_es() const {
return _es;
}

83
fft.h Normal file
View File

@ -0,0 +1,83 @@
#ifndef MENGA_FFT
#define MENGA_FFT
#include <valarray>
#include <vector>
#include "common.h"
#include "cyclequeue.h"
typedef std::complex<float> Complex;
typedef std::valarray<Complex> CArray;
namespace Mengu {
namespace dsp {
using namespace std;
class SFT {
// Slow O(N^2) Fourier Transform
public:
// SFT();
// out-of-place ft
vector<Complex> perform(const vector<Complex> &input);
private:
vector<Complex> _es;
// vector<Complex> _sines;
};
class FFTBuffer; // foward declaration
class FFT {
public:
// Fast Fourier transform that uses lookup tables and performs onto an established array
// Only works for arrays of a declared size. Designed to be cached and use many times
// To use it for different sample rates/lengths, create a new FFT of a different size
// The size needs to be a power of 2 in order for this to work properly
FFT(uint32_t size);
~FFT();
// Both arrays must be at least as long as _size
void transform(const Complex *input, Complex *output) const;
void transform(const CycleQueue<Complex> &input, Complex *output) const; // to avoid having to reallocate array
//unlike 'transform' this edits *input to avoid unnecissary memory allocation
void inverse_transform(const Complex *input, Complex *output) const;
const Complex *get_es() const;
uint32_t size() const { return _size; }
private:
Complex *_es;
uint32_t _size;
uint32_t _fft_size; // size of arrays used in fft computations. must be a power of 2
// cache buffers used in intermediate calculation
Complex *_inp_vec;
Complex *_out_vec;
void _transform_rec(const Complex *input, Complex *output, const uint32_t N, const uint32_t stride) const;
void _transform_rec(const CycleQueue<Complex> &input, uint32_t inp_ind, Complex *output, const uint32_t N, const uint32_t stride) const;
friend class FFTBuffer;
public:
};
class FFTBuffer {
// class designed to do update a n FFT of a signal in real time
public:
void push_signal(const Complex *x, const uint32_t &size) {}
// copies the last 'size'
void pop_transformed_signal(const Complex *output, const uint32_t &size) {}
private:
CycleQueue<Complex> _buffer;
};
};
};
#endif

31
filter.cpp Normal file
View File

@ -0,0 +1,31 @@
#include "filter.h"
#include "common.h"
#include <complex>
#include <cstdint>
using namespace Mengu;
using namespace dsp;
Complex Mengu::dsp::quad_filter_trans(Complex z, float a1, float a2, float b0, float b1, float b2) {
Complex z1 = std::conj(z); //z^-1
Complex z2 = z1 * z1;
return (b0 + b1 * z1 + b2 * z2) / (1.0f + a1 * z1 + a2 * z2);
}
BiquadFilter::BiquadFilter(float pa1, float pa2, float pb0, float pb1, float pb2) :
a1(pa1), a2(pa2), b0(pb0), b1(pb1), b2(pb2) {}
void BiquadFilter::transform(const float *input, float *output, uint32_t size) {
for (uint32_t i = 0; i < size; i++) {
float m = input[i] - a1 * _last_ms[_last_offset] - a2 * _last_ms[(_last_offset + 1) % 2];
output[i] = b0 * m + b1 * _last_ms[_last_offset] + b2 * _last_ms[(_last_offset + 1) % 2];
// push the last m
_last_offset = (_last_offset + 1) % 2;
_last_ms[_last_offset] = m;
}
}
void BiquadFilter::reset() {
_last_ms[0] = _last_ms[1] = 0.0f;
}

44
filter.h Normal file
View File

@ -0,0 +1,44 @@
/**
* @file filter.h
* @author your name (you@domain.com)
* @brief Filter functions
* @version 0.1
* @date 2023-06-10
*
* @copyright Copyright (c) 2023
*
*/
#ifndef MENGU_FILTER
#define MENGU_FILTER
#include "common.h"
#include <cstdint>
namespace Mengu {
namespace dsp {
// the transfer function for a quadratic filter with denominator coefficients a1, a2 and numerator cofficients b0, b1, b2
// assumes z is on a unit circle
Complex quad_filter_trans(Complex z, float a1, float a2, float b0, float b1, float b2);
class BiquadFilter {
public:
float a1;
float a2;
float b0;
float b1;
float b2;
BiquadFilter(float pa1, float pa2, float pb0, float pb1, float pb2);
void transform(const float *input, float *output, uint32_t size);
void reset();
private:
// store the last 2 raw inputs and intermediaries
float _last_ms[2] {0};
uint32_t _last_offset = 0;
};
}
}
#endif

147
formantshifter.cpp Normal file
View File

@ -0,0 +1,147 @@
#include "formantshifter.h"
#include "common.h"
#include "effect.h"
#include "interpolation.h"
#include "loudness.h"
#include "mengumath.h"
#include <algorithm>
#include <array>
#include <cmath>
#include <complex>
#include <cstdint>
using namespace Mengu;
using namespace dsp;
LPCFormantShifter::LPCFormantShifter() {
_transformed_buffer.resize(OverlapSize);
}
Effect::InputDomain LPCFormantShifter::get_input_domain() {
return InputDomain::Time;
};
// push new value of signal
void LPCFormantShifter::push_signal(const Complex *input, const uint32_t &size) {
_raw_buffer.extend_back(input, size);
}
// Last value of transformed signal
uint32_t LPCFormantShifter::pop_transformed_signal(Complex *output, const uint32_t &size) {
std::array<Complex, ProcSize> freq_shifted {0};
std::array<Complex, ProcSize> samples {0};
std::array<Complex, ProcSize> shifted_samples {0};
while (_raw_buffer.size() >= ProcSize && _transformed_buffer.size() < size + OverlapSize) {
// Load sample segment. 0 unused frames
_raw_buffer.to_array(samples.data(), ProcSize);
// do the shifty
_lpc.load_sample(samples.data());
_shift_by_env(
_lpc.get_freq_spectrum().data(),
freq_shifted.data(),
_lpc.get_envelope().data(),
_shift_factor
);
_lpc.get_fft().inverse_transform(freq_shifted.data(), shifted_samples.data());
// Make downward shifts not quieter and upward shifts not louder
// Automatically adjusts for the fact that only half of the frequency spectrum is used
_loudness_norm.normalize(shifted_samples.data(), samples.data(), shifted_samples.data());
// copy to output
mix_and_extend(_transformed_buffer, shifted_samples, OverlapSize, hamming_window);
_raw_buffer.pop_front_many(nullptr, HopSize);
}
if (_transformed_buffer.size() < size + OverlapSize) {
// Not enough samples to output anything
return 0;
}
else {
return _transformed_buffer.pop_front_many(output, size);
}
}
// number of samples that can be output given the current pushed signals of the Effect
uint32_t LPCFormantShifter::n_transformed_ready() const {
return _raw_buffer.size();
};
// resets state of effect to make it reading to take in a new sample
void LPCFormantShifter::reset() {
_raw_sample_filter.reset();
_shifted_sample_filter.reset();
}
// The properties that this Effect exposes to be changed by GUI.
// The index that they are put in is considered the props id
std::vector<EffectPropDesc> LPCFormantShifter::get_property_descs() const {
return {
EffectPropDesc {
.type = EffectPropType::Slider,
.name = "Formant Shift",
.desc = "Scales the formant of pushed signals by this amount",
.slider_data = {
.min_value = 0.5,
.max_value = 2,
.scale = Exp,
}
}
};
}
// Sets a property with the specified id the value declared in the payload
void LPCFormantShifter::set_property(uint32_t id, EffectPropPayload data) {
switch (id) {
case 0:
if (data.type == Slider) {
_shift_factor = data.value;
}
break;
default:
break;
}
}
// Gets the value of a property with the specified id
EffectPropPayload LPCFormantShifter::get_property(uint32_t id) const {
return EffectPropPayload {
.type = Slider,
.value = _shift_factor,
};
};
// rescale an array in the freqency domain by the shape of an envelope if it were to be shifted up or down
void LPCFormantShifter::_shift_by_env(const Complex *input,
Complex *output,
const float *envelope,
const float shift_factor) {
for (uint32_t i = 0; i < ProcSize / 2; i++) {
uint32_t shifted_ind = i / shift_factor;
if (shifted_ind < ProcSize / 2) {
float correction = envelope[shifted_ind] / envelope[i];
if (!std::isfinite(correction)) {
output[i] = Complex(0.0f);
}
else {
output[i] = correction * input[i];
}
}
else {
// output[i] = input[size - 1];
output[i] = Complex(0.0f);
}
}
}

75
formantshifter.h Normal file
View File

@ -0,0 +1,75 @@
/**
* @file formantshifter.h
* @author 9exa
* @brief Shifts the formants of a signal in real time.
* @date 2023-06-08
*/
#ifndef MENGU_FORMANT_SHIFTER
#define MENGU_FORMANT_SHIFTER
#include "common.h"
#include "correlation.h"
#include "effect.h"
#include "loudness.h"
#include "vecdeque.h"
#include <cstdint>
namespace Mengu {
namespace dsp {
// shifts formants using lpc envelope estimation
class LPCFormantShifter: public Effect {
public:
LPCFormantShifter();
// tells an EffectChain what type of input the effect expects
virtual InputDomain get_input_domain() override;
// push new value of signal
virtual void push_signal(const Complex *input, const uint32_t &size) override;
// Last value of transformed signal
virtual uint32_t pop_transformed_signal(Complex *output, const uint32_t &size) override;
// number of samples that can be output given the current pushed signals of the Effect
virtual uint32_t n_transformed_ready() const override;
// resets state of effect to make it reading to take in a new sample
virtual void reset() override;
// The properties that this Effect exposes to be changed by GUI.
// The index that they are put in is considered the props id
virtual std::vector<EffectPropDesc> get_property_descs() const override;
// Sets a property with the specified id the value declared in the payload
virtual void set_property(uint32_t id, EffectPropPayload data) override;
// Gets the value of a property with the specified id
virtual EffectPropPayload get_property(uint32_t id) const override;
private:
VecDeque<Complex> _raw_buffer;
VecDeque<Complex> _transformed_buffer;
static constexpr uint32_t ProcSize = 1 << 11;
static constexpr uint32_t HopSize = ProcSize * 4 / 5;
static constexpr uint32_t OverlapSize = ProcSize - HopSize;
LPC<ProcSize, 60> _lpc;
void _shift_by_env(const Complex *input,
Complex *output,
const float *envelope,
const float shift_factor);
float _shift_factor = 1.0f;
// Amplifies the formant_shifted samples so they have the same LUFS loudness as the raw_sample
LoudnessNormalizer<Complex, ProcSize, 2> _loudness_norm;
LUFSFilter _raw_sample_filter;
LUFSFilter _shifted_sample_filter;
};
}
}
#endif

96
interpolation.h Normal file
View File

@ -0,0 +1,96 @@
#ifndef MENGA_INTERPOLATION
#define MENGA_INTERPOLATION
#include "mengumath.h"
#include "fastmath.h"
#include <cmath>
#include <cstdint>
// Helper functions for various interpolation and windowing methods
namespace Mengu {
namespace dsp{
// The hann function
inline float hann(float a0, float x) {
return a0 - (1.0 - a0) * std::cos(MATH_TAU * x);
// equals
// a - (1-a) (c2 -s2)
// a (1 - c2 + s2) - (c2 - s2)
// a (2s2) - 2s2 - 1
// (a - 1) 2 s2
}
// root of the hann function designed to be used once each before and after processing
inline float hann_root(float a0, float x) {
return std::sqrt(2.0 *(a0 - 1)) * std::sin(MATH_PI * x);
}
// the hann function, centered around 0
// inline float hann_centered(float a0, float x) {
// return hann(a0, (x + 0.5f));
// }
// windows are 1 at 1, 0 at 0
inline float hann_window(float x) {
return hann(0.5, 0.5f * x);
}
inline float hamming_window(float x) {
return hann((float)25 / 46, 0.5f * x);
}
// windows are 1 at 1, 0 at 0
inline float hann_window_root(float x) {
return std::sin(0.5 * MATH_PI * x);
}
// does windowing to both ends of a sequence, in place
template<class T>
inline void window_ends(T* a, uint32_t a_size, uint32_t window_size, float window_f (float)) {
window_size = MIN(a_size / 2, window_size);
for (uint32_t i = 0; i < window_size; i++) {
const float w = window_f((float) i / window_size);
a[i] *= w;
a[a_size - 1 - i] *= w;
}
}
// added the tail and head of an array after applying a window function to bothe sides
template<class T>
inline void overlap_add(const T *prev, const T *next, T *output, uint32_t size, float window_f (float)) {
for (uint32_t i = 0; i < size; i++) {
const float w = window_f((float) i / size);
output[i] = w * next[i] + (1- w) * prev[i];
}
}
// overlap and extend after applying a window function first
// overlap_size may be larger than new_data.size(), in which case only new_data.size() is added. the offset is the same
template<class T, class NewT>
inline void mix_and_extend(T &array, const NewT &new_data, const uint32_t &overlap_size, float window_f (float)) {
uint32_t i = 0;
for(; i < MIN(overlap_size, new_data.size()); i++) {
float w = (float) i / overlap_size;
float w1 = window_f(w);
float w2 = window_f(1.0-w);
w2 = 1.0f - w1;
uint32_t ind = array.size() - overlap_size + i;
array[ind] = array[ind] * w2 + new_data[i] * w1;
}
for (; i < new_data.size(); i++) {
array.push_back(new_data[i]);
}
}
}
}
#endif

603
ladspa.h Normal file
View File

@ -0,0 +1,603 @@
/* ladspa.h
Linux Audio Developer's Simple Plugin API Version 1.1[LGPL].
Copyright (C) 2000-2002 Richard W.E. Furse, Paul Barton-Davis,
Stefan Westerfeld.
This library is free software; you can redistribute it and/or
modify it under the terms of the GNU Lesser General Public License
as published by the Free Software Foundation; either version 2.1 of
the License, or (at your option) any later version.
This library 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
Lesser General Public License for more details.
You should have received a copy of the GNU Lesser General Public
License along with this library; if not, write to the Free Software
Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
USA. */
#ifndef LADSPA_INCLUDED
#define LADSPA_INCLUDED
#define LADSPA_VERSION "1.1"
#define LADSPA_VERSION_MAJOR 1
#define LADSPA_VERSION_MINOR 1
#ifdef __cplusplus
extern "C" {
#endif
/*****************************************************************************/
/* Overview:
There is a large number of synthesis packages in use or development
on the Linux platform at this time. This API (`The Linux Audio
Developer's Simple Plugin API') attempts to give programmers the
ability to write simple `plugin' audio processors in C/C++ and link
them dynamically (`plug') into a range of these packages (`hosts').
It should be possible for any host and any plugin to communicate
completely through this interface.
This API is deliberately short and simple. To achieve compatibility
with a range of promising Linux sound synthesis packages it
attempts to find the `greatest common divisor' in their logical
behaviour. Having said this, certain limiting decisions are
implicit, notably the use of a fixed type (LADSPA_Data) for all
data transfer and absence of a parameterised `initialisation'
phase. See below for the LADSPA_Data typedef.
Plugins are expected to distinguish between control and audio
data. Plugins have `ports' that are inputs or outputs for audio or
control data and each plugin is `run' for a `block' corresponding
to a short time interval measured in samples. Audio data is
communicated using arrays of LADSPA_Data, allowing a block of audio
to be processed by the plugin in a single pass. Control data is
communicated using single LADSPA_Data values. Control data has a
single value at the start of a call to the `run()' or `run_adding()'
function, and may be considered to remain this value for its
duration. The plugin may assume that all its input and output ports
have been connected to the relevant data location (see the
`connect_port()' function below) before it is asked to run.
Plugins will reside in shared object files suitable for dynamic
linking by dlopen() and family. The file will provide a number of
`plugin types' that can be used to instantiate actual plugins
(sometimes known as `plugin instances') that can be connected
together to perform tasks.
This API contains very limited error-handling. */
/*****************************************************************************/
/* Fundamental data type passed in and out of plugin. This data type
is used to communicate audio samples and control values. It is
assumed that the plugin will work sensibly given any numeric input
value although it may have a preferred range (see hints below).
For audio it is generally assumed that 1.0f is the `0dB' reference
amplitude and is a `normal' signal level. */
typedef float LADSPA_Data;
/*****************************************************************************/
/* Special Plugin Properties:
Optional features of the plugin type are encapsulated in the
LADSPA_Properties type. This is assembled by ORing individual
properties together. */
typedef int LADSPA_Properties;
/* Property LADSPA_PROPERTY_REALTIME indicates that the plugin has a
real-time dependency (e.g. listens to a MIDI device) and so its
output must not be cached or subject to significant latency. */
#define LADSPA_PROPERTY_REALTIME 0x1
/* Property LADSPA_PROPERTY_INPLACE_BROKEN indicates that the plugin
may cease to work correctly if the host elects to use the same data
location for both input and output (see connect_port()). This
should be avoided as enabling this flag makes it impossible for
hosts to use the plugin to process audio `in-place.' */
#define LADSPA_PROPERTY_INPLACE_BROKEN 0x2
/* Property LADSPA_PROPERTY_HARD_RT_CAPABLE indicates that the plugin
is capable of running not only in a conventional host but also in a
`hard real-time' environment. To qualify for this the plugin must
satisfy all of the following:
(1) The plugin must not use malloc(), free() or other heap memory
management within its run() or run_adding() functions. All new
memory used in run() must be managed via the stack. These
restrictions only apply to the run() function.
(2) The plugin will not attempt to make use of any library
functions with the exceptions of functions in the ANSI standard C
and C maths libraries, which the host is expected to provide.
(3) The plugin will not access files, devices, pipes, sockets, IPC
or any other mechanism that might result in process or thread
blocking.
(4) The plugin will take an amount of time to execute a run() or
run_adding() call approximately of form (A+B*SampleCount) where A
and B depend on the machine and host in use. This amount of time
may not depend on input signals or plugin state. The host is left
the responsibility to perform timings to estimate upper bounds for
A and B. */
#define LADSPA_PROPERTY_HARD_RT_CAPABLE 0x4
#define LADSPA_IS_REALTIME(x) ((x) & LADSPA_PROPERTY_REALTIME)
#define LADSPA_IS_INPLACE_BROKEN(x) ((x) & LADSPA_PROPERTY_INPLACE_BROKEN)
#define LADSPA_IS_HARD_RT_CAPABLE(x) ((x) & LADSPA_PROPERTY_HARD_RT_CAPABLE)
/*****************************************************************************/
/* Plugin Ports:
Plugins have `ports' that are inputs or outputs for audio or
data. Ports can communicate arrays of LADSPA_Data (for audio
inputs/outputs) or single LADSPA_Data values (for control
input/outputs). This information is encapsulated in the
LADSPA_PortDescriptor type which is assembled by ORing individual
properties together.
Note that a port must be an input or an output port but not both
and that a port must be a control or audio port but not both. */
typedef int LADSPA_PortDescriptor;
/* Property LADSPA_PORT_INPUT indicates that the port is an input. */
#define LADSPA_PORT_INPUT 0x1
/* Property LADSPA_PORT_OUTPUT indicates that the port is an output. */
#define LADSPA_PORT_OUTPUT 0x2
/* Property LADSPA_PORT_CONTROL indicates that the port is a control
port. */
#define LADSPA_PORT_CONTROL 0x4
/* Property LADSPA_PORT_AUDIO indicates that the port is a audio
port. */
#define LADSPA_PORT_AUDIO 0x8
#define LADSPA_IS_PORT_INPUT(x) ((x) & LADSPA_PORT_INPUT)
#define LADSPA_IS_PORT_OUTPUT(x) ((x) & LADSPA_PORT_OUTPUT)
#define LADSPA_IS_PORT_CONTROL(x) ((x) & LADSPA_PORT_CONTROL)
#define LADSPA_IS_PORT_AUDIO(x) ((x) & LADSPA_PORT_AUDIO)
/*****************************************************************************/
/* Plugin Port Range Hints:
The host may wish to provide a representation of data entering or
leaving a plugin (e.g. to generate a GUI automatically). To make
this more meaningful, the plugin should provide `hints' to the host
describing the usual values taken by the data.
Note that these are only hints. The host may ignore them and the
plugin must not assume that data supplied to it is meaningful. If
the plugin receives invalid input data it is expected to continue
to run without failure and, where possible, produce a sensible
output (e.g. a high-pass filter given a negative cutoff frequency
might switch to an all-pass mode).
Hints are meaningful for all input and output ports but hints for
input control ports are expected to be particularly useful.
More hint information is encapsulated in the
LADSPA_PortRangeHintDescriptor type which is assembled by ORing
individual hint types together. Hints may require further
LowerBound and UpperBound information.
All the hint information for a particular port is aggregated in the
LADSPA_PortRangeHint structure. */
typedef int LADSPA_PortRangeHintDescriptor;
/* Hint LADSPA_HINT_BOUNDED_BELOW indicates that the LowerBound field
of the LADSPA_PortRangeHint should be considered meaningful. The
value in this field should be considered the (inclusive) lower
bound of the valid range. If LADSPA_HINT_SAMPLE_RATE is also
specified then the value of LowerBound should be multiplied by the
sample rate. */
#define LADSPA_HINT_BOUNDED_BELOW 0x1
/* Hint LADSPA_HINT_BOUNDED_ABOVE indicates that the UpperBound field
of the LADSPA_PortRangeHint should be considered meaningful. The
value in this field should be considered the (inclusive) upper
bound of the valid range. If LADSPA_HINT_SAMPLE_RATE is also
specified then the value of UpperBound should be multiplied by the
sample rate. */
#define LADSPA_HINT_BOUNDED_ABOVE 0x2
/* Hint LADSPA_HINT_TOGGLED indicates that the data item should be
considered a Boolean toggle. Data less than or equal to zero should
be considered `off' or `false,' and data above zero should be
considered `on' or `true.' LADSPA_HINT_TOGGLED may not be used in
conjunction with any other hint except LADSPA_HINT_DEFAULT_0 or
LADSPA_HINT_DEFAULT_1. */
#define LADSPA_HINT_TOGGLED 0x4
/* Hint LADSPA_HINT_SAMPLE_RATE indicates that any bounds specified
should be interpreted as multiples of the sample rate. For
instance, a frequency range from 0Hz to the Nyquist frequency (half
the sample rate) could be requested by this hint in conjunction
with LowerBound = 0 and UpperBound = 0.5. Hosts that support bounds
at all must support this hint to retain meaning. */
#define LADSPA_HINT_SAMPLE_RATE 0x8
/* Hint LADSPA_HINT_LOGARITHMIC indicates that it is likely that the
user will find it more intuitive to view values using a logarithmic
scale. This is particularly useful for frequencies and gains. */
#define LADSPA_HINT_LOGARITHMIC 0x10
/* Hint LADSPA_HINT_INTEGER indicates that a user interface would
probably wish to provide a stepped control taking only integer
values. Any bounds set should be slightly wider than the actual
integer range required to avoid floating point rounding errors. For
instance, the integer set {0,1,2,3} might be described as [-0.1,
3.1]. */
#define LADSPA_HINT_INTEGER 0x20
/* The various LADSPA_HINT_HAS_DEFAULT_* hints indicate a `normal'
value for the port that is sensible as a default. For instance,
this value is suitable for use as an initial value in a user
interface or as a value the host might assign to a control port
when the user has not provided one. Defaults are encoded using a
mask so only one default may be specified for a port. Some of the
hints make use of lower and upper bounds, in which case the
relevant bound or bounds must be available and
LADSPA_HINT_SAMPLE_RATE must be applied as usual. The resulting
default must be rounded if LADSPA_HINT_INTEGER is present. Default
values were introduced in LADSPA v1.1. */
#define LADSPA_HINT_DEFAULT_MASK 0x3C0
/* This default values indicates that no default is provided. */
#define LADSPA_HINT_DEFAULT_NONE 0x0
/* This default hint indicates that the suggested lower bound for the
port should be used. */
#define LADSPA_HINT_DEFAULT_MINIMUM 0x40
/* This default hint indicates that a low value between the suggested
lower and upper bounds should be chosen. For ports with
LADSPA_HINT_LOGARITHMIC, this should be exp(log(lower) * 0.75 +
log(upper) * 0.25). Otherwise, this should be (lower * 0.75 + upper
* 0.25). */
#define LADSPA_HINT_DEFAULT_LOW 0x80
/* This default hint indicates that a middle value between the
suggested lower and upper bounds should be chosen. For ports with
LADSPA_HINT_LOGARITHMIC, this should be exp(log(lower) * 0.5 +
log(upper) * 0.5). Otherwise, this should be (lower * 0.5 + upper *
0.5). */
#define LADSPA_HINT_DEFAULT_MIDDLE 0xC0
/* This default hint indicates that a high value between the suggested
lower and upper bounds should be chosen. For ports with
LADSPA_HINT_LOGARITHMIC, this should be exp(log(lower) * 0.25 +
log(upper) * 0.75). Otherwise, this should be (lower * 0.25 + upper
* 0.75). */
#define LADSPA_HINT_DEFAULT_HIGH 0x100
/* This default hint indicates that the suggested upper bound for the
port should be used. */
#define LADSPA_HINT_DEFAULT_MAXIMUM 0x140
/* This default hint indicates that the number 0 should be used. Note
that this default may be used in conjunction with
LADSPA_HINT_TOGGLED. */
#define LADSPA_HINT_DEFAULT_0 0x200
/* This default hint indicates that the number 1 should be used. Note
that this default may be used in conjunction with
LADSPA_HINT_TOGGLED. */
#define LADSPA_HINT_DEFAULT_1 0x240
/* This default hint indicates that the number 100 should be used. */
#define LADSPA_HINT_DEFAULT_100 0x280
/* This default hint indicates that the Hz frequency of `concert A'
should be used. This will be 440 unless the host uses an unusual
tuning convention, in which case it may be within a few Hz. */
#define LADSPA_HINT_DEFAULT_440 0x2C0
#define LADSPA_IS_HINT_BOUNDED_BELOW(x) ((x) & LADSPA_HINT_BOUNDED_BELOW)
#define LADSPA_IS_HINT_BOUNDED_ABOVE(x) ((x) & LADSPA_HINT_BOUNDED_ABOVE)
#define LADSPA_IS_HINT_TOGGLED(x) ((x) & LADSPA_HINT_TOGGLED)
#define LADSPA_IS_HINT_SAMPLE_RATE(x) ((x) & LADSPA_HINT_SAMPLE_RATE)
#define LADSPA_IS_HINT_LOGARITHMIC(x) ((x) & LADSPA_HINT_LOGARITHMIC)
#define LADSPA_IS_HINT_INTEGER(x) ((x) & LADSPA_HINT_INTEGER)
#define LADSPA_IS_HINT_HAS_DEFAULT(x) ((x) & LADSPA_HINT_DEFAULT_MASK)
#define LADSPA_IS_HINT_DEFAULT_MINIMUM(x) (((x) & LADSPA_HINT_DEFAULT_MASK) \
== LADSPA_HINT_DEFAULT_MINIMUM)
#define LADSPA_IS_HINT_DEFAULT_LOW(x) (((x) & LADSPA_HINT_DEFAULT_MASK) \
== LADSPA_HINT_DEFAULT_LOW)
#define LADSPA_IS_HINT_DEFAULT_MIDDLE(x) (((x) & LADSPA_HINT_DEFAULT_MASK) \
== LADSPA_HINT_DEFAULT_MIDDLE)
#define LADSPA_IS_HINT_DEFAULT_HIGH(x) (((x) & LADSPA_HINT_DEFAULT_MASK) \
== LADSPA_HINT_DEFAULT_HIGH)
#define LADSPA_IS_HINT_DEFAULT_MAXIMUM(x) (((x) & LADSPA_HINT_DEFAULT_MASK) \
== LADSPA_HINT_DEFAULT_MAXIMUM)
#define LADSPA_IS_HINT_DEFAULT_0(x) (((x) & LADSPA_HINT_DEFAULT_MASK) \
== LADSPA_HINT_DEFAULT_0)
#define LADSPA_IS_HINT_DEFAULT_1(x) (((x) & LADSPA_HINT_DEFAULT_MASK) \
== LADSPA_HINT_DEFAULT_1)
#define LADSPA_IS_HINT_DEFAULT_100(x) (((x) & LADSPA_HINT_DEFAULT_MASK) \
== LADSPA_HINT_DEFAULT_100)
#define LADSPA_IS_HINT_DEFAULT_440(x) (((x) & LADSPA_HINT_DEFAULT_MASK) \
== LADSPA_HINT_DEFAULT_440)
typedef struct _LADSPA_PortRangeHint {
/* Hints about the port. */
LADSPA_PortRangeHintDescriptor HintDescriptor;
/* Meaningful when hint LADSPA_HINT_BOUNDED_BELOW is active. When
LADSPA_HINT_SAMPLE_RATE is also active then this value should be
multiplied by the relevant sample rate. */
LADSPA_Data LowerBound;
/* Meaningful when hint LADSPA_HINT_BOUNDED_ABOVE is active. When
LADSPA_HINT_SAMPLE_RATE is also active then this value should be
multiplied by the relevant sample rate. */
LADSPA_Data UpperBound;
} LADSPA_PortRangeHint;
/*****************************************************************************/
/* Plugin Handles:
This plugin handle indicates a particular instance of the plugin
concerned. It is valid to compare this to NULL (0 for C++) but
otherwise the host should not attempt to interpret it. The plugin
may use it to reference internal instance data. */
typedef void * LADSPA_Handle;
/*****************************************************************************/
/* Descriptor for a Type of Plugin:
This structure is used to describe a plugin type. It provides a
number of functions to examine the type, instantiate it, link it to
buffers and workspaces and to run it. */
typedef struct _LADSPA_Descriptor {
/* This numeric identifier indicates the plugin type
uniquely. Plugin programmers may reserve ranges of IDs from a
central body to avoid clashes. Hosts may assume that IDs are
below 0x1000000. */
unsigned long UniqueID;
/* This identifier can be used as a unique, case-sensitive
identifier for the plugin type within the plugin file. Plugin
types should be identified by file and label rather than by index
or plugin name, which may be changed in new plugin
versions. Labels must not contain white-space characters. */
const char * Label;
/* This indicates a number of properties of the plugin. */
LADSPA_Properties Properties;
/* This member points to the null-terminated name of the plugin
(e.g. "Sine Oscillator"). */
const char * Name;
/* This member points to the null-terminated string indicating the
maker of the plugin. This can be an empty string but not NULL. */
const char * Maker;
/* This member points to the null-terminated string indicating any
copyright applying to the plugin. If no Copyright applies the
string "None" should be used. */
const char * Copyright;
/* This indicates the number of ports (input AND output) present on
the plugin. */
unsigned long PortCount;
/* This member indicates an array of port descriptors. Valid indices
vary from 0 to PortCount-1. */
const LADSPA_PortDescriptor * PortDescriptors;
/* This member indicates an array of null-terminated strings
describing ports (e.g. "Frequency (Hz)"). Valid indices vary from
0 to PortCount-1. */
const char * const * PortNames;
/* This member indicates an array of range hints for each port (see
above). Valid indices vary from 0 to PortCount-1. */
const LADSPA_PortRangeHint * PortRangeHints;
/* This may be used by the plugin developer to pass any custom
implementation data into an instantiate call. It must not be used
or interpreted by the host. It is expected that most plugin
writers will not use this facility as LADSPA_Handle should be
used to hold instance data. */
void * ImplementationData;
/* This member is a function pointer that instantiates a plugin. A
handle is returned indicating the new plugin instance. The
instantiation function accepts a sample rate as a parameter. The
plugin descriptor from which this instantiate function was found
must also be passed. This function must return NULL if
instantiation fails.
Note that instance initialisation should generally occur in
activate() rather than here. */
LADSPA_Handle (*instantiate)(const struct _LADSPA_Descriptor * Descriptor,
unsigned long SampleRate);
/* This member is a function pointer that connects a port on an
instantiated plugin to a memory location at which a block of data
for the port will be read/written. The data location is expected
to be an array of LADSPA_Data for audio ports or a single
LADSPA_Data value for control ports. Memory issues will be
managed by the host. The plugin must read/write the data at these
locations every time run() or run_adding() is called and the data
present at the time of this connection call should not be
considered meaningful.
connect_port() may be called more than once for a plugin instance
to allow the host to change the buffers that the plugin is
reading or writing. These calls may be made before or after
activate() or deactivate() calls.
connect_port() must be called at least once for each port before
run() or run_adding() is called. When working with blocks of
LADSPA_Data the plugin should pay careful attention to the block
size passed to the run function as the block allocated may only
just be large enough to contain the block of samples.
Plugin writers should be aware that the host may elect to use the
same buffer for more than one port and even use the same buffer
for both input and output (see LADSPA_PROPERTY_INPLACE_BROKEN).
However, overlapped buffers or use of a single buffer for both
audio and control data may result in unexpected behaviour. */
void (*connect_port)(LADSPA_Handle Instance,
unsigned long Port,
LADSPA_Data * DataLocation);
/* This member is a function pointer that initialises a plugin
instance and activates it for use. This is separated from
instantiate() to aid real-time support and so that hosts can
reinitialise a plugin instance by calling deactivate() and then
activate(). In this case the plugin instance must reset all state
information dependent on the history of the plugin instance
except for any data locations provided by connect_port() and any
gain set by set_run_adding_gain(). If there is nothing for
activate() to do then the plugin writer may provide a NULL rather
than an empty function.
When present, hosts must call this function once before run() (or
run_adding()) is called for the first time. This call should be
made as close to the run() call as possible and indicates to
real-time plugins that they are now live. Plugins should not rely
on a prompt call to run() after activate(). activate() may not be
called again unless deactivate() is called first. Note that
connect_port() may be called before or after a call to
activate(). */
void (*activate)(LADSPA_Handle Instance);
/* This method is a function pointer that runs an instance of a
plugin for a block. Two parameters are required: the first is a
handle to the particular instance to be run and the second
indicates the block size (in samples) for which the plugin
instance may run.
Note that if an activate() function exists then it must be called
before run() or run_adding(). If deactivate() is called for a
plugin instance then the plugin instance may not be reused until
activate() has been called again.
If the plugin has the property LADSPA_PROPERTY_HARD_RT_CAPABLE
then there are various things that the plugin should not do
within the run() or run_adding() functions (see above). */
void (*run)(LADSPA_Handle Instance,
unsigned long SampleCount);
/* This method is a function pointer that runs an instance of a
plugin for a block. This has identical behaviour to run() except
in the way data is output from the plugin. When run() is used,
values are written directly to the memory areas associated with
the output ports. However when run_adding() is called, values
must be added to the values already present in the memory
areas. Furthermore, output values written must be scaled by the
current gain set by set_run_adding_gain() (see below) before
addition.
run_adding() is optional. When it is not provided by a plugin,
this function pointer must be set to NULL. When it is provided,
the function set_run_adding_gain() must be provided also. */
void (*run_adding)(LADSPA_Handle Instance,
unsigned long SampleCount);
/* This method is a function pointer that sets the output gain for
use when run_adding() is called (see above). If this function is
never called the gain is assumed to default to 1. Gain
information should be retained when activate() or deactivate()
are called.
This function should be provided by the plugin if and only if the
run_adding() function is provided. When it is absent this
function pointer must be set to NULL. */
void (*set_run_adding_gain)(LADSPA_Handle Instance,
LADSPA_Data Gain);
/* This is the counterpart to activate() (see above). If there is
nothing for deactivate() to do then the plugin writer may provide
a NULL rather than an empty function.
Hosts must deactivate all activated units after they have been
run() (or run_adding()) for the last time. This call should be
made as close to the last run() call as possible and indicates to
real-time plugins that they are no longer live. Plugins should
not rely on prompt deactivation. Note that connect_port() may be
called before or after a call to deactivate().
Deactivation is not similar to pausing as the plugin instance
will be reinitialised when activate() is called to reuse it. */
void (*deactivate)(LADSPA_Handle Instance);
/* Once an instance of a plugin has been finished with it can be
deleted using the following function. The instance handle passed
ceases to be valid after this call.
If activate() was called for a plugin instance then a
corresponding call to deactivate() must be made before cleanup()
is called. */
void (*cleanup)(LADSPA_Handle Instance);
} LADSPA_Descriptor;
/**********************************************************************/
/* Accessing a Plugin: */
/* The exact mechanism by which plugins are loaded is host-dependent,
however all most hosts will need to know is the name of shared
object file containing the plugin types. To allow multiple hosts to
share plugin types, hosts may wish to check for environment
variable LADSPA_PATH. If present, this should contain a
colon-separated path indicating directories that should be searched
(in order) when loading plugin types.
A plugin programmer must include a function called
"ladspa_descriptor" with the following function prototype within
the shared object file. This function will have C-style linkage (if
you are using C++ this is taken care of by the `extern "C"' clause
at the top of the file).
A host will find the plugin shared object file by one means or
another, find the ladspa_descriptor() function, call it, and
proceed from there.
Plugin types are accessed by index (not ID) using values from 0
upwards. Out of range indexes must result in this function
returning NULL, so the plugin count can be determined by checking
for the least index that results in NULL being returned. */
const LADSPA_Descriptor * ladspa_descriptor(unsigned long Index);
/* Datatype corresponding to the ladspa_descriptor() function. */
typedef const LADSPA_Descriptor *
(*LADSPA_Descriptor_Function)(unsigned long Index);
/**********************************************************************/
#ifdef __cplusplus
}
#endif
#endif /* LADSPA_INCLUDED */
/* EOF */

50
linalg.cpp Normal file
View File

@ -0,0 +1,50 @@
#include "linalg.h"
#include <algorithm>
#include <array>
#include <vector>
std::vector<float> Mengu::dsp::solve_sym_toeplitz(const std::vector<float> &cols, const std::vector<float> &y) {
// good ol dynamic programming to repeatedly do each step of levinson_recursion consequtively
std::vector<float> backward_vec = {1.0f / cols.at(0)}; // y.size assumed small, so use vec
backward_vec.reserve(y.size());
std::vector<float> forward_vec(1);
float error = 0.0f;
std::vector<float> result = {y.at(0) / cols.at(0)};
result.reserve(y.size());
for (int n = 1; n < y.size(); n++) {
// calculate the current error and backward vec
std::copy(backward_vec.crbegin(), backward_vec.crend(), forward_vec.begin());
error = dot(cols.data() + 1, backward_vec.data(), n);
float denom = 1.0f / (1 - error * error);
scalar_mul_inplace(denom, backward_vec.data(), backward_vec.size());
backward_vec.emplace(backward_vec.begin(), 0);
scalar_mul_inplace(-error * denom, forward_vec.data(), forward_vec.size());
forward_vec.push_back(0);
vec_add_inplace(forward_vec, backward_vec);
// calculate this iterations result
auto cols_iter = cols.crbegin() + (cols.size() - n - 1); // the elements n to 1
float result_error = dot(cols_iter, cols.crend() - 1, result.cbegin());
std::vector<float> scaled_backward_vec = scalar_mul(y[n] - result_error, backward_vec);
result.push_back(0);
vec_add_inplace(scaled_backward_vec, result);
}
return result;
}

148
linalg.h Normal file
View File

@ -0,0 +1,148 @@
/**
* @file linalg.h
* @author 9exa
* @brief objects and algorithms to solve linear equations
* @version 0.1
* @date 2023-05-04
*
* @copyright Copyright (c) 2023
*
*/
#ifndef MENGA_LINALG
#define MENGA_LINALG
#include <array>
#include <cstddef>
#include <vector>
#include <iostream>
namespace Mengu {
namespace dsp {
//2-D matrix
// dot_product
// c-style dot product
template<typename T>
inline T dot(const T *a, const T *b, const int size) {
T total = T();
for (int i = 0; i < size; i++) {
total += a[i] * b[i];
}
return total;
}
// do product on iterator
template<class InputIt1, class InputIt2>
inline float dot(InputIt1 afirst, InputIt1 alast, InputIt2 bfirst) {
float total = 0.0f;
while (afirst != alast) {
total += *afirst * *bfirst;
afirst++;
bfirst++;
}
return total;
}
template<typename T>
std::vector<T> scalar_mul(const T lambda, const std::vector<T> &a) {
std::vector<T> scaled;
for (T x: a) {
scaled.push_back(lambda * x);
}
return scaled;
}
template<typename T, size_t N>
std::array<T, N> scalar_mul(const T lambda, const T *a, size_t size) {
std::array<T, N> scaled;
for (size_t i = 0; i < size; i++) {
scaled[i] = lambda * a[i];
}
return scaled;
}
template<typename T>
void scalar_mul_inplace(const T lambda, T *a, const int size) {
for (int i = 0; i < size; i++) {
a[i] = lambda * a[i];
}
}
// add b to a inplace
template<typename T>
void vec_add_inplace(const std::vector<T> &b, std::vector<T> &a) {
for (int i = 0; i < a.size(); i++) {
a[i] += b.at(i);
}
}
template<typename T>
void vec_add_inplace(const T *b, T *a, const int size) {
for (int i = 0; i < size; i++) {
a[i] += b[i];
}
}
// Solves a system of equations where the matrix is a symmetric TOEplitz matrix,
// The matrix is represented as a list of its column values
std::vector<float> solve_sym_toeplitz(const std::vector<float> &cols, const std::vector<float> &y);
template <typename T, size_t N>
std::array<T, N> solve_sym_toeplitz(const std::array<T, N> &cols, const std::array<T, N> &y) {
// good ol dynamic programming to repeatedly do each step of levinson_recursion consequtively
std::array<T, N> backward_vec;
std::array<T, N> forward_vec;
backward_vec[backward_vec.size() - 1] = forward_vec[0] = 1.0 / cols.at(0);
float error = 0.0f;
std::array<T, N> result;
result[0] = {y.at(0) / cols.at(0)};
for (int n = 1; n < N; n++) {
// calculate the current error and backward vec
std::copy(backward_vec.crbegin(), backward_vec.crbegin() + n, forward_vec.begin());
error = dot(cols.data() + 1, backward_vec.data() + (N - n), n);
float denom = 1.0f / (1 - error * error);
scalar_mul_inplace(denom, backward_vec.data() + (N - n), n);
backward_vec[N - 1 - n] = T();
scalar_mul_inplace(-error * denom, forward_vec.data(), n);
forward_vec[n] = T();
vec_add_inplace(forward_vec.data(), backward_vec.data() + (N - n - 1), n + 1);
// calculate this iterations result
auto cols_iter = cols.crbegin() + (cols.size() - n - 1); // the elements n to 0
float result_error = dot(cols_iter, cols.crend() - 1, result.cbegin());
std::array<T, N> scaled_backward_vec = scalar_mul<T, N>(
y[n] - result_error,
backward_vec.data() + (N - n - 1),
n + 1);
result[n] = T();
vec_add_inplace(scaled_backward_vec.data(), result.data(), n + 1);
}
return result;
}
}
}
//
#endif

67
loudness.cpp Normal file
View File

@ -0,0 +1,67 @@
#include "loudness.h"
#include "filter.h"
#include "common.h"
#include "mengumath.h"
#include <complex>
#include <cstdint>
// Coefficients for LUFS_freq filteters
// High shelf filter stage
#define S1_A1 -1.69065929318241
#define S1_A2 0.73248077421585
#define S1_B0 1.53512485958697
#define S1_B1 -2.69169618940638
#define S1_B2 1.19839281085285
// High pass Filter stage
#define S2_A1 -1.99004745483398
#define S2_A2 0.99007225036621
#define S2_B0 1.0
#define S2_B1 -2.0
#define S2_B2 1.0
#define LUFS_DEFAULT_SAMPLE_RATE 48000
using namespace Mengu;
using namespace dsp;
float Mengu::dsp::LUFS_freq(Complex *freqs, uint32_t size, uint32_t sample_rate) {
// Adjust to the custum sample rate
float sample_rate_correction = (float) sample_rate / LUFS_DEFAULT_SAMPLE_RATE * 2;
float total_amp2 = 0.0f; // total squared amplitude of each freq bin
for (uint32_t i = 0; i < size; i++) {
// argument for the transfer functions
float omega = (float) MATH_PI * i / size * sample_rate_correction;
Complex z = std::polar(1.0f, omega);
// filtered frequency
Complex y = freqs[i] * quad_filter_trans(z, S1_A1, S1_A2, S1_B0, S1_B1, S1_B2);
y = y * quad_filter_trans(z, S2_A1, S2_A2, S2_B0, S2_B1, S2_B2);
total_amp2 += std::norm(y);
}
return -0.691 + 10 * Mengu::log10(total_amp2 / size);
}
Complex Mengu::dsp::LUFS_filter_transfer(float freq) {
float omega = (float) MATH_PI * freq / LUFS_DEFAULT_SAMPLE_RATE * 2;
Complex z = std::polar(1.0f, omega);
Complex y = quad_filter_trans(z, S1_A1, S1_A2, S1_B0, S1_B1, S1_B2);
return y * quad_filter_trans(z, S2_A1, S2_A2, S2_B0, S2_B1, S2_B2);
}
LUFSFilter::LUFSFilter():
_high_shelf_filter(S1_A1, S1_A2, S1_B0, S1_B1, S1_B2),
_high_pass_filter(S2_A1, S2_A2, S2_B0, S2_B1, S2_B2) {}
void LUFSFilter::transform(const float *input, float *output, uint32_t size) {
_high_shelf_filter.transform(input, output, size);
_high_pass_filter.transform(output, output, size);
}
void LUFSFilter::reset() {
_high_shelf_filter.reset();
_high_pass_filter.reset();
}

105
loudness.h Normal file
View File

@ -0,0 +1,105 @@
/**
* @file loudness.h
* @author your name (you@domain.com)
* @brief Functions for calculating loudness in a signal
* @version 0.1
* @date 2023-06-10
*
* @copyright Copyright (c) 2023
*
*/
#ifndef MENGU_LOUDNESS
#define MENGU_LOUDNESS
#include "common.h"
#include "filter.h"
#include <cstdint>
#include <stdint.h>
namespace Mengu {
namespace dsp {
// Loudness Units relative to Full Scale of a sample in the frequency domain
// first (positive) half of the frequency spectrum only
// Only supports 1 channel
float LUFS_freq(Complex *freqs, uint32_t size, uint32_t sample_rate = 48000);
// The value of the transfer function associated with the described frequency bin
Complex LUFS_filter_transfer(float freq);
// Performs the filter step associated with LUFS
class LUFSFilter {
public:
LUFSFilter();
void transform(const float *input, float *output, uint32_t size);
void reset();
private:
BiquadFilter _high_shelf_filter; // Stage 1 filter
BiquadFilter _high_pass_filter; // Stage 2 filter
};
// Scales a raw sample to have the same Loudness (in LUFS) as a reference sample
// The raw sample and reference sample are assumed to come from their own persistant samples
template<typename T, uint32_t N, int DefCorrection = 1>
class LoudnessNormalizer {
public:
void normalize(const T *raw_sample, const T *reference_sample, T *output) {
float filtered_raw[N];
float filtered_reference[N];
for (uint32_t i = 0; i < N; i++) {
filtered_raw[i] = _as_float(raw_sample[i]);
filtered_reference[i] = _as_float(reference_sample[i]);
}
// perform filter
_raw_sample_filter.transform(filtered_raw, filtered_raw, N);
_reference_sample_filter.transform(filtered_reference, filtered_reference, N);
// Get the (unormalized) power of each filtered sample
float raw_power = 0.0f;
float reference_power = 0.0f;
for(uint32_t i = 0; i < N; i++) {
raw_power += filtered_raw[i] * filtered_raw[i];
reference_power += filtered_reference[i] * filtered_reference[i];
}
float correction = sqrt(reference_power / raw_power);
// std::cout << correction << ", " << raw_power << ", " << shifted_power << std::endl;
if (!std::isfinite(correction)) {
// just correct for 1/2 freq spectrum sampling if correction is not real
for (uint32_t i = 0; i < N; i++) {
output[i] = raw_sample[i] * (float) DefCorrection;
}
}
else {
// otherwise scale the shifted by the correction
for (uint32_t i = 0; i < N; i++) {
output[i] = raw_sample[i] * correction;
}
}
}
// resets memory on previous raw and reference samples
void reset() {
_raw_sample_filter.reset();
_reference_sample_filter.reset();
}
private:
LUFSFilter _raw_sample_filter;
LUFSFilter _reference_sample_filter;
static float _as_float(T v) {
return Complex(v).real();
}
};
}
}
#endif

250
mengubah_lv2.cpp Normal file
View File

@ -0,0 +1,250 @@
#include <cstdint>
#include "ladspa.h"
#include <array>
#include "common.h"
#include "effect.h"
#include "pitchshifter.h"
#include "formantshifter.h"
#include "timestretcher.h"
#define PORT_COUNT 6
#define ML_INPUT 0
#define ML_OUTPUT 1
#define ML_PITCH_CHANGE 3
#define ML_TIMBRE_CHANGE 5
#define ML_TIMBRE_IND 4
#define ML_PITCH_IND 2
using namespace Mengu;
using namespace dsp;
struct PluginHandler {
std::array<PitchShifter *, 3> pitch_shifters;
std::array<Effect *, 2> formant_shifters;
const float *in_buffer;
float *out_buffer;
const float *pitch_shifter_ind;
const float *pitch_shift;
const float *formant_shifter_ind;
const float *formant_shift;
};
enum PitchShiftInd {
WSOLAPitch = 0,
PhaseVocoderPitch = 1,
PhaseVocoderDoneRightPitch = 2,
};
enum FormantShiftInd {
LPCFormant = 0,
PSOLAFormant = 1,
};
/* internal core methods */
static LADSPA_Handle instantiate(const LADSPA_Descriptor *descriptor, unsigned long sample_rate) {
PluginHandler *plugin = new PluginHandler;
plugin->pitch_shifters = {
new TimeStretchPitchShifter(new WSOLATimeStretcher(), 1),
new TimeStretchPitchShifter(new PhaseVocoderTimeStretcher(), 1),
new TimeStretchPitchShifter(new PhaseVocoderDoneRightTimeStretcher(), 1),
};
plugin->formant_shifters = {
new LPCFormantShifter(),
new TimeStretchPitchShifter(new PSOLATimeStretcher(), 1),
};
return plugin;
}
static void connect_port(LADSPA_Handle instance, unsigned long port, LADSPA_Data *data_location) {
PluginHandler* plugin = (PluginHandler*) instance;
if (plugin == nullptr) return;
switch (port)
{
case 0:
plugin->in_buffer = (const float*) data_location;
break;
case 1:
plugin->out_buffer = (float*) data_location;
break;
case 2:
plugin->pitch_shifter_ind = (const float*) data_location;
case 3:
plugin->pitch_shift = (const float*) data_location;
break;
case 4:
plugin->formant_shifter_ind = (const float*) data_location;
break;
case 5:
plugin->formant_shift = (const float*) data_location;
break;
default:
break;
}
}
static void activate (LADSPA_Handle instance)
{
/* not needed here */
}
static void run (LADSPA_Handle instance, unsigned long sample_count)
{
PluginHandler* plugin = (PluginHandler*) instance;
if (plugin == nullptr) return;
if ((!plugin->in_buffer) || (!plugin->out_buffer) || (!plugin->pitch_shift)) return;
// apply effects
Effect *pitch_shifter = plugin->pitch_shifters[static_cast<PitchShiftInd>(*plugin->pitch_shifter_ind)];
Effect *formant_shifter = plugin->formant_shifters[static_cast<FormantShiftInd>(*plugin->formant_shifter_ind)];
pitch_shifter->set_property(0, EffectPropPayload {
.type = Slider,
.value = *plugin->pitch_shift,
});
formant_shifter->set_property(0, EffectPropPayload {
.type = Slider,
.value = *plugin->formant_shift,
});
static constexpr uint32_t ProcSize = 1 << 10;
Complex cbuffer[ProcSize];
uint32_t num_processed = 0;
while (num_processed < sample_count) {
uint32_t num_this_pass = MIN(ProcSize, sample_count - num_processed);
for (uint32_t i = 0; i < num_this_pass; i++) {
cbuffer[i] = (Complex) plugin->in_buffer[num_processed + i];
}
pitch_shifter->push_signal(cbuffer, num_this_pass);
pitch_shifter->pop_transformed_signal(cbuffer, num_this_pass);
formant_shifter->push_signal(cbuffer, num_this_pass);
formant_shifter->pop_transformed_signal(cbuffer, num_this_pass);
for (uint32_t i = 0; i < num_this_pass; i++) {
plugin->out_buffer[num_processed + i] = cbuffer[i].real();
}
num_processed += num_this_pass;
}
}
static void deactivate (LADSPA_Handle instance)
{
/* not needed here */
}
static void cleanup (LADSPA_Handle instance) {
PluginHandler *plugin = (PluginHandler *) instance;
if (plugin == nullptr) {
return;
}
for (auto pitch_shifter: plugin->pitch_shifters) {
delete pitch_shifter;
}
for (auto formant_shifter: plugin->formant_shifters) {
delete formant_shifter;
}
delete plugin;
}
static LADSPA_Descriptor *descriptor = NULL;
// On plugin load
static void __attribute__ ((constructor)) init() {
LADSPA_PortDescriptor * portDescriptors;
LADSPA_PortRangeHint * portRangeHints;
char ** portNames;
descriptor = (LADSPA_Descriptor *)malloc(sizeof(LADSPA_Descriptor));
if (!descriptor) {
return;
}
descriptor->UniqueID = 109124; // should be unique
descriptor->Label = "pitcher";
descriptor->Name = "pitchplug";
descriptor->Maker = "Rysertio";
descriptor->Copyright = "(c) 2023, Rysertio";
descriptor->Properties = LADSPA_PROPERTY_HARD_RT_CAPABLE;
descriptor->PortCount = PORT_COUNT;
portDescriptors = (LADSPA_PortDescriptor *) calloc(PORT_COUNT, sizeof(LADSPA_PortDescriptor));
portDescriptors[ML_PITCH_CHANGE] = LADSPA_PORT_INPUT | LADSPA_PORT_CONTROL;
portDescriptors[ML_TIMBRE_CHANGE] = LADSPA_PORT_INPUT | LADSPA_PORT_CONTROL;
portDescriptors[ML_INPUT] = LADSPA_PORT_INPUT | LADSPA_PORT_AUDIO;
portDescriptors[ML_OUTPUT] = LADSPA_PORT_OUTPUT | LADSPA_PORT_AUDIO;
portDescriptors[ML_PITCH_IND] = LADSPA_PORT_INPUT | LADSPA_PORT_CONTROL;
portDescriptors[ML_TIMBRE_IND] = LADSPA_PORT_INPUT | LADSPA_PORT_CONTROL;
descriptor->PortDescriptors = portDescriptors;
portNames = (char **) calloc(PORT_COUNT, sizeof(char *));
portNames[ML_PITCH_CHANGE] = "pitch change factor";
portNames[ML_TIMBRE_CHANGE] = "timbre change factor";
portNames[ML_INPUT] = "Input";
portNames[ML_OUTPUT] = "Output";
portNames[ML_PITCH_IND] = "pitch change algorithm";
portNames[ML_TIMBRE_IND] = "timbre change algorithm";
descriptor->PortNames = (const char * const *) portNames;
portRangeHints = (LADSPA_PortRangeHint *) calloc(PORT_COUNT, sizeof(LADSPA_PortRangeHint));
portRangeHints[ML_PITCH_CHANGE].HintDescriptor = LADSPA_HINT_BOUNDED_BELOW | LADSPA_HINT_DEFAULT_1;
portRangeHints[ML_PITCH_CHANGE].LowerBound = 1;
portRangeHints[ML_TIMBRE_CHANGE].HintDescriptor = LADSPA_HINT_BOUNDED_BELOW | LADSPA_HINT_DEFAULT_1;
portRangeHints[ML_TIMBRE_CHANGE].LowerBound = 1;
portRangeHints[ML_INPUT].HintDescriptor = 0;
portRangeHints[ML_OUTPUT].HintDescriptor = 0;
descriptor->PortRangeHints = portRangeHints;
descriptor->instantiate = instantiate;
descriptor->connect_port = connect_port;
descriptor->activate = activate;
descriptor->run = run;
descriptor->run_adding = NULL;
descriptor->set_run_adding_gain = NULL;
descriptor->deactivate = NULL;
descriptor->cleanup = cleanup;
}
// On plugin unload
static void __attribute__ ((destructor)) fini() {
return;
}
// we only have one type of plugin
const LADSPA_Descriptor * ladspa_descriptor(unsigned long index) {
if (index != 0) {
return NULL;
}
return descriptor;
}

164
mengumath.h Normal file
View File

@ -0,0 +1,164 @@
#ifndef MENGA_MATH
#define MENGA_MATH
#include <cmath>
#include <cstdint>
namespace Mengu {
// Copied from Godot
#define MATH_PI 3.1415926535897932
#define MATH_TAU 6.283185307179586
// Make room for our constexpr's below by overriding potential system-specific macros.
#undef ABS
#undef SIGN
#undef MIN
#undef MAX
#undef CLAMP
// Generic ABS function, for math uses please use Math::abs.
template <typename T>
constexpr T ABS(T m_v) {
return m_v < 0 ? -m_v : m_v;
}
template <typename T>
constexpr const T SIGN(const T m_v) {
return m_v == 0 ? 0.0f : (m_v < 0 ? -1.0f : +1.0f);
}
template <typename T, typename T2>
constexpr auto MIN(const T m_a, const T2 m_b) {
return m_a < m_b ? m_a : m_b;
}
template <typename T, typename T2>
constexpr auto MAX(const T m_a, const T2 m_b) {
return m_a > m_b ? m_a : m_b;
}
template <typename T, typename T2, typename T3>
constexpr auto CLAMP(const T m_a, const T2 m_min, const T3 m_max) {
return m_a < m_min ? m_min : (m_a > m_max ? m_max : m_a);
}
//// #####
template <typename T>
inline T min(T a, T b) {
return a < b ? a : b;
}
template <typename T>
inline T max(T a, T b) {
return a > b ? a : b;
}
template <typename T>
inline T sign(T x) {
return static_cast<T>(SIGN(x));
}
template <typename T>
inline T abs(T x) {
return std::abs(x);
}
template <typename T>
inline T pow(T x, T e) {
return std::pow(x, e);
}
template <typename T>
inline T exp2(T x) {
return std::exp2(x);
}
template <typename T>
inline T log2(T x) {
return std::log2(x);
}
template <typename T>
inline T log10(T x) {
return std::log10(x);
}
template <typename T>
inline T logb(T base, T x) {
return std::log2(x) / std::log2(base);
}
template <typename T>
inline T sqrt(T x) {
return std::sqrt(x);
}
template <typename T>
inline T lerp(T min_v, T max_v, float w) {
return min_v + w * (max_v - min_v);
}
inline float inverse_lerp(float min_v, float max_v, float v) {
return (v - min_v) / (max_v - min_v);
}
// linearly interpelates between the shortest path between 2 angles. Does NOT gaurentee th output be within [-pi, pi]
inline float lerp_angle(float min_v, float max_v, float w) {
if ((max_v - min_v) > MATH_PI) { min_v += MATH_TAU; }
else if ((min_v - max_v) > MATH_PI) { max_v += MATH_TAU; }
return lerp(min_v, max_v, w);
}
template <typename T>
inline T modf(T x, T *iptr) {
return std::modf(x, iptr);
}
inline int posmod(int x, int y) {
int v = x % y;
if (v < 0) {
v += y;
}
return v;
}
inline float fposmod(float x, float y) {
float v = std::fmod(x, y);
if (v * y < 0.0f) {
v += y;
}
return v;
}
inline bool is_pow_2(uint32_t x) {
return x != 0 && ((x - 1) & x) == 0;
}
inline uint32_t last_pow_2(uint32_t x) {
x |= (x >> 1);
x |= (x >> 2);
x |= (x >> 4);
x |= (x >> 8);
x |= (x >> 16);
return x - (x >> 1);
}
inline uint32_t next_pow_2(uint32_t x) {
x |= (x >> 1);
x |= (x >> 2);
x |= (x >> 4);
x |= (x >> 8);
x |= (x >> 16);
return (x << 1) & ~x;
}
};
#endif

92536
miniaudio.h Normal file

File diff suppressed because it is too large Load Diff

291
pitchshifter.cpp Normal file
View File

@ -0,0 +1,291 @@
#include "common.h"
#include "effect.h"
#include "timestretcher.h"
#include "sampling.h"
#include "fft.h"
#include "interpolation.h"
#include "vecdeque.h"
#include <algorithm>
#include <array>
#include <cmath>
#include <complex>
#include <cstdint>
#include "pitchshifter.h"
#include "singletons.h"
#include <iterator>
#include "mengumath.h"
#include <iostream>
#include <vector>
using namespace Mengu;
using namespace dsp;
std::vector<EffectPropDesc> PitchShifter::get_property_descs() const {
return {
EffectPropDesc {
.type = EffectPropType::Slider,
.name = "Pitch Shift",
.desc = "Scales the pitch of pushed signals by this amount",
.slider_data = {
.min_value = 0.5,
.max_value = 2,
.scale = Exp,
}
}
};
}
void PitchShifter::set_property(uint32_t id, EffectPropPayload data) {
if (id == 0) {
if (data.type == Slider) {
set_shift_factor(data.value);
}
}
}
EffectPropPayload PitchShifter::get_property(uint32_t id) const {
if (id == 0) {
return EffectPropPayload {
.type = Slider,
.value = _shift_factor,
};
}
return EffectPropPayload {
.type = Slider,
.value = 0.0f,
};
}
PhaseVocoderPitchShifterV2::PhaseVocoderPitchShifterV2() {
_transformed_buffer.resize(ProcSize);
}
PhaseVocoderPitchShifterV2::~PhaseVocoderPitchShifterV2() {}
void PhaseVocoderPitchShifterV2::push_signal(const Complex *input, const uint32_t &size) {
_raw_buffer.extend_back(input, size);
}
uint32_t PhaseVocoderPitchShifterV2::pop_transformed_signal(Complex *output, const uint32_t &size) {
while ((_transformed_buffer.size() < size + OverlapSize) && (_raw_buffer.size() >= ProcSize)) {
std::array<Complex, ProcSize> samples;
_raw_buffer.to_array(samples.data(), ProcSize);
_lpc.load_sample(samples.data());
const std::array<float, ProcSize> &envelope = _lpc.get_envelope();
const std::array<float, ProcSize> &residuals = _lpc.get_residuals();
const std::array<Complex, ProcSize> &frequencies = _lpc.get_freq_spectrum();
// std::array<float, ProcSize / 2> log_envelope{};
std::array<float, ProcSize / 2> log_residuals{};
std::transform(residuals.cbegin(), residuals.cend(), log_residuals.begin(), [] (float f) {
return log2(f);
});
std::array<float, ProcSize / 2> args{};
std::transform(frequencies.cbegin(), frequencies.cend(), args.begin(), [] (Complex freq) {
// shifted to be positive to make lerping easier
return std::arg(freq);
});
/*
// resample in the log-frequency domain
std::array<Complex, ProcSize / 2> new_freq{};
if (_shift_factor > 1.0f) {
linear_resample_no_filter(log_residuals.data(), log_residuals.data(), ProcSize, _shift_factor, -2e5f);
linear_resample_no_filter(args.data(), args.data(), ProcSize, _shift_factor, 0.0f);
}
else {
linear_resample_no_filter(log_residuals.data(), log_residuals.data(), ProcSize, _shift_factor, -2e5f);
linear_resample_no_filter(args.data(), args.data(), ProcSize, _shift_factor, 0.0f);
// 0 the rest
for (uint32_t i = ProcSize * _shift_factor; i < ProcSize; i++) {
log_residuals[i] = -2e5f;
args[i] = 0.0f;
}
}
std::transform(log_residuals.cbegin(), log_residuals.cend(), new_freq.begin(), [] (float f) {
return exp2(f);
});
std::transform(new_freq.cbegin(), new_freq.cend(), envelope.cbegin(), new_freq.begin(),
[] (Complex resid, float env) {
return resid * env;
}
);
// rescale the frequencies
float max_base_freq = 0.0f;
float max_new_freq = 0.0f;
for (uint32_t i = 0; i < ProcSize / 2; i++) {
float freq_amp = std::norm(frequencies[i]);
if (max_base_freq < freq_amp) { max_base_freq = freq_amp; }
if (max_new_freq < new_freq[i].real()) { max_new_freq = new_freq[i].real(); }
}
max_base_freq = std::sqrt(max_base_freq); // rooting delayed until after the loop
std::transform(new_freq.cbegin(), new_freq.cend(), new_freq.begin(),
[max_base_freq, max_new_freq] (Complex freq) {
return freq * max_base_freq / max_new_freq;
}
);
// adjust phases
for (uint32_t i = 0; i < ProcSize / 2; i++) {
// new_freq[i] = std::polar(std::sqrt(std::norm(frequencies[i])), std::arg(frequencies[i]));
new_freq[i] = std::polar(new_freq[i].real(), args[i]);
// new_freq[i] = std::polar(new_freq[i].real(), (float) -MATH_TAU * (_samples_processed * i) / ProcSize);
}
*/
std::array<Complex, ProcSize / 2> new_freq{};
std::array<float, ProcSize / 2> freq_mags{};
std::transform(frequencies.cbegin(), frequencies.cend(), freq_mags.begin(), [] (Complex c) {
return std::sqrt(std::norm(c));
});
if (_shift_factor > 1.0f) {
linear_resample_no_filter(frequencies.data(), new_freq.data(), ProcSize / 2, _shift_factor);
linear_resample_no_filter(freq_mags.data(), freq_mags.data(), ProcSize / 2, _shift_factor);
// linear_resample_no_filter(args.data(), args.data(), ProcSize, _shift_factor);
}
else {
linear_resample_no_filter(frequencies.data(), new_freq.data(), ProcSize / 2, _shift_factor);
linear_resample_no_filter(freq_mags.data(), freq_mags.data(), ProcSize / 2, _shift_factor);
// linear_resample_no_filter(args.data(), args.data(), ProcSize, _shift_factor);
// 0 the rest
for (uint32_t i = ProcSize * _shift_factor; i < ProcSize; i++) {
freq_mags[i] = 0.0f;
// args[i] = 0.0f;
}
}
// std::transform(freq_mags.cbegin(), freq_mags.cend(), new_freq.cbegin(), new_freq.begin(),
// [] (float mag, Complex old_freq) { return old_freq / std::sqrt(std::norm(old_freq)) * mag; }
// );
// scale by formants
// float envelope_max = -2e5;
// float new_freq_max = -2e5;
// float scaled_freq_max = -2e5;
// for (uint32_t i = 0; i < ProcSize / 2; i++) {
// if (envelope_max < envelope[i]) { envelope_max = envelope[i]; }
// if (new_freq_max < std::norm(new_freq[i])) { new_freq_max = std::norm(new_freq[i]); }
// }
// new_freq_max = std::sqrt(new_freq_max);
// for (uint32_t i = 0; i < ProcSize / 2; i++) {
// new_freq[i] = new_freq[i] * envelope[i] / envelope_max;
// if (scaled_freq_max < std::norm(new_freq[i])) { scaled_freq_max = std::norm(new_freq[i]); }
// }
// scaled_freq_max = std::sqrt(scaled_freq_max);
// for (uint32_t i = 0; i < ProcSize / 2; i++) {
// new_freq[i] *= new_freq_max / scaled_freq_max;
// }
_lpc.get_fft().inverse_transform(new_freq.data(), samples.data());
mix_and_extend(_transformed_buffer, samples, OverlapSize, hann_window);
_samples_processed += ProcSize - OverlapSize;
_raw_buffer.pop_front_many(nullptr, ProcSize - OverlapSize);
}
uint32_t n = _transformed_buffer.pop_front_many(output, size);
return n;
}
uint32_t PhaseVocoderPitchShifterV2::n_transformed_ready() const {
return _transformed_buffer.size() - ProcSize;
}
void PhaseVocoderPitchShifterV2::reset() {
_transformed_buffer.resize(ProcSize, Complex(0.0f));
}
TimeStretchPitchShifter::TimeStretchPitchShifter(TimeStretcher *stretcher, uint32_t nchannels):
_stretcher(stretcher),
_resampler(nchannels, 1.0f) {
_stretcher->set_stretch_factor(1.0f);
_shift_factor = 1.0f;
// Complex zeros[MinResampleInputSize] = {Complex()};
// _stretcher.push_signal(zeros, MinResampleInputSize);
// _pitch_shifting_stretcher.push_signal(zeros, MinResampleInputSize);
// _raw_buffer.resize(MinResampleInputSize * 2, 0);
}
TimeStretchPitchShifter::~TimeStretchPitchShifter() {
delete _stretcher;
}
void TimeStretchPitchShifter::push_signal(const Complex *input, const uint32_t &size) {
// _raw_buffer.extend_back(input, size);
_stretcher->push_signal(input, size);
// _pitch_shifting_stretcher.push_signal(input, size);
}
uint32_t TimeStretchPitchShifter::pop_transformed_signal(Complex *output, const uint32_t &size) {
// do transform, eagerly
// resample the time-stretch pitch shifted samples
// while (_pitch_shifting_stretcher.n_transformed_ready() >= MinResampleInputSize) {
bool can_still_process = true;
while (can_still_process && n_transformed_ready() < size) {
const uint32_t desired_stretched_size = size * _shift_factor;
std::vector<Complex> stretched(desired_stretched_size);
const uint32_t actually_stretched = _stretcher->pop_transformed_signal(stretched.data(), desired_stretched_size);
stretched.resize(actually_stretched);
can_still_process = actually_stretched > 0;
std::vector<Complex> unstretched = _resampler.resample(stretched);
_transformed_buffer.extend_back(unstretched.data(), unstretched.size());
}
uint32_t n = _transformed_buffer.pop_front_many(output, size);
// std::cout << "n " << n << std::endl;
// compensate for drift
if (n_transformed_ready() > IncreaseResampleThreshold) {
static constexpr float ResampleHalflife = IncreaseResampleThreshold * 1.5;
float resample_factor = exp2(-(int)(n_transformed_ready() - IncreaseResampleThreshold) / ResampleHalflife);
_stretcher->set_stretch_factor(_shift_factor * resample_factor);
}
else if (n_transformed_ready() < StandardResampleThreshold) {
_stretcher->set_stretch_factor(_shift_factor);
}
for (uint32_t i = n; i < size; i++) {
output[i] = 0;
}
return n;
}
uint32_t TimeStretchPitchShifter::n_transformed_ready() const {
return _transformed_buffer.size();
}
void TimeStretchPitchShifter::reset() {
_raw_buffer.resize(0);
_transformed_buffer.resize(0);
_stretcher->reset();
}
void TimeStretchPitchShifter::set_shift_factor(const float &shift_factor) {
_shift_factor = shift_factor;
_stretcher->set_stretch_factor(shift_factor);
_resampler.set_stretch_factor(shift_factor);
}

137
pitchshifter.h Normal file
View File

@ -0,0 +1,137 @@
/*
* (Real time) Pitch shifting object implementations
*/
#ifndef MENGA_PITCH_SHIFTER
#define MENGA_PITCH_SHIFTER
#include "correlation.h"
#include "effect.h"
#include "sampling.h"
#include "timestretcher.h"
#include "fft.h"
#include <cstdint>
#include "common.h"
#include "cyclequeue.h"
#include "vecdeque.h"
#include <vector>
namespace Mengu {
namespace dsp {
class PitchShifter: public Effect {
public:
~PitchShifter() {}
virtual InputDomain get_input_domain() override {
return InputDomain::Time;
}
virtual void set_shift_factor(const float &factor) {
_shift_factor = factor;
};
virtual std::vector<EffectPropDesc> get_property_descs() const override;
virtual void set_property(uint32_t id, EffectPropPayload data) override;
virtual EffectPropPayload get_property(uint32_t id) const override;
protected:
float _shift_factor = 1.0f;
};
// Shifter in the frequency domain that uses LPC to estimate formant preservation.
class PhaseVocoderPitchShifterV2: public PitchShifter {
public:
PhaseVocoderPitchShifterV2();
~PhaseVocoderPitchShifterV2();
virtual void push_signal(const Complex *input, const uint32_t &size) override;
virtual uint32_t pop_transformed_signal(Complex *output, const uint32_t &size) override;
virtual uint32_t n_transformed_ready() const override;
virtual void reset() override;
private:
// raw time-domain data
VecDeque<Complex> _raw_buffer;
// time domain data after pitch_shift
VecDeque<Complex> _transformed_buffer;
static constexpr uint32_t ProcSize = 1 << 9;
static constexpr uint32_t OverlapSize = 1 << 6;
LPC<ProcSize, 30> _lpc;
uint32_t _samples_processed;
};
// Shifts by resampling a time stretcher
class TimeStretchPitchShifter: public PitchShifter {
public:
TimeStretchPitchShifter(TimeStretcher *stretcher, uint32_t nchannels);
~TimeStretchPitchShifter();
virtual void push_signal(const Complex *input, const uint32_t &size) override;
virtual uint32_t pop_transformed_signal(Complex *output, const uint32_t &size) override;
virtual uint32_t n_transformed_ready() const override;
virtual void reset() override;
virtual void set_shift_factor(const float &factor) override;
protected:
float _formant_shift = 1.0f;
private:
static constexpr uint32_t MinResampleInputSize = 1 << 10;
// Size of transformed buffer before we increase resampling to compensate for drift
static constexpr uint32_t IncreaseResampleThreshold = 5000;
static constexpr uint32_t StandardResampleThreshold = 3000;
// raw time-domain data
VecDeque<Complex> _raw_buffer;
// time domain data after pitch_shift
VecDeque<Complex> _transformed_buffer;
LinearResampler _resampler;
TimeStretcher *_stretcher;
};
/*
class PTimeStretchPitchShifter: public PitchShifter {
public:
PTimeStretchPitchShifter(uint32_t nchannels);
~PTimeStretchPitchShifter();
virtual void push_signal(const Complex *input, const uint32_t &size) override;
virtual uint32_t pop_transformed_signal(Complex *output, const uint32_t &size) override;
virtual uint32_t n_transformed_ready() override;
virtual void reset() override;
virtual void set_shift_factor(const float &factor) override;
private:
static constexpr uint32_t MinResampleInputSize = 1 << 10;
// Size of transformed buffer before we increase resampling to compensate for drift
static constexpr uint32_t IncreaseResampleThreshold = 10000;
static constexpr uint32_t StandardResampleThreshold = 4000;
// raw time-domain data
VecDeque<Complex> _raw_buffer;
// time domain data after pitch_shift
VecDeque<Complex> _transformed_buffer;
LinearResampler _resampler;
PSOLATimeStretcher _pitch_shifting_stretcher;
};
*/
}; // namespace dsp
}; // namespace Mengu
#endif

78
sampling.cpp Normal file
View File

@ -0,0 +1,78 @@
#include "sampling.h"
#include "common.h"
#include "mengumath.h"
#include "miniaudio.h"
#include <algorithm>
#include <cmath>
#include <cstdint>
#include <iterator>
#include <math.h>
#include <vector>
#include <iostream>
using namespace Mengu;
using namespace dsp;
constexpr uint32_t DefaultSampleRate = 44100;
// void Mengu::dsp::linear_resample_no_filter(const float *input, float *output, uint32_t size, float shift_factor, float def_val)
LinearResampler::LinearResampler(uint32_t nchannels, float stretch_factor) :
_stretch_factor(stretch_factor),
_nchannels(nchannels) {
ma_linear_resampler_config resampler_config = ma_linear_resampler_config_init(
ma_format_f32,
nchannels,
DefaultSampleRate,
DefaultSampleRate * stretch_factor
);
ma_result result = ma_linear_resampler_init(&resampler_config, nullptr, &_resampler);
if (result != MA_SUCCESS) {
std::string err_msg ("Could not create resampler");
err_msg += std::to_string(result);
throw std::runtime_error(err_msg);
}
}
LinearResampler::~LinearResampler() {
ma_linear_resampler_uninit(&_resampler, nullptr);
}
void LinearResampler::set_stretch_factor(float stretch_factor) {
_stretch_factor = stretch_factor;
ma_linear_resampler_set_rate_ratio(&_resampler, stretch_factor);
}
std::vector<Complex> LinearResampler::resample(const std::vector<Complex> &samples) {
std::vector<float> fsamples;
std::transform(samples.cbegin(), samples.cend(), std::back_inserter(fsamples),
[] (Complex f) { return f.real(); }
);
ma_uint64 input_size = fsamples.size();
ma_uint64 output_size;
ma_linear_resampler_get_expected_output_frame_count(&_resampler, samples.size(), &output_size);
std::vector<float> foutput(output_size);
ma_result result = ma_linear_resampler_process_pcm_frames(&_resampler, fsamples.data(), &input_size, foutput.data(), &output_size);
if (result != MA_SUCCESS) {
std::string err_msg ("Could not perform resample");
err_msg += std::to_string(result);
throw std::runtime_error(err_msg);
}
std::vector<Complex> output;
std::transform(foutput.cbegin(), foutput.cend(), std::back_inserter(output),
[] (float f) { return Complex(f); }
);
return output;
}

85
sampling.h Normal file
View File

@ -0,0 +1,85 @@
/**
* @file sampling.h
* @author 9exa
* @brief Objects that resamples signal frames
* @version 0.1
* @date 2023-05-08
*
* @copyright Copyright (c) 2023
*
*/
#ifndef MENGA_SAMPLING
#define MENGA_SAMPLING
#include "miniaudio.h"
#include "mengumath.h"
#include "common.h"
#include <cstdint>
#include <vector>
namespace Mengu {
namespace dsp {
// Performs linear resampling without any filtering. Useful for resampling frequency domain directly
template<class T>
inline void linear_resample_no_filter(const T *input,
T *output,
uint32_t size,
float shift_factor,
T interp(T a, T b, float w) = &lerp,
T def_val = T()) {
float shift_rep = 1 / shift_factor;
// interpelate the positiobs in the input
if (shift_factor > 1.0f) {
// iterate from the back of the array to make the resample work in-place
for (uint32_t j = 0; j < size; j++) {
float inp_index = j * shift_rep;
uint32_t lower = uint32_t(inp_index);
float remainder = inp_index - lower;
output[j] = interp(input[lower], input[lower + 1], remainder);
}
}
else {
int j = size - 1;
int resample_start = (int) (size * shift_factor);
while (j >= resample_start) {
output[j] = def_val;
j -= 1;
}
while (j >= 0) {
int inp_index = j * shift_rep;
int lower = int(inp_index);
float remainder = inp_index - lower;
output[j] = interp(input[lower], input[lower + 1], remainder);
j -= 1;
}
}
}
// Just a wrapper around miniaudios resampler API, which uses linear resampling
// with a 4th order lowpass filter by default
class LinearResampler {
public:
LinearResampler(uint32_t nchannels, float stretch_factor);
~LinearResampler();
void set_stretch_factor(float stretch_factor);
std::vector<Complex> resample(const std::vector<Complex> &samples);
private:
ma_linear_resampler _resampler;
float _stretch_factor;
uint32_t _nchannels;
};
}
}
#endif

42
singletons.h Normal file
View File

@ -0,0 +1,42 @@
/*
Store large objects used by multiple dsp objects in one place (i.e. FFT Transforms)
*/
#ifndef MENGA_SINGLETONS
#define MENGA_SINGLETONS
#include <unordered_map>
#include "common.h"
#include "fft.h"
namespace Mengu {
namespace dsp {
struct Singletons {
private:
static Singletons *_singleton;
std::unordered_map<uint32_t, FFT> _ffts;
public:
static Singletons *get_singleton() {
if (_singleton = nullptr) {
_singleton = new Singletons();
}
return _singleton;
}
// getters are gaurenteed to return a valid object. If none exists one will be created
// ffts initialized to process predetermined length of signal
const FFT &get_fft(uint32_t size) {
if (_ffts.find(size) != _ffts.end()) {
_ffts.emplace(size, FFT(size));
}
return _ffts.at(size);
}
};
}
}
#endif

809
timestretcher.cpp Normal file
View File

@ -0,0 +1,809 @@
#include "timestretcher.h"
#include "common.h"
#include "correlation.h"
#include "effect.h"
#include "interpolation.h"
#include "linalg.h"
#include "mengumath.h"
#include "vecdeque.h"
#include <algorithm>
#include <array>
#include <cassert>
#include <cmath>
#include <complex>
#include <cstdint>
#include <iostream>
#include <numeric>
#include <vector>
using namespace Mengu;
using namespace dsp;
void TimeStretcher::set_stretch_factor(const float &scale) {
_stretch_factor = scale;
}
Effect::InputDomain TimeStretcher::get_input_domain() {
return InputDomain::Time;
}
std::vector<EffectPropDesc> TimeStretcher::get_property_descs() const {
return {
EffectPropDesc {
.type = EffectPropType::Slider,
.name = "stretch_factor",
.desc = "Scales the length of pushed signals by this amount",
.slider_data = {
.min_value = 0.5,
.max_value = 2,
.scale = Exp,
}
}
};
};
void TimeStretcher::set_property(uint32_t id, EffectPropPayload data) {
switch (id) {
case 0:
if (data.type == Slider) {
set_stretch_factor(data.value);
}
break;
};
}
EffectPropPayload TimeStretcher::get_property(uint32_t id) const {
switch (id) {
case 0:
return EffectPropPayload {
.type = Slider,
.value = _stretch_factor,
};
};
return EffectPropPayload {};
}
static float _wrap_phase(float phase) {
return fposmod(phase + MATH_PI, MATH_TAU) - MATH_PI;
}
// difference between 2 unwrapped phases. since, phases are preiodic, picks the closest one to the estimate
static float _phase_diff(float next, float prev, float est = 0.0f) {
return _wrap_phase(next - prev - est) + est;
}
PhaseVocoderTimeStretcher::PhaseVocoderTimeStretcher(bool preserve_formants) {
_preserve_formants = preserve_formants;
reset();
}
void PhaseVocoderTimeStretcher::push_signal(const Complex *input, const uint32_t &size) {
_raw_buffer.extend_back(input, size);
}
uint32_t PhaseVocoderTimeStretcher::pop_transformed_signal(Complex *output, const uint32_t &size) {
while (n_transformed_ready() < size && _raw_buffer.size() >= WindowSize) {
std::array<Complex, WindowSize> sample;
_raw_buffer.to_array(sample.data(), WindowSize);
_load_new_freq_window(sample);
std::array<Complex, WindowSize / 2> curr_freqs;
std::copy(_lpc.get_freq_spectrum().cbegin(), _lpc.get_freq_spectrum().cbegin() + WindowSize / 2, curr_freqs.begin());
float analysis_hop_sizef = SynthesisHopSize / _stretch_factor;
_stretched_sample_truncated += std::modf(analysis_hop_sizef, &analysis_hop_sizef);
uint32_t analysis_hop_size = (uint32_t) analysis_hop_sizef;
if (_stretched_sample_truncated) {
analysis_hop_size += 1;
_stretched_sample_truncated -= 1;
}
std::array<float, WindowSize / 2> amplitudes = _calc_scaled_magnitudes();
std::array<float, WindowSize / 2> new_phases = _calc_scaled_phases(curr_freqs, analysis_hop_size);
// std::array<Complex, WindowSize> new_samples = _calc_new_samples(amplitudes.data(), curr_raw_phases.data());
std::array<Complex, WindowSize> new_samples = _calc_new_samples(sample, amplitudes.data(), new_phases.data());
mix_and_extend(_transformed_buffer, new_samples, WindowSize - SynthesisHopSize, hann_window);
// _transformed_buffer.extend_back(new_samples.data(), WindowSize);
std::transform(new_phases.cbegin(), new_phases.cend(), _last_scaled_phases.begin(),
[] (float unwrapped) { return _wrap_phase(unwrapped); }
);
_raw_buffer.pop_front_many(nullptr, analysis_hop_size);
}
return _transformed_buffer.pop_front_many(output, MIN(n_transformed_ready(), size));
}
uint32_t PhaseVocoderTimeStretcher::n_transformed_ready() const {
return MAX(SynthesisHopSize, _transformed_buffer.size()) - (SynthesisHopSize);
}
void PhaseVocoderTimeStretcher::reset() {
_prev_raw_mag2s.fill(0.0f);
_prev_raw_phases.fill(0.0f);
_last_scaled_phases.fill(0.0);
_raw_buffer.resize(0);
_transformed_buffer.resize(SynthesisHopSize, 0);
_stretched_sample_truncated = 0.0;
}
// void PhaseVocoderTimeStretcher::set_stretch_factor(const float &stretch_factor) {
// _stretch_factor = stretch_factor;
// // _last_scaled_phases.fill(0.0f);
// };
std::array<float, PhaseVocoderTimeStretcher::WindowSize / 2> PhaseVocoderTimeStretcher::_calc_scaled_magnitudes() {
std::array<float, WindowSize / 2> mags;
if (_preserve_formants) {
const std::array<float, WindowSize> &envelope = _lpc.get_envelope();
const std::array<float, WindowSize> &residuals = _lpc.get_residuals();
for (uint32_t i = 0; i < WindowSize / 2; i++) {
const uint32_t stretched_ind = i * _stretch_factor;
if (stretched_ind < WindowSize / 2 && std::isfinite(residuals[i] * envelope[stretched_ind])) {
mags[i] = residuals[i] * envelope[stretched_ind];
}
else {
mags[i] = 0.0f;
}
}
}
else {
const std::array<Complex, WindowSize> &freqs = _lpc.get_freq_spectrum();
for (uint32_t i = 0; i < WindowSize / 2; i++) {
mags[i] = sqrt(std::norm(freqs[i]));
}
}
return mags;
}
std::array<float, PhaseVocoderTimeStretcher::WindowSize / 2> PhaseVocoderTimeStretcher::_calc_scaled_phases(
const std::array<Complex, WindowSize / 2> &curr_freqs,
const uint32_t hopsize) {
std::array<float, WindowSize / 2> curr_mag2s;
std::array<float, WindowSize / 2> curr_phases;
std::transform(curr_freqs.cbegin(), curr_freqs.cend(), curr_mag2s.begin(),
[] (Complex freq) { return std::norm(freq); }
);
std::transform(curr_freqs.cbegin(), curr_freqs.cend(), curr_phases.begin(),
[] (Complex freq) { return std::arg(freq); }
);
std::array<float, WindowSize / 2> phase_deltas;
float stretch_factor = (float) SynthesisHopSize / hopsize;
for (int i = 0; i < WindowSize / 2; i++) {
const float est = i * MATH_TAU * ((float) hopsize / WindowSize);// / stretch_factor;
phase_deltas[i] = _phase_diff(curr_phases[i], _prev_raw_phases[i], est);
}
std::array<float, WindowSize / 2> new_phases;
for (int i = 0; i < WindowSize / 2; i++) {
new_phases[i] = _last_scaled_phases[i] + (phase_deltas[i] * stretch_factor);
}
_replace_prev_freqs(curr_mag2s, curr_phases);
return new_phases;
}
std::array<Complex, PhaseVocoderTimeStretcher::WindowSize / 2> PhaseVocoderTimeStretcher::_load_new_freq_window(const std::array<Complex, WindowSize> &sample) {
_lpc.load_sample(sample.data());
std::array<Complex, WindowSize / 2> new_freqs;
for (uint32_t i = 0; i < WindowSize / 2; i++) {
new_freqs[i] = _lpc.get_freq_spectrum()[i];
}
return new_freqs;
}
void PhaseVocoderTimeStretcher::_replace_prev_freqs(const std::array<float, WindowSize / 2> &curr_mag2s, const std::array<float, WindowSize / 2> &curr_phases) {
std::copy(curr_mag2s.cbegin(), curr_mag2s.cend(), _prev_raw_mag2s.begin());
std::copy(curr_phases.cbegin(), curr_phases.cend(), _prev_raw_phases.begin());
}
std::array<Complex, PhaseVocoderTimeStretcher::WindowSize> PhaseVocoderTimeStretcher::_calc_new_samples(
const std::array<Complex, WindowSize> &raw_samples,
const float *amplitudes,
const float *phases) {
std::array<Complex, WindowSize> freqs = {0.0};
for (uint32_t i = 0; i < WindowSize / 2; i++) {
freqs[i] = 2.0f * std::polar(amplitudes[i], phases[i]);
}
std::array<Complex, WindowSize> new_samples;
_lpc.get_fft().inverse_transform(freqs.data(), new_samples.data());
// Make shifted as loud as raw samples
_loudness_norm.normalize(new_samples.data(), raw_samples.data(), new_samples.data());
return new_samples;
}
PhaseVocoderDoneRightTimeStretcher::PhaseVocoderDoneRightTimeStretcher(bool preserve_formants) : PhaseVocoderTimeStretcher(preserve_formants) {}
std::array<float, PhaseVocoderTimeStretcher::WindowSize / 2> PhaseVocoderDoneRightTimeStretcher::_calc_scaled_phases(
const std::array<Complex, WindowSize / 2> &curr_freqs,
const uint32_t hopsize) {
std::array<float, WindowSize / 2> curr_mag2s;
std::array<float, WindowSize / 2> curr_phases;
std::transform(curr_freqs.cbegin(), curr_freqs.cend(), curr_mag2s.begin(),
[] (Complex freq) { return std::norm(freq); }
);
std::transform(curr_freqs.cbegin(), curr_freqs.cend(), curr_phases.begin(),
[] (Complex freq) { return std::arg(freq); }
);
float stretch_factor = (float) SynthesisHopSize / hopsize;
std::array<float, WindowSize / 2> time_phase_deltas;
for (int i = 0; i < WindowSize / 2; i++) {
const float est = i * MATH_TAU * ((float) hopsize / WindowSize);
time_phase_deltas[i] = _phase_diff(curr_phases[i], _prev_raw_phases[i], est);
}
std::array<float, WindowSize / 2> freq_phase_deltas;
for (uint32_t i = 1; i < WindowSize / 2 - 1; i++) {
const float up_delta = _wrap_phase(curr_phases[i + 1] - curr_phases[i]);
const float down_delta = _wrap_phase(curr_phases[i] - curr_phases[i - 1]);
freq_phase_deltas[i] = 0.5f * (up_delta + down_delta);
}
freq_phase_deltas[0] = _wrap_phase(curr_phases[1] - curr_phases[0]);
freq_phase_deltas[WindowSize / 2 - 1] = _wrap_phase(curr_phases[WindowSize / 2 - 1] - curr_phases[WindowSize / 2 - 2]);
std::array<float, WindowSize / 2> new_phases = _propagate_phase_gradients (
time_phase_deltas,
freq_phase_deltas,
_last_scaled_phases,
_prev_raw_mag2s,
curr_mag2s,
stretch_factor,
1e-3f
);
// for (int i = 0; i < WindowSize / 2; i++) {
// new_phases[i] = _last_scaled_phases[i] + (time_phase_deltas[i] * stretch_factor);
// }
_replace_prev_freqs(curr_mag2s, curr_phases);
return new_phases;
}
std::array<float, PhaseVocoderTimeStretcher::WindowSize / 2> PhaseVocoderDoneRightTimeStretcher::_propagate_phase_gradients(
const std::array<float, WindowSize / 2> &phase_time_deltas,
const std::array<float, WindowSize / 2> &phase_freq_deltas,
const std::array<float, WindowSize / 2> &last_stretched_phases,
const std::array<float, WindowSize / 2> &prev_freq_mags,
const std::array<float, WindowSize / 2> &next_freq_mags,
const float stretch_factor,
const float tolerance) {
// sort indexes to the next frequencies based on the magnitude of that bin, in descending order
// also, bins with magnitude under the threshold do not propagate in the frequency domain
float max_mag = 0.0f;
for (uint32_t i = 0; i < WindowSize / 2; i++) {
max_mag = MAX(MAX(max_mag, prev_freq_mags[i]), next_freq_mags[i]);
}
const float abs_tol = max_mag * (tolerance * tolerance);
// Used to store indexs to each frequencing in the propagation_queue
struct FreqBin {
enum {
Prev,
Next,
};
uint8_t frame;
uint32_t bin;
};
const auto freq_bin_cmp = [&prev_freq_mags, &next_freq_mags] (FreqBin a, FreqBin b) {
float a_mag = (a.frame == FreqBin::Prev) ? prev_freq_mags[a.bin] : next_freq_mags[a.bin];
float b_mag = (b.frame == FreqBin::Prev) ? prev_freq_mags[b.bin] : next_freq_mags[b.bin];
return a_mag < b_mag;
};
// Max Heap of the frequency bins to propagate. (Listen I REALLY don't want allocate dynamically)
std::array<FreqBin, WindowSize * 3> propagation_queue;
uint32_t propagation_queue_size = WindowSize / 2;
for (uint32_t i = 0; i < WindowSize / 2; i++) {
propagation_queue[i] = FreqBin { .frame = FreqBin::Prev, .bin = i };
}
std::make_heap(propagation_queue.begin(), propagation_queue.begin() + propagation_queue_size, freq_bin_cmp);
// Set of frequency bins to propagate to
bool can_recieve_propagation[WindowSize / 2];
uint32_t n_can_recieve_propagation = WindowSize / 2;
for (uint32_t i = 0; i < WindowSize / 2; i++) {
if ((next_freq_mags[i] < abs_tol)) {
can_recieve_propagation[i] = false;
n_can_recieve_propagation -= 1;
}
else {
can_recieve_propagation[i] = true;
}
}
// perform propagation in all dimension
std::array<float, WindowSize / 2> new_phases {0};
while (n_can_recieve_propagation > 0) {
std::pop_heap(propagation_queue.begin(), propagation_queue.begin() + propagation_queue_size, freq_bin_cmp);
FreqBin next_bin = propagation_queue[propagation_queue_size - 1];
propagation_queue_size -= 1;
const uint32_t freq_ind = next_bin.bin;
if (next_bin.frame == FreqBin::Prev) {
if (can_recieve_propagation[freq_ind]) {
new_phases[freq_ind] = last_stretched_phases[freq_ind] + (phase_time_deltas[freq_ind] * stretch_factor);
//remove from set
can_recieve_propagation[freq_ind] = false;
n_can_recieve_propagation -= 1;
// push to the heap
propagation_queue[propagation_queue_size] = FreqBin {.frame = FreqBin::Next, .bin = freq_ind};
propagation_queue_size += 1;
std::push_heap(propagation_queue.begin(), propagation_queue.begin() + propagation_queue_size, freq_bin_cmp);
}
}
else {
if (freq_ind > 0 && can_recieve_propagation[freq_ind - 1]) {
const uint32_t freq_down = freq_ind - 1;
new_phases[freq_down] = new_phases[freq_ind] - (0.5 * (phase_freq_deltas[freq_down] + phase_freq_deltas[freq_ind]) * stretch_factor);
//remove from set
can_recieve_propagation[freq_down] = false;
n_can_recieve_propagation -= 1;
// push to the heap
propagation_queue[propagation_queue_size] = FreqBin {.frame = FreqBin::Next, .bin = freq_down};
propagation_queue_size += 1;
std::push_heap(propagation_queue.begin(), propagation_queue.begin() + propagation_queue_size, freq_bin_cmp);
}
if ((freq_ind < WindowSize / 2 - 1) && can_recieve_propagation[freq_ind + 1]) {
const uint32_t freq_up = freq_ind + 1;
new_phases[freq_up] = new_phases[freq_ind] + (0.5 * (phase_freq_deltas[freq_up] + phase_freq_deltas[freq_ind]) * stretch_factor);
//remove from set
can_recieve_propagation[freq_up] = false;
n_can_recieve_propagation -= 1;
// push to the heap
propagation_queue[propagation_queue_size] = FreqBin {.frame = FreqBin::Next, .bin = freq_up};
propagation_queue_size += 1;
std::push_heap(propagation_queue.begin(), propagation_queue.begin() + propagation_queue_size, freq_bin_cmp);
}
}
}
// std::array<float, WindowSize / 2> new_phases {0};
// auto freq_ind_cmp = [next_freq_mags] (uint32_t a, uint32_t b) { return next_freq_mags[a] > next_freq_mags[b]; };
// std::array<uint32_t, WindowSize / 2> freq_inds;
// for (uint32_t i = 0; i < WindowSize / 2; i++) { freq_inds[i] = i;}
// std::sort(freq_inds.begin(), freq_inds.end(), freq_ind_cmp);
// std::array<float, WindowSize / 2> prop_source_mag {0.0f};
// for (uint32_t i = 0; i < WindowSize / 2; i++) {
// new_phases[i] = last_stretched_phases[i] + phase_time_deltas[i] * stretch_factor;
// prop_source_mag[i] = prev_freq_mags[i];
// }
// for (const uint32_t freq_ind: freq_inds) {
// if (freq_ind > 0 && prop_source_mag[freq_ind - 1] < next_freq_mags[freq_ind]) {
// const uint32_t freq_down = freq_ind - 1;
// new_phases[freq_down] = new_phases[freq_ind] - 0.5 * stretch_factor * (phase_freq_deltas[freq_down] + phase_freq_deltas[freq_ind]);
// prop_source_mag[freq_down] = next_freq_mags[freq_ind];
// }
// if (freq_ind < WindowSize / 2 && prop_source_mag[freq_ind + 1] < next_freq_mags[freq_ind]) {
// const uint32_t freq_up = freq_ind + 1;
// new_phases[freq_up] = new_phases[freq_ind] - 0.5 * stretch_factor * (phase_freq_deltas[freq_up] + phase_freq_deltas[freq_ind]);
// prop_source_mag[freq_up] = next_freq_mags[freq_ind];
// }
// prop_source_mag[freq_ind] = 2e10;
// }
return new_phases;
}
OLATimeStretcher::OLATimeStretcher(uint32_t w_size) {
_window_size = w_size;
_overlap = _window_size / 5;
_selection_window = _window_size / 2;
_transformed_buffer.resize(_window_size);
}
template<class T>
static void mix_into_extend_by_pointer(const Complex *new_data,
T &output,
const uint32_t &window_size,
const uint32_t &overlap_size) {
uint32_t i = 0;
for (; i < overlap_size; i++) {
const float w = hann_window((float) i / overlap_size);
const uint32_t output_ind = output.size() - overlap_size + i;
output[output_ind] = lerp(new_data[i], output[output_ind], w);
}
for (; i < window_size; i++) {
output.push_back(new_data[i]);
}
}
void OLATimeStretcher::push_signal(const Complex *input, const uint32_t &size) {
_raw_buffer.extend_back(input, size);
}
uint32_t OLATimeStretcher::pop_transformed_signal(Complex *output, const uint32_t &size) {
// Do stretchy
// Theoretical interval between samples per process
const uint32_t sample_skip = (_window_size - _overlap) / _stretch_factor;
// based on the example given by https://www.surina.net/article/time-and-pitch-scaling.html
const uint32_t length_for_process = MAX(sample_skip, _window_size + _selection_window);
while (_raw_buffer.size() > length_for_process) {
std::vector<Complex> new_data(_window_size + _selection_window);
_raw_buffer.to_array(new_data.data(), _window_size + _selection_window);
std::vector<Complex> prev_tail(_overlap);
_transformed_buffer.pop_back_many(prev_tail.data(), _overlap);
// find best start for overlap
uint32_t overlap_ind = find_max_correlation_quad(prev_tail.data(), new_data.data(), _overlap, _selection_window);
mix_into_extend_by_pointer(new_data.data() + overlap_ind, prev_tail, _window_size, _overlap);
for (auto v: prev_tail) {
_transformed_buffer.push_back(v);
}
_raw_buffer.pop_front_many(nullptr, sample_skip);
}
uint32_t n = _transformed_buffer.pop_front_many(output, MIN(size, n_transformed_ready()));
for (uint32_t i = n; i < size; i++ ) {
output[i] = 0;
}
return n;
}
uint32_t OLATimeStretcher::n_transformed_ready() const {
return _transformed_buffer.size() - _overlap;
}
void OLATimeStretcher::reset() {
_raw_buffer.resize(0);
_transformed_buffer.resize(_window_size, 0);
}
WSOLATimeStretcher::WSOLATimeStretcher() {
// _transformed_buffer.resize(MaxBackWindowOverlap);
}
void WSOLATimeStretcher::push_signal(const Complex *input, const uint32_t &size) {
_raw_buffer.extend_back(input, size);
// perform the stretchy
// if (_raw_buffer.size() > SampleProcSize) {
// std::array<Complex, SampleProcSize> samples;
// _raw_buffer.to_array(samples.data(), SampleProcSize);
// // do the stretchy
// uint32_t frames_used = _stretch_sample_and_add(samples.data());
// // delete everything used
// _raw_buffer.pop_front_many(nullptr, frames_used);
// _last_overlap_start -= frames_used;
// _next_overlap_start -= frames_used;
// }
}
uint32_t WSOLATimeStretcher::pop_transformed_signal(Complex *output, const uint32_t &size) {
// lazily perform the stretchy
while ((_raw_buffer.size() > SampleProcSize) && (n_transformed_ready() < size)) {
std::array<Complex, SampleProcSize> samples;
_raw_buffer.to_array(samples.data(), SampleProcSize);
// do the stretchy
uint32_t frames_used = _stretch_sample_and_add(samples.data());
// delete everything used
_raw_buffer.pop_front_many(nullptr, frames_used);
_last_overlap_start -= frames_used;
_next_overlap_start -= frames_used;
}
uint32_t n = _transformed_buffer.pop_front_many(output, MIN(size, n_transformed_ready()));
for (uint32_t i = n; i < size; i++) {
output[i] = 0;
}
return n;
}
uint32_t WSOLATimeStretcher::n_transformed_ready() const {
return _transformed_buffer.size();
}
void WSOLATimeStretcher::reset() {
_raw_buffer.resize(0);
_transformed_buffer.resize(0);
}
uint32_t WSOLATimeStretcher::_stretch_sample_and_add(const Complex *sample) {
// base overlap
const uint32_t overlap_size = WindowSize / 4;
// search forward for better overlap point
const uint32_t search_window = WindowSize / 5;
const uint32_t flat_duration = WindowSize - 2 * overlap_size;
// consider te amount of frames skiped through truncation
double sample_skipd;
double dropped_per_window = std::modf((WindowSize - overlap_size) / _stretch_factor, &sample_skipd);
const uint32_t sample_skip = sample_skipd;
while (_next_overlap_start + WindowSize < SampleProcSize) {
// Find insertion that best fits the new sample
uint32_t prev_not_overlapped = find_max_correlation(
sample + _next_overlap_start,
sample + _last_overlap_start,
overlap_size,
search_window
);
uint32_t actual_last_overlap = _last_overlap_start + prev_not_overlapped;
const uint32_t actually_overlapped = overlap_size - prev_not_overlapped;
Complex overlap_buffer[SampleProcSize / 2];
overlap_add(
sample + actual_last_overlap,
sample + _next_overlap_start,
overlap_buffer,
actually_overlapped,
hamming_window
);
// append new data
_transformed_buffer.extend_back(sample + _last_overlap_start, prev_not_overlapped);
_transformed_buffer.extend_back(overlap_buffer, actually_overlapped);
// take 2 * prev_not_overlapped from both sides of the flat duration
_transformed_buffer.extend_back(sample + _next_overlap_start + actually_overlapped, flat_duration);
// set for next cycle
_last_overlap_start = _next_overlap_start + actually_overlapped + flat_duration;
_next_overlap_start = _next_overlap_start + sample_skip;
_stretched_sample_truncated += dropped_per_window;
if (_stretched_sample_truncated > 1.0) {
_next_overlap_start += 1;
_stretched_sample_truncated -= 1.0;
}
}
return MIN(_last_overlap_start, _next_overlap_start);
}
PSOLATimeStretcher::PSOLATimeStretcher() {
_transformed_buffer.resize(MaxBackWindowOverlap);
}
void PSOLATimeStretcher::push_signal(const Complex *input, const uint32_t &size) {
_raw_buffer.extend_back(input, size);
}
uint32_t PSOLATimeStretcher::pop_transformed_signal(Complex *output, const uint32_t &size) {
// lazily perform the stretchy
while ((_raw_buffer.size() > SampleProcSize) && (size > n_transformed_ready())) {
std::array<Complex, SampleProcSize> samples;
std::array<Complex, SampleProcSize> windowed_samples;
_raw_buffer.to_array(samples.data(), SampleProcSize);
_raw_buffer.to_array(windowed_samples.data(), SampleProcSize);
window_ends(windowed_samples.data(), SampleProcSize, SampleProcSize / 10, hann_window);
int est_freq = _est_fund_frequency(windowed_samples.data());
uint32_t est_period = (1.0 / est_freq) * SampleProcSize / 2;
std::vector<uint32_t> est_peaks = _find_upcoming_peaks(samples.data(), est_period);
_stretch_peaks_and_add(samples.data(), est_peaks);
// delete everything used
_raw_buffer.pop_front_many(nullptr, est_peaks.back() + 1); // +1 because est_peaks are the INDEX of the peaks, not the num of frames used
}
uint32_t n = _transformed_buffer.pop_front_many(output, MIN(size, n_transformed_ready()));
for (uint32_t i = n; i < size; i++) {
output[i] = 0;
}
return n;
}
uint32_t PSOLATimeStretcher::n_transformed_ready() const {
return _transformed_buffer.size() - MaxBackWindowOverlap;
}
void PSOLATimeStretcher::reset() {
_transformed_buffer.resize(MaxBackWindowOverlap, 0);
_raw_buffer.resize(0);
}
int PSOLATimeStretcher::_est_fund_frequency(const Complex *samples) {
// Todo: move all this to a PitchDetecter object or something
_lpc.load_sample(samples);
const std::array<float, SampleProcSize> residuals = _lpc.get_residuals();
// find candidate from residual peaks (only use positive half of the spectrum)
std::vector<float> srhs = calc_srhs(residuals.data(), residuals.size() / 2, MinFreqInd, MaxFreqInd, 10);
float max_srhs = -10e32;
int max_pitch_ind = MaxFreqInd;
for (int i = 0; i < srhs.size(); i++) {
if (srhs[i] > max_srhs) {
max_pitch_ind = MinFreqInd + i;
max_srhs = srhs[i];
}
}
return max_pitch_ind;
}
std::vector<uint32_t> PSOLATimeStretcher::_find_upcoming_peaks(const Complex *samples, const uint32_t est_period) {
// assume that peaks are around est_period apart, but give some sllack as pitches change slightly
const uint32_t search_start = 0.8 * est_period;
const uint32_t search_end = 1.2 * est_period;
std::vector<uint32_t> peaks;
uint32_t last_peak = 0;
while ((last_peak + est_period) < SampleProcSize) {
uint32_t peak = last_peak + search_start;
float peak_size = 0.0f;
for (uint32_t i = last_peak + search_start; i < MIN(SampleProcSize, last_peak + search_end); i++) {
if (samples[i].real() > peak_size) {
peak_size = samples[i].real();
peak = i;
}
}
peaks.push_back(peak);
last_peak = peak;
}
return peaks;
}
// // overlap and extend without applying a window function first
// // overlap_size may be larger than new_data.size(), in which case only new_data.size() is added. the offset is the same
template<class T, class NewT>
static void _mix_and_extend_no_window(T &array, const NewT new_data, const uint32_t &overlap_size) {
uint32_t i = 0;
for(; i < MIN(overlap_size, new_data.size()); i++) {
uint32_t ind = array.size() - overlap_size + i;
array[ind] = array[ind] + new_data[i];
}
for (; i < new_data.size(); i++) {
array.push_back(new_data[i]);
}
}
void PSOLATimeStretcher::_stretch_peaks_and_add(const Complex *samples, const std::vector<uint32_t> &est_peaks) {
// accumulate by collecting the left halves and right halves if each peak sepeartly
std::vector<std::vector<Complex>> right_windows;
std::vector<std::vector<Complex>> left_windows;
uint32_t last_peak = 0;
for (uint32_t next_peak: est_peaks) {
std::vector<Complex> left_window;
std::vector<Complex> right_window;
for (uint32_t i = last_peak; i < next_peak; i++) {
float w = (float) (i - last_peak) / (next_peak - last_peak);
w = hann_window(w);
left_window.emplace_back(w * samples[i]);
right_window.emplace_back((1.0f - w) * samples[i]);
}
left_windows.push_back(std::move(left_window));
right_windows.push_back(std::move(right_window));
last_peak = next_peak;
}
// arrange the peaks together
last_peak = 0;
for (uint32_t i = 0; i < left_windows.size(); i++) {
const std::vector<Complex> left_window = std::move(left_windows[i]);
const std::vector<Complex> right_window = std::move(right_windows[i]);
// use overlapsize of previous period to make the right window continuous with the previous left window
_mix_and_extend_no_window(_transformed_buffer, right_window, _next_right_window_overlap);
// calculate the frames overlap taking the truncated part of previous processes into account
int curr_period = est_peaks[i] - last_peak;
double overlap_sized;
_stretched_sample_truncated += std::abs(std::modf((2.0 - _stretch_factor) * curr_period, &overlap_sized));
uint32_t overlap_size = std::abs(overlap_sized);
if (_stretched_sample_truncated > 1.0f) {
overlap_size += 1;
_stretched_sample_truncated -= 1.0f;
}
if (_stretch_factor >= 1.0f) {
if (_stretch_factor > 2.0f) {
// pad the gap between each side with 0s
uint32_t padding_size = std::move(overlap_size);
for (uint32_t j = 0; j < padding_size; j++) {
_transformed_buffer.push_back(0);
}
_transformed_buffer.extend_back(left_window.data(), left_window.size());
}
else {
// some parts are overlapped
// mix it with the left window of the next
_mix_and_extend_no_window(_transformed_buffer, left_window, overlap_size);
}
// the end of the _transform buffer is the end of the left window, which is equivalenting to setting both of these values to 0
_next_right_window_overlap = 0;
}
else {
// since the left peak ends before the right peak ends, we have to overlap both the right and left peaks into the transformed data
_mix_and_extend_no_window(_transformed_buffer, left_window, overlap_size);
_next_right_window_overlap = overlap_size - curr_period;
}
last_peak = est_peaks[i];
}
}

226
timestretcher.h Normal file
View File

@ -0,0 +1,226 @@
#ifndef MENGA_TIME_STRETCHER
#define MENGA_TIME_STRETCHER
#include "common.h"
#include "correlation.h"
#include "effect.h"
#include "fft.h"
#include "loudness.h"
#include "vecdeque.h"
#include <array>
#include <cstdint>
#include <vector>
namespace Mengu {
namespace dsp {
class TimeStretcher: public Effect {
public:
virtual ~TimeStretcher() {}
virtual InputDomain get_input_domain() override;
virtual void set_stretch_factor(const float &scale);
virtual std::vector<EffectPropDesc> get_property_descs() const override;
virtual void set_property(uint32_t id, EffectPropPayload data) override;
virtual EffectPropPayload get_property(uint32_t id) const override;
protected:
float _stretch_factor = 1.0f;
// keep track of resampled size error due to rounding errors
double _stretched_sample_truncated = 0.0;
};
// Classic timeshifter be scale the phases of frequency bins in the time dimension
class PhaseVocoderTimeStretcher: public TimeStretcher {
public:
PhaseVocoderTimeStretcher(bool _preserve_formants = false);
virtual void push_signal(const Complex *input, const uint32_t &size) override;
virtual uint32_t pop_transformed_signal(Complex *output, const uint32_t &size) override;
virtual uint32_t n_transformed_ready() const override;
// virtual void set_stretch_factor(const float &stretch_factor) override;
virtual void reset() override;
protected:
static constexpr uint32_t WindowSize = 1 << 9;
static constexpr uint32_t SynthesisHopSize = 400;
static constexpr uint32_t NStoredWindows = 2;
std::array<float, WindowSize / 2> _prev_raw_mag2s;
std::array<float, WindowSize / 2> _prev_raw_phases;
// phases of the last transformed samples
std::array<float, WindowSize / 2> _last_scaled_phases;
VecDeque<Complex> _raw_buffer;
VecDeque<Complex> _transformed_buffer;
virtual std::array<float, WindowSize / 2> _calc_scaled_magnitudes();
// expects time_deltas to be scaled
virtual std::array<float, WindowSize / 2> _calc_scaled_phases(const std::array<Complex, WindowSize / 2> &curr_freqs, const uint32_t hopsize);
std::array<Complex, WindowSize / 2> _load_new_freq_window(const std::array<Complex, WindowSize> &sample);
void _replace_prev_freqs(const std::array<float, WindowSize / 2> &curr_mags, const std::array<float, WindowSize / 2> &curr_phases);
private:
// FFT _fft;
// for finding the envelope and doing fft
LPC<WindowSize, 50> _lpc;
bool _preserve_formants = false;
// Used to make sure the percieved loudness of the sample is preserved
LoudnessNormalizer<Complex, WindowSize, 1> _loudness_norm;
std::array<Complex, WindowSize> _calc_new_samples(const std::array<Complex, WindowSize> &raw_samples, const float *amplitudes, const float *phase_deltas);
};
// 'Phase vocoder done right' implementation
class PhaseVocoderDoneRightTimeStretcher: public PhaseVocoderTimeStretcher {
public:
PhaseVocoderDoneRightTimeStretcher(bool _preserve_formants = false);
protected:
virtual std::array<float, WindowSize / 2> _calc_scaled_phases(const std::array<Complex, WindowSize / 2> &curr_phases, const uint32_t hopsize) override;
private:
// Phase propagation algorithm as described in the paper
std::array<float, WindowSize / 2> _propagate_phase_gradients(const std::array<float, WindowSize / 2> &phase_time_deltas,
const std::array<float, WindowSize / 2> &phase_freq_deltas,
const std::array<float, WindowSize / 2> &last_stretched_phases,
const std::array<float, WindowSize / 2> &prev_mag2s,
const std::array<float, WindowSize / 2> &next_mag2s,
const float stretch_factor,
const float tolerance = 1e-5f);
};
// Syncronised OverLap and Add time stretcher with fixed window size
class OLATimeStretcher: public TimeStretcher {
public:
OLATimeStretcher(uint32_t w_size);
static inline const float MinScale = 0.05;
// std::vector<Complex> stretch_and_overlap_window2(const Complex *input);
virtual void push_signal(const Complex *input, const uint32_t &size) override;
virtual uint32_t pop_transformed_signal(Complex *output, const uint32_t &size) override;
virtual uint32_t n_transformed_ready() const override;
virtual void reset() override;
private:
uint32_t _window_size;
uint32_t _overlap;
uint32_t _selection_window;
// uint32_t _sample_skip;
float _desired_extension = 0.0f;
VecDeque<Complex> _raw_buffer;
VecDeque<Complex> _transformed_buffer;
};
// timestrech where extensions are added to best match wave form with. Fixed window size. not fixed soutput size (may be slightly longer)
class WSOLATimeStretcher: public TimeStretcher {
public:
WSOLATimeStretcher();
virtual void push_signal(const Complex *input, const uint32_t &size) override;
virtual uint32_t pop_transformed_signal(Complex *output, const uint32_t &size) override;
virtual uint32_t n_transformed_ready() const override;
virtual void reset() override;
private:
VecDeque<Complex> _raw_buffer;
VecDeque<Complex> _transformed_buffer;
// length of the the input buffer must be before transforming and size of arrays in intermediate calculations. Should catch up to 1000hz
static constexpr uint32_t SampleProcSize = 1 << 11;
// length of a window
static constexpr uint32_t WindowSize = 1 << 9;
// stretches the sample, and adds it to the transform buffer, tje position of each window is based on the autocorrelation
// returns how many frames were used and can be discarded
uint32_t _stretch_sample_and_add(const Complex *sample);
// Basically, the beginning of the next overlap(_left_tail_start) can be before the beggining of the last overlap(_right_tail_start)
// and the size of the overlap can change on each process. So store both.
uint32_t _last_overlap_start = 0;
uint32_t _next_overlap_start = 0;
};
// Formant preserving, pitch changing time stretch
// Inherits a lot of algorithmic functionality from SOLA
class PSOLATimeStretcher: public TimeStretcher {
public:
PSOLATimeStretcher();
virtual void push_signal(const Complex *input, const uint32_t &size) override;
virtual uint32_t pop_transformed_signal(Complex *output, const uint32_t &size) override;
virtual uint32_t n_transformed_ready() const override;
virtual void reset() override;
private:
VecDeque<Complex> _raw_buffer;
VecDeque<Complex> _transformed_buffer;
static constexpr uint32_t SampleProcSize = 1 << 11;
// basically how accurate the reconstruction will be. Assumes ~44kHz input, good for operation in up to 16kHz
static constexpr uint32_t LPCSize = (uint32_t) 1.25 * 16;
static constexpr uint32_t MinFreqHz = 50;
static constexpr uint32_t MaxFreqHz = 800;
static constexpr uint32_t InputSampleRate = 44100;
static constexpr uint32_t MinFreqInd = MAX(MinFreqHz * SampleProcSize / InputSampleRate, 2);
static constexpr uint32_t MaxFreqInd = MAX(MaxFreqHz * SampleProcSize / InputSampleRate, MinFreqInd * 2);
// the amount of frames in _transformed_buffer that need to remain in case of overlapping future samples
static constexpr uint32_t MaxBackWindowOverlap = (1.0f / MinFreqHz * InputSampleRate) + 1;
// for finding fundemental frequency
// FFT _fft;
LPC<SampleProcSize, LPCSize> _lpc;
// returns the index of the fundemental frequency (pitch) in an fft of samples
int _est_fund_frequency(const Complex *sample);
// store the last estimated pitches. the median will be used
VecDeque<int> _last_periods;
//used to meld windows in the same sample of different length (peaks are not uniformly spaced)
uint32_t _next_right_window_overlap = 0;
// estimate the peaks in the upcoming sample
std::vector<uint32_t> _find_upcoming_peaks(const Complex *samples, const uint32_t est_period);
// stretches the sample, and adds it to the transform buffer;
void _stretch_peaks_and_add(const Complex *samples, const std::vector<uint32_t> &est_peaks);
};
}
}
#endif

243
vecdeque.h Normal file
View File

@ -0,0 +1,243 @@
/**
* @file vecdeque.h
* @author 9exa
* @brief An dynamically-sized array where elements can be added on either end. Useful for queues
* Supposed to be like Rusts VecDeque, implemented with a ring_buffer
* @version 0.1
* @date 2023-04-22
*
* @copyright Copyright (c) 2023
*
*/
#ifndef MENGA_VECDEQUE
#define MENGA_VECDEQUE
#include <algorithm>
#include <cstdint>
#include <iostream>
#include "mengumath.h"
namespace Mengu {
template<class T>
class VecDeque {
private:
T *_data = nullptr;
uint32_t _front = 0;
// uint32_t _end = 0;
uint32_t _size = 0;
uint32_t _capacity = 0;
public:
VecDeque() {}
VecDeque(const VecDeque &from) {
_front = from._front;
_size = from._size;
_capacity = from._capacity;
if (_capacity > 0) {
_data = new T[_capacity];
for (uint32_t i = 0; i < _capacity; i++) {
_data[i] = from._data[i];
}
}
}
~VecDeque() {
if (_data != nullptr) {
delete[] _data;
}
}
inline uint32_t size() const {
return _size;
}
inline uint32_t capacity() const {
return _capacity;
}
inline const T *data() const {
return _data;
}
void reserve(const uint32_t &new_cap) {
if (_capacity < new_cap) {
T *new_data = new T[new_cap];
// copy and initialise new array
uint32_t i = 0;
if (_data != nullptr) {
for (; i < _size; i++) {
new_data[i] = std::move(_data[(_front + i) % _capacity]);
}
delete[] _data;
_front = 0;
}
_data = new_data;
_capacity = new_cap;
if (_size > _capacity) {
resize(_capacity);
}
}
}
void resize(const uint32_t &new_size) {
if (_capacity < new_size) {
uint32_t new_cap = MAX(_capacity, 1);
while (new_cap < new_size) {
new_cap = new_cap << 1; // if you leave out the new_cap = the optimizer just skips this loop.
// Which makes this infinite loop bug hard to spot
}
reserve(new_cap);
}
if (new_size > _size) {
// expand from back
if ( _size > 0) {
for (uint32_t i = _size; i < new_size; i++) {
_data[(_front + i) % _capacity] = _data[(_front + _size - 1) % _capacity];
}
}
else {
for (uint32_t i = _size; i < new_size; i++) {
_data[(_front + i) % _capacity] = T();
}
}
}
_size = new_size;
}
void resize(const uint32_t &new_size, const T &x) {
resize(new_size);
for (uint32_t i = 0; i < new_size; i++) {
_data[(_front + i) % _capacity] = x;
}
}
inline void push_back(const T &x) {
resize(_size + 1);
_data[(_front + _size - 1) % _capacity] = x;
}
inline void push_front(const T &x) {
resize(_size + 1);
_front = (_front == 0) ? _capacity - 1 : _front - 1;
_data[_front] = x;
}
// adds elements to the end of the queue
inline void extend_back(const T *array, const uint32_t n) {
uint32_t old_size = _size;
resize(_size + n);
for (uint32_t i = 0; i < n; i++) {
_data[(_front + old_size + i) % _capacity] = array[i];
}
}
// moves at most n elements from the front of the queue to the array output. outputs memory must be validated elsewhere
inline uint32_t pop_front_many(T *output, uint32_t n) {
n = MIN(n, _size);
if (output != nullptr) {
for (uint32_t i = 0; i < n; i++) {
output[i] = _data[(_front + i) % _capacity];
}
}
if (n != 0){
_front = (_front + n) % _capacity;
resize(_size - n);
}
return n;
}
inline uint32_t pop_back_many(T *output, uint32_t n) {
n = MIN(n, _size);
if (output != nullptr) {
for (uint32_t i = 0; i < n; i++) {
output[i] = _data[(_front + _size - n + i) % _capacity];
}
}
resize(_size - n);
return n;
}
void make_contiguous() {
if (_size > 0 && _front > 0) {
T *new_data = new T[_capacity];
for (uint32_t i = 0; i < _size; i++) {
new_data[i] = std::move(_data[(_front + i) % _capacity]);
}
delete[] _data;
_data = new_data;
}
}
//// Operators
inline T &operator[](int i) {
if (_size == 0) {
resize(1);
}
return _data[(_front + posmod(i, _size)) % _capacity];
}
inline const T &operator[](int i) const {
if (_size == 0) {
resize(1);
}
return _data[(_front + posmod(i, _size)) % _capacity];
}
VecDeque &operator=(const VecDeque &from) {
resize(from._size);
_front = 0;
for (uint32_t i = 0; i < _size; i++) {
_data[i] = from._data[(from._front + _size % from._capacity)];
}
return *this;
};
// converts the first 'size' items into a contiguous array. -1 does the whole queue
uint32_t to_array(T *out, int size = -1) {
if (size == -1) {
size = _size;
}
size = MIN(size, _size);
for (uint32_t i = 0; i < size; i++) {
out[i] = _data[(i + _front) % _capacity];
}
return size;
}
// writes the last 'size' items into a contiguous array
uint32_t to_array_back(T *out, int size) {
size = MIN(size, _size);
for (uint32_t i = 0; i < size; i++) {
out[i] = _data[(i + _front + _size - size) % _capacity];
}
return size;
}
};
}
#endif