circuitpython-ulab/docs/ulab-manual.ipynb
2020-02-26 18:05:49 +01:00

4772 lines
150 KiB
Text

{
"cells": [
{
"cell_type": "code",
"execution_count": 1,
"metadata": {
"ExecuteTime": {
"end_time": "2019-11-04T18:52:11.229037Z",
"start_time": "2019-11-04T18:52:10.081797Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Populating the interactive namespace from numpy and matplotlib\n"
]
}
],
"source": [
"%pylab inline"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {
"ExecuteTime": {
"end_time": "2020-02-26T17:04:16.562607Z",
"start_time": "2020-02-26T17:04:16.502531Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Overwriting manual/source/conf.py\n"
]
}
],
"source": [
"%%writefile manual/source/conf.py\n",
"# Configuration file for the Sphinx documentation builder.\n",
"#\n",
"# This file only contains a selection of the most common options. For a full\n",
"# list see the documentation:\n",
"# http://www.sphinx-doc.org/en/master/config\n",
"\n",
"# -- Path setup --------------------------------------------------------------\n",
"\n",
"# If extensions (or modules to document with autodoc) are in another directory,\n",
"# add these directories to sys.path here. If the directory is relative to the\n",
"# documentation root, use os.path.abspath to make it absolute, like shown here.\n",
"#\n",
"# import os\n",
"# import sys\n",
"# sys.path.insert(0, os.path.abspath('.'))\n",
"\n",
"\n",
"# -- Project information -----------------------------------------------------\n",
"\n",
"project = 'micropython-ulab'\n",
"copyright = '2019, Zoltán Vörös'\n",
"author = 'Zoltán Vörös'\n",
"\n",
"# The full version, including alpha/beta/rc tags\n",
"release = '0.32.2'\n",
"\n",
"\n",
"# -- General configuration ---------------------------------------------------\n",
"\n",
"# Add any Sphinx extension module names here, as strings. They can be\n",
"# extensions coming with Sphinx (named 'sphinx.ext.*') or your custom\n",
"# ones.\n",
"extensions = [\n",
"]\n",
"\n",
"# Add any paths that contain templates here, relative to this directory.\n",
"templates_path = ['_templates']\n",
"\n",
"# List of patterns, relative to source directory, that match files and\n",
"# directories to ignore when looking for source files.\n",
"# This pattern also affects html_static_path and html_extra_path.\n",
"exclude_patterns = []\n",
"\n",
"\n",
"# -- Options for HTML output -------------------------------------------------\n",
"\n",
"# The theme to use for HTML and HTML Help pages. See the documentation for\n",
"# a list of builtin themes.\n",
"#\n",
"html_theme = 'sphinx_rtd_theme'\n",
"\n",
"# Add any paths that contain custom static files (such as style sheets) here,\n",
"# relative to this directory. They are copied after the builtin static files,\n",
"# so a file named \"default.css\" will overwrite the builtin \"default.css\".\n",
"html_static_path = ['_static']\n",
"\n",
"master_doc = 'index'\n",
"\n",
"author=u'Zoltán Vörös'\n",
"copyright=author\n",
"language='en'\n",
"\n",
"latex_documents = [\n",
"(master_doc, 'ulab-manual.tex', 'Micropython ulab documentation', \n",
"'Zoltán Vörös', 'manual'),\n",
"]\n",
" "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Notebook conversion"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {
"ExecuteTime": {
"end_time": "2020-02-26T17:04:49.527515Z",
"start_time": "2020-02-26T17:04:21.416456Z"
}
},
"outputs": [],
"source": [
"import nbformat as nb\n",
"import nbformat.v4.nbbase as nb4\n",
"from nbconvert import RSTExporter\n",
"\n",
"def convert_notebook(node, fn):\n",
" (rst, resources) = rstexporter.from_notebook_node(notebook)\n",
" with open(fn, 'w') as fout:\n",
" fout.write(rst)\n",
" \n",
"rstexporter = RSTExporter()\n",
"rstexporter.template_file = './templates/manual.tpl'\n",
"\n",
"source = nb.read('ulab-manual.ipynb', nb.NO_CONVERT)\n",
"append_cell = False\n",
"\n",
"notebook = nb4.new_notebook()\n",
"for j, cell in enumerate(source['cells']):\n",
" if cell['cell_type'] == 'markdown':\n",
" # skip everything before Introduction\n",
" if cell['source'].split('\\n')[0].startswith('# Introduction'):\n",
" append_cell = True\n",
" if append_cell:\n",
" notebook.cells.append(cell)\n",
" \n",
"convert_notebook(notebook,'./manual/source/ulab.rst')"
]
},
{
"cell_type": "code",
"execution_count": 516,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-19T19:46:35.995470Z",
"start_time": "2019-10-19T19:46:35.989079Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Overwriting ./templates/manual.tpl\n"
]
}
],
"source": [
"%%writefile ./templates/manual.tpl\n",
"\n",
"{%- extends 'display_priority.tpl' -%}\n",
"\n",
"\n",
"{% block in_prompt %}\n",
"{% endblock in_prompt %}\n",
"\n",
"{% block output_prompt %}\n",
"{% endblock output_prompt %}\n",
"\n",
"{% block input scoped%}\n",
"\n",
"{%- if cell.source.split('\\n')[0].startswith('%%micropython') -%}\n",
".. code::\n",
" \n",
"{{ '\\n'.join(['# code to be run in micropython'] + cell.source.strip().split('\\n')[1:]) | indent}}\n",
"\n",
"{%- else -%}\n",
".. code::\n",
"\n",
"{{ '\\n'.join(['# code to be run in CPython\\n'] + cell.source.strip().split('\\n')) | indent}}\n",
"{%- endif -%}\n",
"{% endblock input %}\n",
"\n",
"{% block error %}\n",
"::\n",
"\n",
"{{ super() }}\n",
"{% endblock error %}\n",
"\n",
"{% block traceback_line %}\n",
"{{ line | indent | strip_ansi }}\n",
"{% endblock traceback_line %}\n",
"\n",
"{% block execute_result %}\n",
"{% block data_priority scoped %}\n",
"{{ super() }}\n",
"{% endblock %}\n",
"{% endblock execute_result %}\n",
"\n",
"{% block stream %}\n",
".. parsed-literal::\n",
"\n",
"{{ output.text | indent }}\n",
"{% endblock stream %}\n",
"\n",
"{% block data_svg %}\n",
".. image:: {{ output.metadata.filenames['image/svg+xml'] | urlencode }}\n",
"{% endblock data_svg %}\n",
"\n",
"{% block data_png %}\n",
".. image:: {{ output.metadata.filenames['image/png'] | urlencode }}\n",
"{%- set width=output | get_metadata('width', 'image/png') -%}\n",
"{%- if width is not none %}\n",
" :width: {{ width }}px\n",
"{%- endif %}\n",
"{%- set height=output | get_metadata('height', 'image/png') -%}\n",
"{%- if height is not none %}\n",
" :height: {{ height }}px\n",
"{%- endif %}\n",
"{% endblock data_png %}\n",
"\n",
"{% block data_jpg %}\n",
".. image:: {{ output.metadata.filenames['image/jpeg'] | urlencode }}\n",
"{%- set width=output | get_metadata('width', 'image/jpeg') -%}\n",
"{%- if width is not none %}\n",
" :width: {{ width }}px\n",
"{%- endif %}\n",
"{%- set height=output | get_metadata('height', 'image/jpeg') -%}\n",
"{%- if height is not none %}\n",
" :height: {{ height }}px\n",
"{%- endif %}\n",
"{% endblock data_jpg %}\n",
"\n",
"{% block data_markdown %}\n",
"{{ output.data['text/markdown'] | convert_pandoc(\"markdown\", \"rst\") }}\n",
"{% endblock data_markdown %}\n",
"\n",
"{% block data_latex %}\n",
".. math::\n",
"\n",
"{{ output.data['text/latex'] | strip_dollars | indent }}\n",
"{% endblock data_latex %}\n",
"\n",
"{% block data_text scoped %}\n",
".. parsed-literal::\n",
"\n",
"{{ output.data['text/plain'] | indent }}\n",
"{% endblock data_text %}\n",
"\n",
"{% block data_html scoped %}\n",
".. raw:: html\n",
"\n",
"{{ output.data['text/html'] | indent }}\n",
"{% endblock data_html %}\n",
"\n",
"{% block markdowncell scoped %}\n",
"{{ cell.source | convert_pandoc(\"markdown\", \"rst\") }}\n",
"{% endblock markdowncell %}\n",
"\n",
"{%- block rawcell scoped -%}\n",
"{%- if cell.metadata.get('raw_mimetype', '').lower() in resources.get('raw_mimetypes', ['']) %}\n",
"{{cell.source}}\n",
"{% endif -%}\n",
"{%- endblock rawcell -%}\n",
"\n",
"{% block headingcell scoped %}\n",
"{{ (\"#\" * cell.level + cell.source) | replace('\\n', ' ') | convert_pandoc(\"markdown\", \"rst\") }}\n",
"{% endblock headingcell %}\n",
"\n",
"{% block unknowncell scoped %}\n",
"unknown type {{cell.type}}\n",
"{% endblock unknowncell %}\n"
]
},
{
"cell_type": "code",
"execution_count": 17,
"metadata": {
"ExecuteTime": {
"end_time": "2020-02-16T14:53:49.098172Z",
"start_time": "2020-02-16T14:53:49.093201Z"
}
},
"outputs": [],
"source": [
"from IPython.core.magic import Magics, magics_class, line_cell_magic\n",
"from IPython.core.magic import cell_magic, register_cell_magic, register_line_magic\n",
"from IPython.core.magic_arguments import argument, magic_arguments, parse_argstring\n",
"import subprocess\n",
"import os"
]
},
{
"cell_type": "code",
"execution_count": 18,
"metadata": {
"ExecuteTime": {
"end_time": "2020-02-16T14:53:53.396267Z",
"start_time": "2020-02-16T14:53:53.375754Z"
}
},
"outputs": [],
"source": [
"@magics_class\n",
"class PyboardMagic(Magics):\n",
" @cell_magic\n",
" @magic_arguments()\n",
" @argument('-skip')\n",
" @argument('-unix')\n",
" @argument('-pyboard')\n",
" @argument('-file')\n",
" @argument('-data')\n",
" @argument('-time')\n",
" @argument('-memory')\n",
" def micropython(self, line='', cell=None):\n",
" args = parse_argstring(self.micropython, line)\n",
" if args.skip: # doesn't care about the cell's content\n",
" print('skipped execution')\n",
" return None # do not parse the rest\n",
" if args.unix: # tests the code on the unix port. Note that this works on unix only\n",
" with open('/dev/shm/micropython.py', 'w') as fout:\n",
" fout.write(cell)\n",
" proc = subprocess.Popen([\"../../micropython/ports/unix/micropython\", \"/dev/shm/micropython.py\"], \n",
" stdout=subprocess.PIPE, stderr=subprocess.PIPE)\n",
" print(proc.stdout.read().decode(\"utf-8\"))\n",
" print(proc.stderr.read().decode(\"utf-8\"))\n",
" return None\n",
" if args.file: # can be used to copy the cell content onto the pyboard's flash\n",
" spaces = \" \"\n",
" try:\n",
" with open(args.file, 'w') as fout:\n",
" fout.write(cell.replace('\\t', spaces))\n",
" printf('written cell to {}'.format(args.file))\n",
" except:\n",
" print('Failed to write to disc!')\n",
" return None # do not parse the rest\n",
" if args.data: # can be used to load data from the pyboard directly into kernel space\n",
" message = pyb.exec(cell)\n",
" if len(message) == 0:\n",
" print('pyboard >>>')\n",
" else:\n",
" print(message.decode('utf-8'))\n",
" # register new variable in user namespace\n",
" self.shell.user_ns[args.data] = string_to_matrix(message.decode(\"utf-8\"))\n",
" \n",
" if args.time: # measures the time of executions\n",
" pyb.exec('import utime')\n",
" message = pyb.exec('t = utime.ticks_us()\\n' + cell + '\\ndelta = utime.ticks_diff(utime.ticks_us(), t)' + \n",
" \"\\nprint('execution time: {:d} us'.format(delta))\")\n",
" print(message.decode('utf-8'))\n",
" \n",
" if args.memory: # prints out memory information \n",
" message = pyb.exec('from micropython import mem_info\\nprint(mem_info())\\n')\n",
" print(\"memory before execution:\\n========================\\n\", message.decode('utf-8'))\n",
" message = pyb.exec(cell)\n",
" print(\">>> \", message.decode('utf-8'))\n",
" message = pyb.exec('print(mem_info())')\n",
" print(\"memory after execution:\\n========================\\n\", message.decode('utf-8'))\n",
"\n",
" if args.pyboard:\n",
" message = pyb.exec(cell)\n",
" print(message.decode('utf-8'))\n",
"\n",
"ip = get_ipython()\n",
"ip.register_magics(PyboardMagic)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## pyboard"
]
},
{
"cell_type": "code",
"execution_count": 111,
"metadata": {
"ExecuteTime": {
"end_time": "2020-02-16T18:36:59.172039Z",
"start_time": "2020-02-16T18:36:59.144651Z"
}
},
"outputs": [],
"source": [
"import pyboard\n",
"pyb = pyboard.Pyboard('/dev/ttyACM0')\n",
"pyb.enter_raw_repl()"
]
},
{
"cell_type": "code",
"execution_count": 115,
"metadata": {
"ExecuteTime": {
"end_time": "2020-02-16T18:50:42.907664Z",
"start_time": "2020-02-16T18:50:42.903709Z"
}
},
"outputs": [],
"source": [
"pyb.exit_raw_repl()\n",
"pyb.close()"
]
},
{
"cell_type": "code",
"execution_count": 521,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-20T06:48:05.984879Z",
"start_time": "2019-10-20T06:48:05.619747Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"\n"
]
}
],
"source": [
"%%micropython -pyboard 1\n",
"\n",
"import utime\n",
"import ulab as np\n",
"\n",
"def timeit(n=1000):\n",
" def wrapper(f, *args, **kwargs):\n",
" func_name = str(f).split(' ')[1]\n",
" def new_func(*args, **kwargs):\n",
" run_times = np.zeros(n, dtype=np.uint16)\n",
" for i in range(n):\n",
" t = utime.ticks_us()\n",
" result = f(*args, **kwargs)\n",
" run_times[i] = utime.ticks_diff(utime.ticks_us(), t)\n",
" print('{}() execution times based on {} cycles'.format(func_name, n, (delta2-delta1)/n))\n",
" print('\\tbest: %d us'%np.min(run_times))\n",
" print('\\tworst: %d us'%np.max(run_times))\n",
" print('\\taverage: %d us'%np.mean(run_times))\n",
" print('\\tdeviation: +/-%.3f us'%np.std(run_times)) \n",
" return result\n",
" return new_func\n",
" return wrapper\n",
"\n",
"def timeit(f, *args, **kwargs):\n",
" func_name = str(f).split(' ')[1]\n",
" def new_func(*args, **kwargs):\n",
" t = utime.ticks_us()\n",
" result = f(*args, **kwargs)\n",
" print('execution time: ', utime.ticks_diff(utime.ticks_us(), t), ' us')\n",
" return result\n",
" return new_func"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Introduction"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"In the [last chapter](https://micropython-usermod.readthedocs.io/en/latest/usermods_15.html) of the usermod documentation, I mentioned that I have another story, for another day. The day has come, so here is my story.\n",
"\n",
"## Enter ulab\n",
"\n",
"`ulab` is a numpy-like module for `micropython`, meant to simplify and speed up common mathematical operations on arrays. The primary goal was to implement a small subset of numpy that might be useful in the context of a microcontroller. This means low-level data processing of linear (array) and two-dimensional (matrix) data.\n",
"\n",
"## Purpose\n",
"\n",
"Of course, the first question that one has to answer is, why on Earth one would need a fast math library on a microcontroller. After all, it is not expected that heavy number crunching is going to take place on bare metal. It is not meant to. On a PC, the main reason for writing fast code is the sheer amount of data that one wants to process. On a microcontroller, the data volume is probably small, but it might lead to catastrophic system failure, if these data are not processed in time, because the microcontroller is supposed to interact with the outside world in a timely fashion. In fact, this latter objective was the initiator of this project: I needed the Fourier transform of a signal coming from the ADC of the pyboard, and all available options were simply too slow. \n",
"\n",
"In addition to speed, another issue that one has to keep in mind when working with embedded systems is the amount of available RAM: I believe, everything here could be implemented in pure python with relatively little effort, but the price we would have to pay for that is not only speed, but RAM, too. python code, if is not frozen, and compiled into the firmware, has to be compiled at runtime, which is not exactly a cheap process. On top of that, if numbers are stored in a list or tuple, which would be the high-level container, then they occupy 8 bytes, no matter, whether they are all smaller than 100, or larger than one hundred million. This is obviously a waste of resources in an environment, where resources are scarce. \n",
"\n",
"Finally, there is a reason for using `micropython` in the first place. Namely, that a microcontroller can be programmed in a very elegant, and *pythonic* way. But if it is so, why should we not extend this idea to other tasks and concepts that might come up in this context? If there was no other reason than this *elegance*, I would find that convincing enough.\n",
"\n",
"Based on the above-mentioned considerations, all functions in `ulab` are implemented in a way that \n",
"\n",
"1. conforms to `numpy` as much as possible\n",
"2. is so frugal with RAM as possible,\n",
"3. and yet, fast. Much faster than pure python.\n",
"\n",
"The main points of `ulab` are \n",
"\n",
"- compact, iterable and slicable containers of numerical data in 1, and 2 dimensions (arrays and matrices). These containers support all the relevant unary and binary operators (e.g., `len`, ==, +, *, etc.)\n",
"- vectorised computations on micropython iterables and numerical arrays/matrices (in `numpy`-speak, universal functions)\n",
"- basic linear algebra routines (matrix inversion, multiplication, reshaping, transposition, determinant, and eigenvalues)\n",
"- polynomial fits to numerical data\n",
"- fast Fourier transforms\n",
"\n",
"At the time of writing this manual (for version 0.33.2), the library adds approximately 30 kB of extra compiled code to the micropython (pyboard.v.11) firmware. However, if you are tight with flash space, you can easily shave off a couple of kB. See the section on [customising ulab](#Custom_builds).\n",
"\n",
"## Resources and legal matters\n",
"\n",
"The source code of the module can be found under https://github.com/v923z/micropython-ulab/tree/master/code. The source of this user manual is under https://github.com/v923z/micropython-ulab/tree/master/docs, while the technical details of the implementation are discussed at great length in https://github.com/v923z/micropython-ulab/tree/master/docs/ulab.ipynb. If you want an even thorougher explanation on why the various constructs of the implementation work, and work in that particular way, you can read more on the subject under https://micropython-usermod.readthedocs.io/en/latest/, where I demonstrate, what you have to do, if you want to make a C object behave in a *pythonic* way. \n",
"\n",
"The MIT licence applies to all material. \n",
"\n",
"## Friendly request\n",
"\n",
"If you use `ulab`, and bump into a bug, or think that a particular function is missing, or its behaviour does not conform to `numpy`, please, raise a [ulab issue](#https://github.com/v923z/micropython-ulab/issues) on github, so that the community can profit from your experiences. \n",
"\n",
"Even better, if you find the project useful, and think that it could be made better, faster, tighter, and shinier, please, consider contributing, and issue a pull request with the implementation of your improvements and new features. `ulab` can only become successful, if it offers what the community needs.\n",
"\n",
"These last comments apply to the documentation, too. If, in your opinion, the documentation is obscure, misleading, or not detailed enough, please, let me know, so that *we* can fix it.\n",
"\n",
"## Differences between micropython-ulab and circuitpython-ulab\n",
"\n",
"`ulab` has originally been developed for `micropython`, but has since been integrated into a number of its flavours. Most of these flavours are simply forks of `micropython` itself, with some additional functionality. One of the notable exceptions is `circuitpython`, which has slightly diverged at the core level, and this has some minor consequences. Some of these concern the C implementation details only, which all have been sorted out with the generous and enthusiastic support of Jeff Epler from [Adafruit Industries](http://www.adafruit.com).\n",
"\n",
"There are, however, a couple of instances, where the usage in the two environments is slightly different at the python level. These are how the packges can be imported, and how the class properties can be accessed. In both cases, the `circuitpython` implementation results in `numpy`-conform code. `numpy`-compatibility in `micropython` will be implemented as soon as `micropython` itself has the required tools. Till then we have to live with a workaround, which I will point out at the relevant places."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Customising `ulab`\n",
"\n",
"`ulab` implements a great number of functions, which are organised in sub-modules. E.g., functions related to Fourier transforms are located in the `ulab.fft` sub-module, so you would import `fft` as\n",
"\n",
"```python\n",
"import ulab\n",
"from ulab import fft\n",
"```\n",
"by which point you can get the FFT of your data by calling `fft.fft(...)`. \n",
"\n",
"The idea of such grouping of functions and methods is to provide a means for granularity: It is quite possible that you do not need all functions in a particular application. If you want to save some flash space, you can easily exclude arbitrary sub-modules from the firmware. The [ulab.h](https://github.com/v923z/micropython-ulab/blob/master/code/ulab.h) header file contains a pre-processor flag for each sub-module. The default setting is 1 for each of them. Setting them to 0 removes the module from the compiled firmware. \n",
"\n",
"The first couple of lines of the file look like this\n",
"\n",
"```c\n",
"// vectorise (all functions) takes approx. 3 kB of flash space\n",
"#define ULAB_VECTORISE_MODULE (1)\n",
"\n",
"// linalg adds around 6 kB\n",
"#define ULAB_LINALG_MODULE (1)\n",
"\n",
"// poly is approx. 2.5 kB\n",
"#define ULAB_POLY_MODULE (1)\n",
"```\n",
"\n",
"In order to simplify navigation in the header, each flag begins with `ULAB_`, and continues with the name of the sub-module. This name is also the `.c` file, where the sub-module is implemented. So, e.g., the linear algebra routines can be found in `linalg.c`, and the corresponding compiler flag is `ULAB_LINALG_MODULE`. Each section displays a hint as to how much space you can save by un-setting the flag.\n",
"\n",
"At first, having to import everything in this way might appear to be overly complicated, but there is a very good reason behind all this: you can find out at the time of importing, whether a function or sub-module is part of your `ulab` firmware, or not. The alternative, namely, that you do not have to import anything beyond `ulab`, could prove catastrophic: you would learn only at run time (at the moment of calling the function in your code) that a particular function is not in the firmware, and that is most probably too late.\n",
"\n",
"Except for `fft`, the standard sub-modules, `vector`, `linalg`, `numerical`, `and poly`all `numpy`-compatible. User-defined functions that accept `ndarray`s as their argument should be implemented in the `extra` sub-module, or its sub-modules. Hints as to how to do that can be found in the section [Extending ulab](#Extending-ulab)."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Supported functions and methods\n",
"\n",
"`ulab` supports a number of array operators, which are listed here. I tried to follow the specifications of the `numpy` interface as closely as possible, though, it was not always practical to implement verbatim behaviour. The differences, if any, are in each case small (e.g., a function cannot take all possible keyword arguments), and should not hinder everyday use. In the list below, a single asterisk denotes slight deviations from `numpy`'s nomenclature, and a double asterisk denotes those cases, where a bit more caution should be exercised, though this usually means functions that are not supported by `numpy`.\n",
"\n",
"The detailed discussion of the various functions always contains a link to the corresponding `numpy` documentation. However, before going down the rabbit hole, the module also defines a constant, the version, which can always be queried as "
]
},
{
"cell_type": "code",
"execution_count": 209,
"metadata": {
"ExecuteTime": {
"end_time": "2019-11-01T13:19:05.158219Z",
"start_time": "2019-11-01T13:19:05.140455Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"you are running ulab version 0.24\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"print('you are running ulab version', np.__version__)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"If you find a bug, please, include this number in your report!"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Basic ndarray operations\n",
"\n",
"[Unary operators](#Unary-operators)\n",
"\n",
"[Binary operators](#Binary-operators)\n",
"\n",
"[Indexing and slicing](#Slicing-and-indexing)\n",
"\n",
"[ndarray iterators](#Iterating-over-arrays)\n",
"\n",
"[Comparison operators<sup>*</sup>](#Comparison-operators)\n",
"\n",
"[Universal functions](#Universal-functions) (also support function calls on general iterables)\n",
"\n",
"\n",
"## Methods of ndarrays\n",
"\n",
"[.shape<sup>*</sup>](#.shape)\n",
"\n",
"[size<sup>*</sup>](#size)\n",
"\n",
"[itemsize<sup>*</sup>](#itemsize)\n",
"\n",
"[.reshape](#.reshape)\n",
"\n",
"[.transpose](#.transpose)\n",
"\n",
"[.flatten<sup>**</sup>](#.flatten)\n",
"\n",
"## Matrix methods\n",
"\n",
"[inv](#inv)\n",
"\n",
"[dot](#dot)\n",
"\n",
"[det](#det)\n",
"\n",
"[roll](#roll)\n",
"\n",
"[flip](#flip)\n",
"\n",
"## Array initialisation functions\n",
"\n",
"[eye](#eye)\n",
"\n",
"[ones](#ones,-zeros)\n",
"\n",
"[zeros](#ones,-zeros)\n",
"\n",
"[linspace](#linspace)\n",
"\n",
"## Statistical and other properties of arrays\n",
"\n",
"[min](#min,-argmin,-max,-argmax)\n",
"\n",
"[argmin](#min,-argmin,-max,-argmax)\n",
"\n",
"[max](#min,-argmin,-max,-argmax)\n",
"\n",
"[argmax](#min,-argmin,-max,-argmax)\n",
"\n",
"[sum](#sum,-std,-mean)\n",
"\n",
"[std](#sum,-std,-mean)\n",
"\n",
"[mean](#sum,-std,-mean)\n",
"\n",
"[diff](#diff)\n",
"\n",
"[sort](#sort)\n",
"\n",
"[argsort](#argsort)\n",
"\n",
"## Manipulation of polynomials\n",
"\n",
"[polyval](#polyval)\n",
"\n",
"[polyfit](#polyfit)\n",
"\n",
"## FFT routines\n",
"\n",
"[fft<sup>**</sup>](#fft)\n",
"\n",
"[ifft<sup>**</sup>](#ifft)\n",
"\n",
"[spectrum<sup>**</sup>](#spectrum)\n",
"\n",
"## Filter functions\n",
"\n",
"[convolve](#convolve)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# ndarray, the basic container\n",
"\n",
"The `ndarray` is the underlying container of numerical data. It is derived from micropython's own `array` object, but has a great number of extra features starting with how it can be initialised, which operations can be done on it, and which functions can accept it as an argument. One important property of an `ndarray` is that it is also a proper `micropython` iterable.\n",
"\n",
"Since the `ndarray` is a binary container, it is also compact, meaning that it takes only a couple of bytes of extra RAM in addition to what is required for storing the numbers themselves. `ndarray`s are also type-aware, i.e., one can save RAM by specifying a data type, and using the smallest reasonable one. Five such types are defined, namely `uint8`, `int8`, which occupy a single byte of memory per datum, `uint16`, and `int16`, which occupy two bytes per datum, and `float`, which occupies four or eight bytes per datum. The precision/size of the `float` type depends on the definition of `mp_float_t`. Some platforms, e.g., the PYBD, implement `double`s, but some, e.g., the pyboard.v.11, don't. You can find out, what type of float your particular platform implements by looking at the output of the [.itemsize](#.itemsize) class property.\n",
"\n",
"On the following pages, we will see how one can work with `ndarray`s. Those familiar with `numpy` should find that the nomenclature and naming conventions of `numpy` are adhered to as closely as possible. I will point out the few differences, where necessary.\n",
"\n",
"For the sake of comparison, in addition to the `ulab` code snippets, sometimes the equivalent `numpy` code is also presented. You can find out, where the snippet is supposed to run by looking at its first line, the header.\n",
"\n",
"Hint: you can easily port existing `numpy` code, if you `import ulab as np`."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Initialising an array\n",
"\n",
"A new array can be created by passing either a standard micropython iterable, or another `ndarray` into the constructor."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Initialising by passing iterables\n",
"\n",
"If the iterable is one-dimensional, i.e., one whose elements are numbers, then a row vector will be created and returned. If the iterable is two-dimensional, i.e., one whose elements are again iterables, a matrix will be created. If the lengths of the iterables is not consistent, a `ValueError` will be raised. Iterables of different types can be mixed in the initialisation function. \n",
"\n",
"If the `dtype` keyword with the possible `uint8/int8/uint16/int16/float` values is supplied, the new `ndarray` will have that type, otherwise, it assumes `float` as default. "
]
},
{
"cell_type": "code",
"execution_count": 174,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-11T17:43:37.591149Z",
"start_time": "2019-10-11T17:43:37.571853Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a:\t [1, 2, 3, 4, 5, 6, 7, 8]\n",
"b:\t array([1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0], dtype=float)\n",
"\n",
"c:\t array([[0, 1, 2, 3, 4],\n",
"\t [20, 21, 22, 23, 24],\n",
"\t [44, 55, 66, 77, 88]], dtype=uint8)\n",
"\n",
"Traceback (most recent call last):\n",
" File \"/dev/shm/micropython.py\", line 15, in <module>\n",
"ValueError: iterables are not of the same length\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = [1, 2, 3, 4, 5, 6, 7, 8]\n",
"b = np.array(a)\n",
"\n",
"print(\"a:\\t\", a)\n",
"print(\"b:\\t\", b)\n",
"\n",
"# a two-dimensional array with mixed-type initialisers\n",
"c = np.array([range(5), range(20, 25, 1), [44, 55, 66, 77, 88]], dtype=np.uint8)\n",
"print(\"\\nc:\\t\", c)\n",
"\n",
"# and now we throw an exception\n",
"d = np.array([range(5), range(10), [44, 55, 66, 77, 88]], dtype=np.uint8)\n",
"print(\"\\nd:\\t\", d)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"`ndarray`s are pretty-printed, i.e., if the length is larger than 10, then only the first and last three entries will be printed. Also note that, as opposed to `numpy`, the printout always contains the `dtype`."
]
},
{
"cell_type": "code",
"execution_count": 448,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-19T12:09:32.898039Z",
"start_time": "2019-10-19T12:09:32.877872Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a:\t array([0.0, 1.0, 2.0, ..., 197.0, 198.0, 199.0], dtype=float)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array(range(200))\n",
"print(\"a:\\t\", a)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Initialising by passing arrays\n",
"\n",
"An `ndarray` can be initialised by supplying another array. This statement is almost trivial, since `ndarray`s are iterables themselves, though it should be pointed out that initialising through arrays should be faster, because simply a new copy is created, without inspection, iteration etc."
]
},
{
"cell_type": "code",
"execution_count": 175,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-11T17:46:37.653805Z",
"start_time": "2019-10-11T17:46:37.641837Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a:\t [1, 2, 3, 4, 5, 6, 7, 8]\n",
"b:\t array([1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0], dtype=float)\n",
"\n",
"c:\t array([1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0], dtype=float)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = [1, 2, 3, 4, 5, 6, 7, 8]\n",
"b = np.array(a)\n",
"c = np.array(b)\n",
"\n",
"print(\"a:\\t\", a)\n",
"print(\"b:\\t\", b)\n",
"print(\"\\nc:\\t\", c)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Methods of ndarrays"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### .shape\n",
"\n",
"The `.shape` method (property) returns a 2-tuple with the number of rows, and columns. \n",
"\n",
"**WARNING:** In `circuitpython`, you can call the method as a property, i.e., "
]
},
{
"cell_type": "code",
"execution_count": 12,
"metadata": {
"ExecuteTime": {
"end_time": "2020-02-11T19:01:40.377272Z",
"start_time": "2020-02-11T19:01:40.364448Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a:\n",
" array([1, 2, 3, 4], dtype=int8)\n",
"shape of a: (1, 4)\n",
"\n",
"b:\n",
" array([[1, 2],\n",
"\t [3, 4]], dtype=int8)\n",
"shape of b: (2, 2)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array([1, 2, 3, 4], dtype=np.int8)\n",
"print(\"a:\\n\", a)\n",
"print(\"shape of a:\", a.shape)\n",
"\n",
"b= np.array([[1, 2], [3, 4]], dtype=np.int8)\n",
"print(\"\\nb:\\n\", b)\n",
"print(\"shape of b:\", b.shape)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"**WARNING:** On the other hand, since properties are not implemented in `micropython`, there you would call the method as a function, i.e., "
]
},
{
"cell_type": "code",
"execution_count": 12,
"metadata": {
"ExecuteTime": {
"end_time": "2020-02-11T19:01:40.377272Z",
"start_time": "2020-02-11T19:01:40.364448Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a:\n",
" array([1, 2, 3, 4], dtype=int8)\n",
"shape of a: (1, 4)\n",
"\n",
"b:\n",
" array([[1, 2],\n",
"\t [3, 4]], dtype=int8)\n",
"shape of b: (2, 2)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array([1, 2, 3, 4], dtype=np.int8)\n",
"print(\"a:\\n\", a)\n",
"print(\"shape of a:\", a.shape)\n",
"\n",
"b= np.array([[1, 2], [3, 4]], dtype=np.int8)\n",
"print(\"\\nb:\\n\", b)\n",
"print(\"shape of b:\", b.shape())"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### .size\n",
"\n",
"The `.size` method (property) returns an integer with the number of elements in the array. \n",
"\n",
"**WARNING:** In `circuitpython`, the `numpy` nomenclature applies, i.e., "
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {
"ExecuteTime": {
"end_time": "2020-02-11T06:32:22.721112Z",
"start_time": "2020-02-11T06:32:22.713111Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a:\n",
" array([1, 2, 3], dtype=int8)\n",
"size of a: 3\n",
"\n",
"b:\n",
" array([[1, 2],\n",
"\t [3, 4]], dtype=int8)\n",
"size of b: 4\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array([1, 2, 3], dtype=np.int8)\n",
"print(\"a:\\n\", a)\n",
"print(\"size of a:\", a.size)\n",
"\n",
"b= np.array([[1, 2], [3, 4]], dtype=np.int8)\n",
"print(\"\\nb:\\n\", b)\n",
"print(\"size of b:\", b.size)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"**WARNING:** In `micropython`, `size` is a method, i.e., "
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {
"ExecuteTime": {
"end_time": "2020-02-11T06:32:22.721112Z",
"start_time": "2020-02-11T06:32:22.713111Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a:\n",
" array([1, 2, 3], dtype=int8)\n",
"size of a: 3\n",
"\n",
"b:\n",
" array([[1, 2],\n",
"\t [3, 4]], dtype=int8)\n",
"size of b: 4\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array([1, 2, 3], dtype=np.int8)\n",
"print(\"a:\\n\", a)\n",
"print(\"size of a:\", a.size)\n",
"\n",
"b= np.array([[1, 2], [3, 4]], dtype=np.int8)\n",
"print(\"\\nb:\\n\", b)\n",
"print(\"size of b:\", b.size())"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### .itemsize\n",
"\n",
"The `.itemsize` method (property) returns an integer with the siz enumber of elements in the array.\n",
"\n",
"**WARNING:** In `circuitpython`:"
]
},
{
"cell_type": "code",
"execution_count": 13,
"metadata": {
"ExecuteTime": {
"end_time": "2020-02-11T19:05:04.296601Z",
"start_time": "2020-02-11T19:05:04.280669Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a:\n",
" array([1, 2, 3], dtype=int8)\n",
"itemsize of a: 1\n",
"\n",
"b:\n",
" array([[1.0, 2.0],\n",
"\t [3.0, 4.0]], dtype=float)\n",
"itemsize of b: 8\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array([1, 2, 3], dtype=np.int8)\n",
"print(\"a:\\n\", a)\n",
"print(\"itemsize of a:\", a.itemsize)\n",
"\n",
"b= np.array([[1, 2], [3, 4]], dtype=np.float)\n",
"print(\"\\nb:\\n\", b)\n",
"print(\"itemsize of b:\", b.itemsize)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"**WARNING:** In `micropython`:"
]
},
{
"cell_type": "code",
"execution_count": 13,
"metadata": {
"ExecuteTime": {
"end_time": "2020-02-11T19:05:04.296601Z",
"start_time": "2020-02-11T19:05:04.280669Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a:\n",
" array([1, 2, 3], dtype=int8)\n",
"itemsize of a: 1\n",
"\n",
"b:\n",
" array([[1.0, 2.0],\n",
"\t [3.0, 4.0]], dtype=float)\n",
"itemsize of b: 8\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array([1, 2, 3], dtype=np.int8)\n",
"print(\"a:\\n\", a)\n",
"print(\"itemsize of a:\", a.itemsize)\n",
"\n",
"b= np.array([[1, 2], [3, 4]], dtype=np.float)\n",
"print(\"\\nb:\\n\", b)\n",
"print(\"itemsize of b:\", b.itemsize())"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### .reshape\n",
"\n",
"numpy: https://docs.scipy.org/doc/numpy/reference/generated/numpy.reshape.html\n",
"\n",
"`reshape` re-writes the shape properties of an `ndarray`, but the array will not be modified in any other way. The function takes a single 2-tuple with two integers as its argument. The 2-tuple should specify the desired number of rows and columns. If the new shape is not consistent with the old, a `ValueError` exception will be raised."
]
},
{
"cell_type": "code",
"execution_count": 379,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-18T13:17:12.656602Z",
"start_time": "2019-10-18T13:17:12.641007Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a (4 by 4): array([[1, 2, 3, 4],\n",
"\t [5, 6, 7, 8],\n",
"\t [9, 10, 11, 12],\n",
"\t [13, 14, 15, 16]], dtype=uint8)\n",
"a (2 by 8): array([[1, 2, 3, 4, 5, 6, 7, 8],\n",
"\t [9, 10, 11, 12, 13, 14, 15, 16]], dtype=uint8)\n",
"a (1 by 16): array([1, 2, 3, ..., 14, 15, 16], dtype=uint8)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array([[1, 2, 3, 4], [5, 6, 7, 8], [9, 10, 11, 12], [13, 14, 15, 16]], dtype=np.uint8)\n",
"print('a (4 by 4):', a)\n",
"print('a (2 by 8):', a.reshape((2, 8)))\n",
"print('a (1 by 16):', a.reshape((1, 16)))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### .flatten\n",
"\n",
"numpy: https://docs.scipy.org/doc/numpy/reference/generated/numpy.ndarray.flatten.htm\n",
"\n",
"`.flatten` returns the flattened array. The array can be flattened in `C` style (i.e., moving horizontally in the matrix), or in `fortran` style (i.e., moving vertically in the matrix). The `C`-style flattening is the default, and it is also fast, because this is just a verbatim copy of the contents."
]
},
{
"cell_type": "code",
"execution_count": 335,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-16T18:04:12.241634Z",
"start_time": "2019-10-16T18:04:12.226130Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a: \t\t array([1, 2, 3, 4], dtype=int8)\n",
"a flattened: \t array([1, 2, 3, 4], dtype=int8)\n",
"\n",
"b: array([[1, 2, 3],\n",
"\t [4, 5, 6]], dtype=int8)\n",
"b flattened (C): \t array([1, 2, 3, 4, 5, 6], dtype=int8)\n",
"b flattened (F): \t array([1, 4, 2, 5, 3, 6], dtype=int8)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array([1, 2, 3, 4], dtype=np.int8)\n",
"print(\"a: \\t\\t\", a)\n",
"print(\"a flattened: \\t\", a.flatten())\n",
"\n",
"b = np.array([[1, 2, 3], [4, 5, 6]], dtype=np.int8)\n",
"print(\"\\nb:\", b)\n",
"\n",
"print(\"b flattened (C): \\t\", b.flatten())\n",
"print(\"b flattened (F): \\t\", b.flatten(order='F'))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### .transpose\n",
"\n",
"numpy: https://docs.scipy.org/doc/numpy/reference/generated/numpy.transpose.html\n",
"\n",
"Note that only square matrices can be transposed in place, and in general, an internal copy of the matrix is required. If RAM is a concern, plan accordingly!"
]
},
{
"cell_type": "code",
"execution_count": 384,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-19T08:39:11.844987Z",
"start_time": "2019-10-19T08:39:11.828099Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a:\n",
" array([[1, 2, 3],\n",
"\t [4, 5, 6],\n",
"\t [7, 8, 9],\n",
"\t [10, 11, 12]], dtype=uint8)\n",
"shape of a: (4, 3)\n",
"\n",
"transpose of a:\n",
" array([[1, 4, 7, 10],\n",
"\t [2, 5, 8, 11],\n",
"\t [3, 6, 9, 12]], dtype=uint8)\n",
"shape of a: (3, 4)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9], [10, 11, 12]], dtype=np.uint8)\n",
"print('a:\\n', a)\n",
"print('shape of a:', a.shape())\n",
"a.transpose()\n",
"print('\\ntranspose of a:\\n', a)\n",
"print('shape of a:', a.shape())"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### .sort\n",
"\n",
"numpy: https://docs.scipy.org/doc/numpy/reference/generated/numpy.sort.html\n",
"\n",
"In-place sorting of an `ndarray`. For a more detailed exposition, see [sort](#sort)."
]
},
{
"cell_type": "code",
"execution_count": 14,
"metadata": {
"ExecuteTime": {
"end_time": "2019-11-05T06:23:39.251388Z",
"start_time": "2019-11-05T06:23:39.235989Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"\n",
"a:\n",
" array([[1, 12, 3, 0],\n",
"\t [5, 3, 4, 1],\n",
"\t [9, 11, 1, 8],\n",
"\t [7, 10, 0, 1]], dtype=uint8)\n",
"\n",
"a sorted along vertical axis:\n",
" array([[1, 3, 0, 0],\n",
"\t [5, 10, 1, 1],\n",
"\t [7, 11, 3, 1],\n",
"\t [9, 12, 4, 8]], dtype=uint8)\n",
"\n",
"a sorted along horizontal axis:\n",
" array([[0, 1, 3, 12],\n",
"\t [1, 3, 4, 5],\n",
"\t [1, 8, 9, 11],\n",
"\t [0, 1, 7, 10]], dtype=uint8)\n",
"\n",
"flattened a sorted:\n",
" array([0, 0, 1, ..., 10, 11, 12], dtype=uint8)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array([[1, 12, 3, 0], [5, 3, 4, 1], [9, 11, 1, 8], [7, 10, 0, 1]], dtype=np.uint8)\n",
"print('\\na:\\n', a)\n",
"a.sort(axis=0)\n",
"print('\\na sorted along vertical axis:\\n', a)\n",
"\n",
"a = np.array([[1, 12, 3, 0], [5, 3, 4, 1], [9, 11, 1, 8], [7, 10, 0, 1]], dtype=np.uint8)\n",
"a.sort(a, axis=1)\n",
"print('\\na sorted along horizontal axis:\\n', a)\n",
"\n",
"a = np.array([[1, 12, 3, 0], [5, 3, 4, 1], [9, 11, 1, 8], [7, 10, 0, 1]], dtype=np.uint8)\n",
"a.sort(a, axis=None)\n",
"print('\\nflattened a sorted:\\n', a)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Unary operators\n",
"\n",
"With the exception of `len`, which returns a single number, all unary operators manipulate the underlying data element-wise. "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### len\n",
"\n",
"This operator takes a single argument, and returns either the length (for row vectors), or the number of rows (for matrices) of its argument."
]
},
{
"cell_type": "code",
"execution_count": 105,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-11T13:23:54.932077Z",
"start_time": "2019-10-11T13:23:54.911337Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a:\t array([1, 2, 3, 4, 5], dtype=uint8)\n",
"length of a: 5\n",
"shape of a: (1, 5)\n",
"\n",
"b:\t array([[0, 1, 2, 3, 4],\n",
"\t [0, 1, 2, 3, 4],\n",
"\t [0, 1, 2, 3, 4],\n",
"\t [0, 1, 2, 3, 4]], dtype=uint8)\n",
"length of b: 4\n",
"shape of b: (4, 5)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array([1, 2, 3, 4, 5], dtype=np.uint8)\n",
"b = np.array([range(5), range(5), range(5), range(5)], dtype=np.uint8)\n",
"\n",
"print(\"a:\\t\", a)\n",
"print(\"length of a: \", len(a))\n",
"print(\"shape of a: \", a.shape())\n",
"print(\"\\nb:\\t\", b)\n",
"print(\"length of b: \", len(b))\n",
"print(\"shape of b: \", b.shape())"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
" The number returned by `len` is also the length of the iterations, when the array supplies the elements for an iteration (see later)."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### invert\n",
"\n",
"The function function is defined for integer data types (`uint8`, `int8`, `uint16`, and `int16`) only, takes a single argument, and returns the element-by-element, bit-wise inverse of the array. If a `float` is supplied, the function raises a `ValueError` exception.\n",
"\n",
"With signed integers (`int8`, and `int16`), the results might be unexpected, as in the example below:"
]
},
{
"cell_type": "code",
"execution_count": 98,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-11T13:16:16.754210Z",
"start_time": "2019-10-11T13:16:16.735618Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a:\t\t array([0, -1, -100], dtype=int8)\n",
"inverse of a:\t array([-1, 0, 99], dtype=int8)\n",
"\n",
"a:\t\t array([0, 1, 254, 255], dtype=uint8)\n",
"inverse of a:\t array([255, 254, 1, 0], dtype=uint8)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array([0, -1, -100], dtype=np.int8)\n",
"print(\"a:\\t\\t\", a)\n",
"print(\"inverse of a:\\t\", ~a)\n",
"\n",
"a = np.array([0, 1, 254, 255], dtype=np.uint8)\n",
"print(\"\\na:\\t\\t\", a)\n",
"print(\"inverse of a:\\t\", ~a)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### abs\n",
"\n",
"This function takes a single argument, and returns the element-by-element absolute value of the array. When the data type is unsigned (`uint8`, or `uint16`), a copy of the array will be returned immediately, and no calculation takes place."
]
},
{
"cell_type": "code",
"execution_count": 73,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-11T13:05:43.926821Z",
"start_time": "2019-10-11T13:05:43.912629Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a:\t\t\t array([0, -1, -100], dtype=int8)\n",
"absolute value of a:\t array([0, 1, 100], dtype=int8)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array([0, -1, -100], dtype=np.int8)\n",
"print(\"a:\\t\\t\\t \", a)\n",
"print(\"absolute value of a:\\t \", abs(a))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### neg\n",
"\n",
"This operator takes a single argument, and changes the sign of each element in the array. Unsigned values are wrapped. "
]
},
{
"cell_type": "code",
"execution_count": 99,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-11T13:17:00.946009Z",
"start_time": "2019-10-11T13:17:00.927264Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a:\t\t array([10, -1, 1], dtype=int8)\n",
"negative of a:\t array([-10, 1, -1], dtype=int8)\n",
"\n",
"b:\t\t array([0, 100, 200], dtype=uint8)\n",
"negative of b:\t array([0, 156, 56], dtype=uint8)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array([10, -1, 1], dtype=np.int8)\n",
"print(\"a:\\t\\t\", a)\n",
"print(\"negative of a:\\t\", -a)\n",
"\n",
"b = np.array([0, 100, 200], dtype=np.uint8)\n",
"print(\"\\nb:\\t\\t\", b)\n",
"print(\"negative of b:\\t\", -b)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### pos\n",
"\n",
"This function takes a single argument, and simply returns a copy of the array."
]
},
{
"cell_type": "code",
"execution_count": 85,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-11T13:09:15.965662Z",
"start_time": "2019-10-11T13:09:15.945461Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a:\t\t array([10, -1, 1], dtype=int8)\n",
"positive of a:\t array([10, -1, 1], dtype=int8)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array([10, -1, 1], dtype=np.int8)\n",
"print(\"a:\\t\\t\", a)\n",
"print(\"positive of a:\\t\", +a)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Binary operators\n",
"\n",
"All binary operators work element-wise. This also means that the operands either must have the same shape, or one of them must be a scalar.\n",
"\n",
"**WARNING:** `numpy` also allows operations between a matrix, and a row vector, if the row vector has exactly as many elements, as many columns the matrix has. This feature will be added in future versions of `ulab`."
]
},
{
"cell_type": "code",
"execution_count": 188,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-11T18:59:35.810060Z",
"start_time": "2019-10-11T18:59:35.804683Z"
}
},
"outputs": [
{
"data": {
"text/plain": [
"array([[11, 22, 33],\n",
" [14, 25, 36],\n",
" [17, 28, 36]])"
]
},
"execution_count": 188,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"a = array([[1, 2, 3], [4, 5, 6], [7, 8, 6]])\n",
"b = array([10, 20, 30])\n",
"a+b"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Upcasting\n",
"\n",
"Binary operations require special attention, because two arrays with different typecodes can be the operands of an operation, in which case it is not trivial, what the typecode of the result is. This decision on the result's typecode is called upcasting. Since the number of typecodes in `ulab` is significantly smaller than in `numpy`, we have to define new upcasting rules. Where possible, I followed `numpy`'s conventions. \n",
"\n",
"`ulab` observes the following upcasting rules:\n",
"\n",
"1. Operations with two `ndarray`s of the same `dtype` preserve their `dtype`, even when the results overflow.\n",
"\n",
"2. if either of the operands is a float, the result is automatically a float\n",
"\n",
"3. When the right hand side of a binary operator is a micropython variable, `mp_obj_int`, or `mp_obj_float`, then the result will be promoted to `dtype` `float`. This is necessary, because a micropython integer can be 31 bites wide. Other micropython types (e.g., lists, tuples, etc.) raise a `TypeError` exception. \n",
"\n",
"4. \n",
" \n",
"| left hand side | right hand side | ulab result | numpy result |\n",
"|----------------|-----------------|-------------|--------------|\n",
"|`uint8` |`int8` |`int16` |`int16` |\n",
"|`uint8` |`int16` |`int16` |`int16` |\n",
"|`uint8` |`uint16` |`uint16` |`uint16` |\n",
"|`int8` |`int16` |`int16` |`int16` | \n",
"|`int8` |`uint16` |`uint16` |`int32` |\n",
"|`uint16` |`int16` |`float` |`int32` |\n",
" \n",
"Note that the last two operations are promoted to `int32` in `numpy`.\n",
" \n",
"**WARNING:** Due to the lower number of available data types, the upcasting rules of `ulab` are slightly different to those of `numpy`. Watch out for this, when porting code!\n",
"\n",
"Upcasting can be seen in action in the following snippet:"
]
},
{
"cell_type": "code",
"execution_count": 212,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-11T19:28:50.688741Z",
"start_time": "2019-10-11T19:28:50.674459Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a:\t array([1, 2, 3, 4], dtype=uint8)\n",
"b:\t array([1, 2, 3, 4], dtype=int8)\n",
"a+b:\t array([2, 4, 6, 8], dtype=int16)\n",
"\n",
"a:\t array([1, 2, 3, 4], dtype=uint8)\n",
"c:\t array([1.0, 2.0, 3.0, 4.0], dtype=float)\n",
"a*c:\t array([1.0, 4.0, 9.0, 16.0], dtype=float)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array([1, 2, 3, 4], dtype=np.uint8)\n",
"b = np.array([1, 2, 3, 4], dtype=np.int8)\n",
"print(\"a:\\t\", a)\n",
"print(\"b:\\t\", b)\n",
"print(\"a+b:\\t\", a+b)\n",
"\n",
"c = np.array([1, 2, 3, 4], dtype=np.float)\n",
"print(\"\\na:\\t\", a)\n",
"print(\"c:\\t\", c)\n",
"print(\"a*c:\\t\", a*c)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"**WARNING:** If a binary operation involves an `ndarray` and a micropython type (integer, or float), then the array must be on the left hand side. "
]
},
{
"cell_type": "code",
"execution_count": 181,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-11T18:38:30.369710Z",
"start_time": "2019-10-11T18:38:30.354603Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a:\t array([1, 2, 3, 4], dtype=uint8)\n",
"b:\t 12\n",
"a+b:\t array([13, 14, 15, 16], dtype=uint8)\n",
"\n",
"Traceback (most recent call last):\n",
" File \"/dev/shm/micropython.py\", line 12, in <module>\n",
"TypeError: unsupported types for __add__: 'int', 'ndarray'\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"# this is going to work\n",
"a = np.array([1, 2, 3, 4], dtype=np.uint8)\n",
"b = 12\n",
"print(\"a:\\t\", a)\n",
"print(\"b:\\t\", b)\n",
"print(\"a+b:\\t\", a+b)\n",
"\n",
"# but this will spectacularly fail\n",
"print(\"b+a:\\t\", b+a)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The reason for this lies in how micropython resolves binary operators, and this means that a fix can only be implemented, if micropython itself changes the corresponding function(s). Till then, keep `ndarray`s on the left hand side. "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Benchmarks\n",
"\n",
"The following snippet compares the performance of binary operations to a possible implementation in python. For the time measurement, we will take the following snippet from the micropython manual:"
]
},
{
"cell_type": "code",
"execution_count": 522,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-20T06:49:27.584150Z",
"start_time": "2019-10-20T06:49:27.551381Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"\n"
]
}
],
"source": [
"%%micropython -pyboard 1\n",
"\n",
"def timeit(f, *args, **kwargs):\n",
" func_name = str(f).split(' ')[1]\n",
" def new_func(*args, **kwargs):\n",
" t = utime.ticks_us()\n",
" result = f(*args, **kwargs)\n",
" print('execution time: ', utime.ticks_diff(utime.ticks_us(), t), ' us')\n",
" return result\n",
" return new_func"
]
},
{
"cell_type": "code",
"execution_count": 490,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-19T13:23:45.432395Z",
"start_time": "2019-10-19T13:23:45.344021Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"python add:\n",
"execution time: 10051 us\n",
"\n",
"python multiply:\n",
"execution time: 14175 us\n",
"\n",
"ulab add:\n",
"execution time: 222 us\n",
"\n",
"ulab multiply:\n",
"execution time: 213 us\n",
"\n"
]
}
],
"source": [
"%%micropython -pyboard 1\n",
"\n",
"import ulab as np\n",
"\n",
"@timeit\n",
"def py_add(a, b):\n",
" return [a[i]+b[i] for i in range(1000)]\n",
"\n",
"@timeit\n",
"def py_multiply(a, b):\n",
" return [a[i]*b[i] for i in range(1000)]\n",
"\n",
"@timeit\n",
"def ulab_add(a, b):\n",
" return a + b\n",
"\n",
"@timeit\n",
"def ulab_multiply(a, b):\n",
" return a * b\n",
"\n",
"a = [0.0]*1000\n",
"b = range(1000)\n",
"\n",
"print('python add:')\n",
"py_add(a, b)\n",
"\n",
"print('\\npython multiply:')\n",
"py_multiply(a, b)\n",
"\n",
"a = np.linspace(0, 10, num=1000)\n",
"b = np.ones(1000)\n",
"\n",
"print('\\nulab add:')\n",
"ulab_add(a, b)\n",
"\n",
"print('\\nulab multiply:')\n",
"ulab_multiply(a, b)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"I do not claim that the python implementation above is perfect, and certainly, there is much room for improvement. However, the factor of 50 difference in execution time is very spectacular. This is nothing but a consequence of the fact that the `ulab` functions run `C` code, with very little python overhead. The factor of 50 appears to be quite universal: the FFT routine obeys similar scaling (see [Speed of FFTs](#Speed-of-FFTs)), and this number came up with font rendering, too: [fast font rendering on graphical displays](https://forum.micropython.org/viewtopic.php?f=15&t=5815&p=33362&hilit=ufont#p33383)."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Comparison operators\n",
"\n",
"The smaller than, greater than, smaller or equal, and greater or equal operators return a vector of Booleans indicating the positions (`True`), where the condition is satisfied. "
]
},
{
"cell_type": "code",
"execution_count": 118,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-11T16:24:11.562136Z",
"start_time": "2019-10-11T16:24:11.548252Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"[True, True, True, True, False, False, False, False]\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array([1, 2, 3, 4, 5, 6, 7, 8], dtype=np.uint8)\n",
"print(a < 5)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"**WARNING:** Note that `numpy` returns an array of Booleans. For most use cases this fact should not make a difference. "
]
},
{
"cell_type": "code",
"execution_count": 119,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-11T16:25:30.888882Z",
"start_time": "2019-10-11T16:25:30.865641Z"
}
},
"outputs": [
{
"data": {
"text/plain": [
"array([ True, True, True, True, False, False, False, False])"
]
},
"execution_count": 119,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"a = array([1, 2, 3, 4, 5, 6, 7, 8])\n",
"a < 5"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"These operators work with matrices, too, in which case a list of lists of Booleans will be returned:"
]
},
{
"cell_type": "code",
"execution_count": 122,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-11T16:28:07.876371Z",
"start_time": "2019-10-11T16:28:07.859304Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"array([[0, 1, 2, 3, 4],\n",
"\t [1, 2, 3, 4, 5],\n",
"\t [2, 3, 4, 5, 6]], dtype=uint8)\n",
"[[True, True, True, True, True], [True, True, True, True, False], [True, True, True, False, False]]\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array([range(0, 5, 1), range(1, 6, 1), range(2, 7, 1)], dtype=np.uint8)\n",
"print(a)\n",
"print(a < 5)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Iterating over arrays\n",
"\n",
"`ndarray`s are iterable, which means that their elements can also be accessed as can the elements of a list, tuple, etc. If the array is one-dimensional, the iterator returns scalars, otherwise a new one-dimensional `ndarray`, which is simply a copy of the corresponding row of the matrix, i.e, its data type will be inherited."
]
},
{
"cell_type": "code",
"execution_count": 104,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-11T13:23:32.570237Z",
"start_time": "2019-10-11T13:23:32.550470Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a:\t array([1, 2, 3, 4, 5], dtype=uint8)\n",
"element 0 in a: 1\n",
"element 1 in a: 2\n",
"element 2 in a: 3\n",
"element 3 in a: 4\n",
"element 4 in a: 5\n",
"\n",
"b:\t array([[0, 1, 2, 3, 4],\n",
"\t [10, 11, 12, 13, 14],\n",
"\t [20, 21, 22, 23, 24],\n",
"\t [30, 31, 32, 33, 34]], dtype=uint8)\n",
"element 0 in b: array([0, 1, 2, 3, 4], dtype=uint8)\n",
"element 1 in b: array([10, 11, 12, 13, 14], dtype=uint8)\n",
"element 2 in b: array([20, 21, 22, 23, 24], dtype=uint8)\n",
"element 3 in b: array([30, 31, 32, 33, 34], dtype=uint8)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array([1, 2, 3, 4, 5], dtype=np.uint8)\n",
"b = np.array([range(5), range(10, 15, 1), range(20, 25, 1), range(30, 35, 1)], dtype=np.uint8)\n",
"\n",
"print(\"a:\\t\", a)\n",
"\n",
"for i, _a in enumerate(a):\n",
" print(\"element %d in a:\"%i, _a)\n",
" \n",
"print(\"\\nb:\\t\", b)\n",
"\n",
"for i, _b in enumerate(b):\n",
" print(\"element %d in b:\"%i, _b)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Slicing and indexing\n",
"\n",
"Copies of sub-arrays can be created by indexing, and slicing."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Indexing\n",
"\n",
"The simplest form of indexing is specifying a single integer between the square brackets as in "
]
},
{
"cell_type": "code",
"execution_count": 116,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-11T15:58:24.845241Z",
"start_time": "2019-10-11T15:58:24.821490Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a:\t\t\t\t\t\t array([0, 1, 2, ..., 7, 8, 9], dtype=uint8)\n",
"the first, and first from right element of a:\t 0 9\n",
"the second, and second from right element of a:\t 1 8\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array(range(10), dtype=np.uint8)\n",
"print(\"a:\\t\\t\\t\\t\\t\\t\", a)\n",
"print(\"the first, and first from right element of a:\\t\", a[0], a[-1])\n",
"print(\"the second, and second from right element of a:\\t\", a[1], a[-2])"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Indices are (not necessarily non-negative) integers, or a list of Booleans. By using a Boolean list, we can select those elements of an array that satisfy a specific condition. At the moment, such indexing is defined for row vectors only, for matrices the function raises a `ValueError` exception, though this will be rectified in a future version of `ulab`."
]
},
{
"cell_type": "code",
"execution_count": 183,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-11T18:51:18.226269Z",
"start_time": "2019-10-11T18:51:18.208328Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a:\t array([0.0, 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0], dtype=float)\n",
"a < 5:\t array([0.0, 1.0, 2.0, 3.0, 4.0], dtype=float)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array(range(9), dtype=np.float)\n",
"print(\"a:\\t\", a)\n",
"print(\"a < 5:\\t\", a[a < 5])"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Indexing with Boolean arrays can take more complicated expressions. This is a very concise way of comparing two vectors, e.g.:"
]
},
{
"cell_type": "code",
"execution_count": 523,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-20T06:50:17.527745Z",
"start_time": "2019-10-20T06:50:17.466717Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a:\t array([0, 1, 2, 3, 4, 5, 6, 7, 8], dtype=uint8)\n",
"\n",
"a**2:\t array([0, 1, 4, 9, 16, 25, 36, 49, 64], dtype=uint8)\n",
"\n",
"b:\t array([4, 4, 4, 3, 3, 3, 13, 13, 13], dtype=uint8)\n",
"\n",
"100*sin(b):\t array([-75.68025, -75.68025, -75.68025, 14.112, 14.112, 14.112, 42.01671, 42.01671, 42.01671], dtype=float)\n",
"\n",
"a[a*a > np.sin(b)*100.0]:\t array([0, 1, 2, 4, 5, 7, 8], dtype=uint8)\n",
"\n"
]
}
],
"source": [
"%%micropython -pyboard 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array(range(9), dtype=np.uint8)\n",
"b = np.array([4, 4, 4, 3, 3, 3, 13, 13, 13], dtype=np.uint8)\n",
"print(\"a:\\t\", a)\n",
"print(\"\\na**2:\\t\", a*a)\n",
"print(\"\\nb:\\t\", b)\n",
"print(\"\\n100*sin(b):\\t\", np.sin(b)*100.0)\n",
"print(\"\\na[a*a > np.sin(b)*100.0]:\\t\", a[a*a > np.sin(b)*100.0])"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Slicing and assigning to slices\n",
"\n",
"You can also generate sub-arrays by specifying slices as the index of an array. Slices are special python objects of the form \n",
"\n",
"```python\n",
"slice = start:end:stop\n",
"```\n",
"where `start`, `end`, and `stop` are (not necessarily non-negative) integers. Not all of these three numbers must be specified in an index, in fact, all three of them can be missing. The interpreter takes care of filling in the missing values. (Note that slices cannot be defined in this way, only there, where an index is expected.) For a good explanation on how slices work in python, you can read the stackoverflow question https://stackoverflow.com/questions/509211/understanding-slice-notation.\n",
"\n",
"Slices work on both axes:"
]
},
{
"cell_type": "code",
"execution_count": 129,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-29T19:50:23.551363Z",
"start_time": "2019-10-29T19:50:23.528032Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a:\n",
" array([[1, 2, 3],\n",
"\t [4, 5, 6],\n",
"\t [7, 8, 9]], dtype=uint8)\n",
"\n",
"a[0]:\n",
" array([1, 2, 3], dtype=uint8)\n",
"\n",
"a[0,:2]:\n",
" array([1, 2], dtype=uint8)\n",
"\n",
"a[:,0]:\n",
" array([1, 4, 7], dtype=uint8)\n",
"\n",
"a[-1]:\n",
" array([7, 8, 9], dtype=uint8)\n",
"\n",
"a[::1]:\n",
" array([[7, 8, 9],\n",
"\t [4, 5, 6]], dtype=uint8)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]], dtype=np.uint8)\n",
"print('a:\\n', a)\n",
"\n",
"# the first row\n",
"print('\\na[0]:\\n', a[0])\n",
"\n",
"# the first two elements of the first row\n",
"print('\\na[0,:2]:\\n', a[0,:2])\n",
"\n",
"# the zeroth element in each row (also known as the zeroth column)\n",
"print('\\na[:,0]:\\n', a[:,0])\n",
"\n",
"# the last but one row\n",
"print('\\na[-1]:\\n', a[-1])\n",
"\n",
"# the last two rows backwards\n",
"print('\\na[::1]:\\n', a[::-1])"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Assignment to slices can be done for the whole slice, per row, and per column. A couple of examples should make these statements clearer:"
]
},
{
"cell_type": "code",
"execution_count": 132,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-29T19:54:32.232516Z",
"start_time": "2019-10-29T19:54:32.217379Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a:\n",
" array([[0, 0, 0],\n",
"\t [0, 0, 0],\n",
"\t [0, 0, 0]], dtype=uint8)\n",
"\n",
"a[0] = 1\n",
" array([[1, 1, 1],\n",
"\t [0, 0, 0],\n",
"\t [0, 0, 0]], dtype=uint8)\n",
"\n",
"a[0] = np.array([1, 2, 3])\n",
" array([[1, 2, 179],\n",
"\t [0, 0, 0],\n",
"\t [0, 0, 0]], dtype=uint8)\n",
"\n",
"a[:,0]:\n",
" array([[1, 2, 3],\n",
"\t [0, 0, 3],\n",
"\t [0, 0, 3]], dtype=uint8)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"zero_list = [0, 0, 0]\n",
"a = np.array([zero_list, zero_list, zero_list], dtype=np.uint8)\n",
"print('a:\\n', a)\n",
"\n",
"# assigning to the whole row\n",
"a[0] = 1\n",
"print('\\na[0] = 1\\n', a)\n",
"\n",
"# assigning to the whole row\n",
"a[0] = np.array([1, 2, -333], dtype=np.float)\n",
"print('\\na[0] = np.array([1, 2, 3])\\n', a)\n",
"\n",
"# assigning to a column\n",
"a[:,2] = 3.0\n",
"print('\\na[:,0]:\\n', a)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Universal functions\n",
"\n",
"Standard mathematical functions can be calculated on any scalar-valued iterable (ranges, lists, tuples containing numbers), and on `ndarray`s without having to change the call signature. In all cases the functions return a new `ndarray` of typecode `float` (since these functions usually generate float values, anyway). The functions execute faster with `ndarray` arguments than with iterables, because the values of the input vector can be extracted faster. \n",
"\n",
"At present, the following functions are supported:\n",
"\n",
"`acos`, `acosh`, `asin`, `asinh`, `atan`, `atanh`, `ceil`, `cos`, `erf`, `erfc`, `exp`, `expm1`, `floor`, `tgamma`, `lgamma`, `log`, `log10`, `log2`, `sin`, `sinh`, `sqrt`, `tan`, `tanh`.\n",
"\n",
"These functions are applied element-wise to the arguments, thus, e.g., the exponential of a matrix cannot be calculated in this way."
]
},
{
"cell_type": "code",
"execution_count": 525,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-20T06:50:55.762508Z",
"start_time": "2019-10-20T06:50:55.655715Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a:\t range(0, 9)\n",
"exp(a):\t array([1.0, 2.718282, 7.389056, 20.08554, 54.59816, 148.4132, 403.4288, 1096.633, 2980.958], dtype=float)\n",
"\n",
"b:\t array([0.0, 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0], dtype=float)\n",
"exp(b):\t array([1.0, 2.718282, 7.389056, 20.08554, 54.59816, 148.4132, 403.4288, 1096.633, 2980.958], dtype=float)\n",
"\n",
"c:\t array([[1.0, 2.0, 3.0],\n",
"\t [4.0, 5.0, 6.0],\n",
"\t [7.0, 8.0, 9.0]], dtype=float)\n",
"exp(c):\t array([[2.718282, 7.389056, 20.08554],\n",
"\t [54.59816, 148.4132, 403.4288],\n",
"\t [1096.633, 2980.958, 8103.084]], dtype=float)\n",
"\n"
]
}
],
"source": [
"%%micropython -pyboard 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = range(9)\n",
"b = np.array(a)\n",
"\n",
"# works with ranges, lists, tuples etc.\n",
"print('a:\\t', a)\n",
"print('exp(a):\\t', np.exp(a))\n",
"\n",
"# with 1D arrays\n",
"print('\\nb:\\t', b)\n",
"print('exp(b):\\t', np.exp(b))\n",
"\n",
"# as well as with matrices\n",
"c = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]])\n",
"print('\\nc:\\t', c)\n",
"print('exp(c):\\t', np.exp(c))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Computation expenses\n",
"\n",
"The overhead for calculating with micropython iterables is quite significant: for the 1000 samples below, the difference is more than 800 microseconds, because internally the function has to create the `ndarray` for the output, has to fetch the iterable's items of unknown type, and then convert them to floats. All these steps are skipped for `ndarray`s, because these pieces of information are already known. "
]
},
{
"cell_type": "code",
"execution_count": 526,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-20T06:51:19.483011Z",
"start_time": "2019-10-20T06:51:19.444764Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"execution time: 1259 us\n",
"execution time: 408 us\n",
"\n"
]
}
],
"source": [
"%%micropython -pyboard 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = [0]*1000\n",
"b = np.array(a)\n",
"\n",
"@timeit\n",
"def measure_run_time(x):\n",
" return np.exp(x)\n",
"\n",
"measure_run_time(a)\n",
"\n",
"measure_run_time(b)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Of course, such a time saving is reasonable only, if the data are already available as an `ndarray`. If one has to initialise the `ndarray` from the list, then there is no gain, because the iterator was simply pushed into the initialisation function."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Numerical"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## linspace\n",
"\n",
"numpy: https://docs.scipy.org/doc/numpy/reference/generated/numpy.linspace.html\n",
"\n",
"This function returns an array, whose elements are uniformly spaced between the `start`, and `stop` points. The number of intervals is determined by the `num` keyword argument, whose default value is 50. With the `endpoint` keyword argument (defaults to `True`) one can include `stop` in the sequence. In addition, the `dtype` keyword can be supplied to force type conversion of the output. The default is `float`. Note that, when `dtype` is of integer type, the sequence is not necessarily evenly spaced. This is not an error, rather a consequence of rounding. (This is also the `numpy` behaviour.)"
]
},
{
"cell_type": "code",
"execution_count": 371,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-18T12:59:53.458138Z",
"start_time": "2019-10-18T12:59:53.438981Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"default sequence:\t array([0.0, 0.2040816396474838, 0.4081632792949677, ..., 9.591833114624023, 9.795914649963379, 9.999996185302734], dtype=float)\n",
"num=5:\t\t\t array([0.0, 2.5, 5.0, 7.5, 10.0], dtype=float)\n",
"num=5:\t\t\t array([0.0, 2.0, 4.0, 6.0, 8.0], dtype=float)\n",
"num=5:\t\t\t array([0, 0, 1, 2, 2, 3, 4], dtype=uint8)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"# generate a sequence with defaults\n",
"print('default sequence:\\t', np.linspace(0, 10))\n",
"\n",
"# num=5\n",
"print('num=5:\\t\\t\\t', np.linspace(0, 10, num=5))\n",
"\n",
"# num=5, endpoint=False\n",
"print('num=5:\\t\\t\\t', np.linspace(0, 10, num=5, endpoint=False))\n",
"\n",
"# num=5, endpoint=False, dtype=uint8\n",
"print('num=5:\\t\\t\\t', np.linspace(0, 5, num=7, endpoint=False, dtype=np.uint8))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## min, argmin, max, argmax\n",
"\n",
"numpy: https://docs.scipy.org/doc/numpy/reference/generated/numpy.min.html\n",
"\n",
"numpy: https://docs.scipy.org/doc/numpy/reference/generated/numpy.argmax.html\n",
"\n",
"numpy: https://docs.scipy.org/doc/numpy/reference/generated/numpy.max.html\n",
"\n",
"numpy: https://docs.scipy.org/doc/numpy/reference/generated/numpy.argmax.html\n",
"\n",
"**WARNING:** Difference to `numpy`: the `out` keyword argument is not implemented.\n",
"\n",
"These functions follow the same pattern, and work with generic iterables, and `ndarray`s. `min`, and `max` return the minimum or maximum of a sequence. If the input array is two-dimensional, the `axis` keyword argument can be supplied, in which case the minimum/maximum along the given axis will be returned. If `axis=None` (this is also the default value), the minimum/maximum of the flattened array will be determined.\n",
"\n",
"`argmin/argmax` return the position (index) of the minimum/maximum in the sequence."
]
},
{
"cell_type": "code",
"execution_count": 375,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-18T13:08:28.113525Z",
"start_time": "2019-10-18T13:08:28.093518Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a: array([1.0, 2.0, 0.0, 1.0, 10.0], dtype=float)\n",
"min of a: 0.0\n",
"argmin of a: 2\n",
"\n",
"b:\n",
" array([[1.0, 2.0, 0.0],\n",
"\t [1.0, 10.0, -1.0]], dtype=float)\n",
"min of b (flattened): -1.0\n",
"min of b (axis=0): array([1.0, 2.0, -1.0], dtype=float)\n",
"min of b (axis=1): array([0.0, -1.0], dtype=float)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array([1, 2, 0, 1, 10])\n",
"print('a:', a)\n",
"print('min of a:', np.min(a))\n",
"print('argmin of a:', np.argmin(a))\n",
"\n",
"b = np.array([[1, 2, 0], [1, 10, -1]])\n",
"print('\\nb:\\n', b)\n",
"print('min of b (flattened):', np.min(b))\n",
"print('min of b (axis=0):', np.min(b, axis=0))\n",
"print('min of b (axis=1):', np.min(b, axis=1))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## sum, std, mean\n",
"\n",
"numpy: https://docs.scipy.org/doc/numpy/reference/generated/numpy.sum.html\n",
"\n",
"numpy: https://docs.scipy.org/doc/numpy/reference/generated/numpy.std.html\n",
"\n",
"numpy: https://docs.scipy.org/doc/numpy/reference/generated/numpy.mean.html\n",
"\n",
"These three functions follow the same pattern: if the axis keyword is not specified, it assumes the default value of `None`, and returns the result of the computation for the flattened array. Otherwise, the calculation is along the given axis."
]
},
{
"cell_type": "code",
"execution_count": 527,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-20T06:51:58.845076Z",
"start_time": "2019-10-20T06:51:58.798730Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a: \n",
" array([[1.0, 2.0, 3.0],\n",
"\t [4.0, 5.0, 6.0],\n",
"\t [7.0, 8.0, 9.0]], dtype=float)\n",
"sum, flat array: 45.0\n",
"mean, horizontal: array([2.0, 5.0, 8.0], dtype=float)\n",
"std, vertical: array([2.44949, 2.44949, 2.44949], dtype=float)\n",
"\n"
]
}
],
"source": [
"%%micropython -pyboard 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]])\n",
"print('a: \\n', a)\n",
"\n",
"print('sum, flat array: ', np.sum(a))\n",
"\n",
"print('mean, horizontal: ', np.mean(a, axis=1))\n",
"\n",
"print('std, vertical: ', np.std(a, axis=0))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## roll\n",
"\n",
"numpy: https://docs.scipy.org/doc/numpy/reference/generated/numpy.roll.html\n",
"\n",
"The roll function shifts the content of a vector by the positions given as the second argument. If the `axis` keyword is supplied, the shift is applied to the given axis."
]
},
{
"cell_type": "code",
"execution_count": 229,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-11T19:39:47.459395Z",
"start_time": "2019-10-11T19:39:47.443691Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a:\t\t\t array([1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0], dtype=float)\n",
"a rolled to the left:\t array([3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 1.0, 2.0], dtype=float)\n",
"a rolled to the right:\t array([1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0], dtype=float)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array([1, 2, 3, 4, 5, 6, 7, 8])\n",
"print(\"a:\\t\\t\\t\", a)\n",
"\n",
"np.roll(a, 2)\n",
"print(\"a rolled to the left:\\t\", a)\n",
"\n",
"# this should be the original vector\n",
"np.roll(a, -2)\n",
"print(\"a rolled to the right:\\t\", a)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Rolling works with matrices, too. If the `axis` keyword is 0, the matrix is rolled along its vertical axis, otherwise, horizontally. \n",
"\n",
"Horizontal rolls are faster, because they require fewer steps, and larger memory chunks are copied, however, they also require more RAM: basically the whole row must be stored internally. Most expensive are the `None` keyword values, because with `axis = None`, the array is flattened first, hence the row's length is the size of the whole matrix.\n",
"\n",
"Vertical rolls require two internal copies of single columns. "
]
},
{
"cell_type": "code",
"execution_count": 268,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-15T17:46:20.051069Z",
"start_time": "2019-10-15T17:46:20.033205Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a:\n",
" array([[1.0, 2.0, 3.0, 4.0],\n",
"\t [5.0, 6.0, 7.0, 8.0]], dtype=float)\n",
"\n",
"a rolled to the left:\n",
" array([[3.0, 4.0, 5.0, 6.0],\n",
"\t [7.0, 8.0, 1.0, 2.0]], dtype=float)\n",
"\n",
"a rolled up:\n",
" array([[6.0, 3.0, 4.0, 5.0],\n",
"\t [2.0, 7.0, 8.0, 1.0]], dtype=float)\n",
"\n",
"a rolled with None:\n",
" array([[3.0, 4.0, 5.0, 2.0],\n",
"\t [7.0, 8.0, 1.0, 6.0]], dtype=float)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array([[1, 2, 3, 4], [5, 6, 7, 8]])\n",
"print(\"a:\\n\", a)\n",
"\n",
"np.roll(a, 2)\n",
"print(\"\\na rolled to the left:\\n\", a)\n",
"\n",
"np.roll(a, -1, axis=1)\n",
"print(\"\\na rolled up:\\n\", a)\n",
"\n",
"np.roll(a, 1, axis=None)\n",
"print(\"\\na rolled with None:\\n\", a)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Simple running weighted average\n",
"\n",
"As a demonstration of the conciseness of `ulab/numpy` operations, we will calculate an exponentially weighted running average of a measurement vector in just a couple of lines. I chose this particular example, because I think that this can indeed be used in real-life applications."
]
},
{
"cell_type": "code",
"execution_count": 230,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-11T20:03:00.713235Z",
"start_time": "2019-10-11T20:03:00.696932Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"array([0.01165623031556606, 0.03168492019176483, 0.08612854033708572, 0.234121635556221, 0.6364086270332336], dtype=float)\n",
"0.2545634508132935\n",
"array([0.0, 0.0, 0.0, 0.0, 2.0], dtype=float)\n",
"0.3482121050357819\n",
"array([0.0, 0.0, 0.0, 2.0, 2.0], dtype=float)\n",
"0.3826635211706161\n",
"array([0.0, 0.0, 2.0, 2.0, 2.0], dtype=float)\n",
"0.3953374892473221\n",
"array([0.0, 2.0, 2.0, 2.0, 2.0], dtype=float)\n",
"0.3999999813735485\n",
"array([2.0, 2.0, 2.0, 2.0, 2.0], dtype=float)\n",
"0.3999999813735485\n",
"array([2.0, 2.0, 2.0, 2.0, 2.0], dtype=float)\n",
"0.3999999813735485\n",
"array([2.0, 2.0, 2.0, 2.0, 2.0], dtype=float)\n",
"0.3999999813735485\n",
"array([2.0, 2.0, 2.0, 2.0, 2.0], dtype=float)\n",
"0.3999999813735485\n",
"array([2.0, 2.0, 2.0, 2.0, 2.0], dtype=float)\n",
"0.3999999813735485\n",
"array([2.0, 2.0, 2.0, 2.0, 2.0], dtype=float)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"def dummy_adc():\n",
" # dummy adc function, so that the results are reproducible\n",
" return 2\n",
" \n",
"n = 10\n",
"# These are the normalised weights; the last entry is the most dominant\n",
"weight = np.exp([1, 2, 3, 4, 5])\n",
"weight = weight/np.sum(weight)\n",
"\n",
"print(weight)\n",
"# initial array of samples\n",
"samples = np.array([0]*n)\n",
"\n",
"for i in range(n):\n",
" # a new datum is inserted on the right hand side. This simply overwrites whatever was in the last slot\n",
" samples[-1] = dummy_adc()\n",
" print(np.mean(samples[-5:]*weight))\n",
" print(samples[-5:])\n",
" # the data are shifted by one position to the left\n",
" np.roll(samples, 1)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## flip\n",
"\n",
"numpy: https://docs.scipy.org/doc/numpy/reference/generated/numpy.flip.html\n",
"\n",
"The `flip` function takes one positional, an `ndarray`, and one keyword argument, `axis = None`, and reverses the order of elements along the given axis. If the keyword argument is `None`, the matrix' entries are flipped along all axes. `flip` returns a new copy of the array."
]
},
{
"cell_type": "code",
"execution_count": 275,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-16T06:35:52.163725Z",
"start_time": "2019-10-16T06:35:52.149231Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a: \t array([1.0, 2.0, 3.0, 4.0, 5.0], dtype=float)\n",
"a flipped:\t array([5.0, 4.0, 3.0, 2.0, 1.0], dtype=float)\n",
"\n",
"a flipped horizontally\n",
" array([[3, 2, 1],\n",
"\t [6, 5, 4],\n",
"\t [9, 8, 7]], dtype=uint8)\n",
"\n",
"a flipped vertically\n",
" array([[7, 8, 9],\n",
"\t [4, 5, 6],\n",
"\t [1, 2, 3]], dtype=uint8)\n",
"\n",
"a flipped horizontally+vertically\n",
" array([[9, 8, 7],\n",
"\t [6, 5, 4],\n",
"\t [3, 2, 1]], dtype=uint8)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array([1, 2, 3, 4, 5])\n",
"print(\"a: \\t\", a)\n",
"print(\"a flipped:\\t\", np.flip(a))\n",
"\n",
"a = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]], dtype=np.uint8)\n",
"print(\"\\na flipped horizontally\\n\", np.flip(a, axis=1))\n",
"print(\"\\na flipped vertically\\n\", np.flip(a, axis=0))\n",
"print(\"\\na flipped horizontally+vertically\\n\", np.flip(a))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## diff\n",
"\n",
"numpy: https://docs.scipy.org/doc/numpy/reference/generated/numpy.diff.html\n",
"\n",
"The `diff` function returns the numerical derivative of the forward scheme, or more accurately, the differences of an `ndarray` along a given axis. The order of derivative can be stipulated with the `n` keyword argument, which should be between 0, and 9. Default is 1. If higher order derivatives are required, they can be gotten by repeated calls to the function. The `axis` keyword argument should be -1 (last axis, in `ulab` equivalent to the second axis, and this also happens to be the default value), 0, or 1. \n",
"\n",
"Beyond the output array, the function requires only a couple of bytes of extra RAM for the differentiation stencil. (The stencil is an `int8` array, one byte longer than `n`. This also explains, why the highest order is 9: the coefficients of a ninth-order stencil all fit in signed bytes, while 10 would require `int16`.) Note that as usual in numerical differentiation (and also in `numpy`), the length of the respective axis will be reduced by `n` after the operation. If `n` is larger than, or equal to the length of the axis, an empty array will be returned.\n",
"\n",
"**WARNING**: the `diff` function does not implement the `prepend` and `append` keywords that can be found in `numpy`. "
]
},
{
"cell_type": "code",
"execution_count": 169,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-31T11:51:02.854338Z",
"start_time": "2019-10-31T11:51:02.838000Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a:\n",
" array([0, 1, 2, 3, 4, 5, 6, 7, 8], dtype=uint8)\n",
"\n",
"first derivative:\n",
" array([1, 1, 1, 1, 1, 1, 1, 1], dtype=uint8)\n",
"\n",
"second derivative:\n",
" array([0, 0, 0, 0, 0, 0, 0], dtype=uint8)\n",
"\n",
"c:\n",
" array([[1.0, 2.0, 3.0, 4.0],\n",
"\t [4.0, 3.0, 2.0, 1.0],\n",
"\t [1.0, 4.0, 9.0, 16.0],\n",
"\t [0.0, 0.0, 0.0, 0.0]], dtype=float)\n",
"\n",
"first derivative, first axis:\n",
" array([[3.0, 1.0, -1.0, -3.0],\n",
"\t [-3.0, 1.0, 7.0, 15.0],\n",
"\t [-1.0, -4.0, -9.0, -16.0]], dtype=float)\n",
"\n",
"first derivative, second axis:\n",
" array([[1.0, 1.0, 1.0],\n",
"\t [-1.0, -1.0, -1.0],\n",
"\t [3.0, 5.0, 7.0],\n",
"\t [0.0, 0.0, 0.0]], dtype=float)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array(range(9), dtype=np.uint8)\n",
"print('a:\\n', a)\n",
"\n",
"print('\\nfirst derivative:\\n', np.diff(a, n=1))\n",
"print('\\nsecond derivative:\\n', np.diff(a, n=2))\n",
"\n",
"c = np.array([[1, 2, 3, 4], [4, 3, 2, 1], [1, 4, 9, 16], [0, 0, 0, 0]])\n",
"print('\\nc:\\n', c)\n",
"print('\\nfirst derivative, first axis:\\n', np.diff(c, axis=0))\n",
"print('\\nfirst derivative, second axis:\\n', np.diff(c, axis=1))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## sort\n",
"\n",
"numpy: https://docs.scipy.org/doc/numpy/reference/generated/numpy.sort.html\n",
"\n",
"The sort function takes an ndarray, and sorts its elements in ascending order along the specified axis using a heap sort algorithm. As opposed to the `.sort()` method discussed earlier, this function creates a copy of its input before sorting, and at the end, returns this copy. Sorting takes place in place, without auxiliary storage. The `axis` keyword argument takes on the possible values of -1 (the last axis, in `ulab` equivalent to the second axis, and this also happens to be the default value), 0, 1, or `None`. The first three cases are identical to those in [diff](#diff), while the last one flattens the array before sorting. \n",
"\n",
"If descending order is required, the result can simply be `flip`ped, see [flip](#flip).\n",
"\n",
"**WARNING:** `numpy` defines the `kind`, and `order` keyword arguments that are not implemented here. The function in `ulab` always uses heap sort, and since `ulab` does not have the concept of data fields, the `order` keyword argument would have no meaning."
]
},
{
"cell_type": "code",
"execution_count": 29,
"metadata": {
"ExecuteTime": {
"end_time": "2019-11-05T16:06:27.536193Z",
"start_time": "2019-11-05T16:06:27.521792Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"\n",
"a:\n",
" array([[1.0, 12.0, 3.0, 0.0],\n",
"\t [5.0, 3.0, 4.0, 1.0],\n",
"\t [9.0, 11.0, 1.0, 8.0],\n",
"\t [7.0, 10.0, 0.0, 1.0]], dtype=float)\n",
"\n",
"a sorted along vertical axis:\n",
" array([[1.0, 3.0, 0.0, 0.0],\n",
"\t [5.0, 10.0, 1.0, 1.0],\n",
"\t [7.0, 11.0, 3.0, 1.0],\n",
"\t [9.0, 12.0, 4.0, 8.0]], dtype=float)\n",
"\n",
"a sorted along horizontal axis:\n",
" array([[0.0, 1.0, 3.0, 12.0],\n",
"\t [1.0, 3.0, 4.0, 5.0],\n",
"\t [1.0, 8.0, 9.0, 11.0],\n",
"\t [0.0, 1.0, 7.0, 10.0]], dtype=float)\n",
"\n",
"flattened a sorted:\n",
" array([0.0, 0.0, 1.0, ..., 10.0, 11.0, 12.0], dtype=float)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array([[1, 12, 3, 0], [5, 3, 4, 1], [9, 11, 1, 8], [7, 10, 0, 1]], dtype=np.float)\n",
"print('\\na:\\n', a)\n",
"b = np.sort(a, axis=0)\n",
"print('\\na sorted along vertical axis:\\n', b)\n",
"\n",
"c = np.sort(a, axis=1)\n",
"print('\\na sorted along horizontal axis:\\n', c)\n",
"\n",
"c = np.sort(a, axis=None)\n",
"print('\\nflattened a sorted:\\n', c)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Heap sort requires $\\sim N\\log N$ operations, and notably, the worst case costs only 20% more time than the average. In order to get an order-of-magnitude estimate, we will take the sine of 1000 uniformly spaced numbers between 0, and two pi, and sort them:"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"%%micropython -pyboard 1\n",
"\n",
"import ulab as np\n",
"\n",
"@timeit\n",
"def sort_time(array):\n",
" return np.sort(array)\n",
"\n",
"b = np.sin(np.linspace(0, 6.28, num=1000))\n",
"print('b: ', b)\n",
"sort_time(b)\n",
"print('\\nb sorted:\\n', b)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## argsort\n",
"\n",
"numpy: https://docs.scipy.org/doc/numpy/reference/generated/numpy.argsort.html\n",
"\n",
"Similarly to [sort](#sort), `argsort` takes a positional, and a keyword argument, and returns an unsigned short index array of type `ndarray` with the same dimensions as the input, or, if `axis=None`, as a row vector with length equal to the number of elements in the input (i.e., the flattened array). The indices in the output sort the input in ascending order. The routine in `argsort` is the same as in `sort`, therefore, the comments on computational expenses (time and RAM) also apply. In particular, since no copy of the original data is required, virtually no RAM beyond the output array is used. \n",
"\n",
"Since the underlying container of the output array is of type `uint16_t`, neither of the output dimensions should be larger than 65535."
]
},
{
"cell_type": "code",
"execution_count": 30,
"metadata": {
"ExecuteTime": {
"end_time": "2019-11-06T06:28:45.719578Z",
"start_time": "2019-11-06T06:28:45.704072Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"\n",
"a:\n",
" array([[1.0, 12.0, 3.0, 0.0],\n",
"\t [5.0, 3.0, 4.0, 1.0],\n",
"\t [9.0, 11.0, 1.0, 8.0],\n",
"\t [7.0, 10.0, 0.0, 1.0]], dtype=float)\n",
"\n",
"a sorted along vertical axis:\n",
" array([[0, 1, 3, 0],\n",
"\t [1, 3, 2, 1],\n",
"\t [3, 2, 0, 3],\n",
"\t [2, 0, 1, 2]], dtype=uint16)\n",
"\n",
"a sorted along horizontal axis:\n",
" array([[3, 0, 2, 1],\n",
"\t [3, 1, 2, 0],\n",
"\t [2, 3, 0, 1],\n",
"\t [2, 3, 0, 1]], dtype=uint16)\n",
"\n",
"flattened a sorted:\n",
" array([3, 14, 0, ..., 13, 9, 1], dtype=uint16)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array([[1, 12, 3, 0], [5, 3, 4, 1], [9, 11, 1, 8], [7, 10, 0, 1]], dtype=np.float)\n",
"print('\\na:\\n', a)\n",
"b = np.argsort(a, axis=0)\n",
"print('\\na sorted along vertical axis:\\n', b)\n",
"\n",
"c = np.argsort(a, axis=1)\n",
"print('\\na sorted along horizontal axis:\\n', c)\n",
"\n",
"c = np.argsort(a, axis=None)\n",
"print('\\nflattened a sorted:\\n', c)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Since during the sorting, only the indices are shuffled, `argsort` does not modify the input array, as one can verify this by the following example:"
]
},
{
"cell_type": "code",
"execution_count": 33,
"metadata": {
"ExecuteTime": {
"end_time": "2019-11-06T16:04:31.653444Z",
"start_time": "2019-11-06T16:04:31.634995Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"\n",
"a:\n",
" array([0, 5, 1, 3, 2, 4], dtype=uint8)\n",
"\n",
"sorting indices:\n",
" array([0, 2, 4, 3, 5, 1], dtype=uint16)\n",
"\n",
"the original array:\n",
" array([0, 5, 1, 3, 2, 4], dtype=uint8)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array([0, 5, 1, 3, 2, 4], dtype=np.uint8)\n",
"print('\\na:\\n', a)\n",
"b = np.argsort(a, axis=1)\n",
"print('\\nsorting indices:\\n', b)\n",
"print('\\nthe original array:\\n', a)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Linalg"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## size\n",
"\n",
"`size` takes a single argument, the axis, whose size is to be returned. Depending on the value of the argument, the following information will be returned:\n",
"\n",
"1. argument is 0: the number of elements of the array\n",
"2. argument is 1: the number of rows\n",
"3. argument is 2: the number of columns"
]
},
{
"cell_type": "code",
"execution_count": 321,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-16T17:47:38.439137Z",
"start_time": "2019-10-16T17:47:38.422003Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a:\n",
" array([1, 2, 3, 4], dtype=int8)\n",
"size of a: 4 , 4\n",
"\n",
"b:\n",
" array([[1, 2],\n",
"\t [3, 4]], dtype=int8)\n",
"size of b: 4 , 2 , 2\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array([1, 2, 3, 4], dtype=np.int8)\n",
"print(\"a:\\n\", a)\n",
"print(\"size of a:\", np.size(a, axis=None), \",\", np.size(a, axis=0))\n",
"\n",
"b= np.array([[1, 2], [3, 4]], dtype=np.int8)\n",
"print(\"\\nb:\\n\", b)\n",
"print(\"size of b:\", np.size(b, axis=None), \",\", np.size(b, axis=0), \",\", np.size(b, axis=1))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## ones, zeros\n",
"\n",
"numpy: https://docs.scipy.org/doc/numpy/reference/generated/numpy.zeros.html\n",
"\n",
"numpy: https://docs.scipy.org/doc/numpy/reference/generated/numpy.ones.html\n",
"\n",
"A couple of special arrays and matrices can easily be initialised by calling one of the `ones`, or `zeros` functions. `ones` and `zeros` follow the same pattern, and have the call signature\n",
"\n",
"```python\n",
"ones(shape, dtype=float)\n",
"zeros(shape, dtype=float)\n",
"```\n",
"where shape is either an integer, or a 2-tuple."
]
},
{
"cell_type": "code",
"execution_count": 14,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-10T17:15:39.366695Z",
"start_time": "2019-10-10T17:15:39.344152Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"array([1, 1, 1, 1, 1, 1], dtype=uint8)\n",
"array([[0.0, 0.0, 0.0, 0.0],\n",
"\t [0.0, 0.0, 0.0, 0.0],\n",
"\t [0.0, 0.0, 0.0, 0.0],\n",
"\t [0.0, 0.0, 0.0, 0.0],\n",
"\t [0.0, 0.0, 0.0, 0.0],\n",
"\t [0.0, 0.0, 0.0, 0.0]], dtype=float)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"print(np.ones(6, dtype=np.uint8))\n",
"print(np.zeros((6, 4)))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## eye\n",
"\n",
"numpy: https://docs.scipy.org/doc/numpy/reference/generated/numpy.eye.html\n",
"\n",
"Another special array method is the `eye` function, whose call signature is \n",
"\n",
"```python\n",
"eye(N, M, k=0, dtype=float)\n",
"```\n",
"where `N` (`M`) specify the dimensions of the matrix (if only `N` is supplied, then we get a square matrix, otherwise one with `M` rows, and `N` columns), and `k` is the shift of the ones (the main diagonal corresponds to `k=0`). Here are a couple of examples."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### With a single argument"
]
},
{
"cell_type": "code",
"execution_count": 23,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-10T17:22:21.351732Z",
"start_time": "2019-10-10T17:22:21.336206Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"array([[1.0, 0.0, 0.0, 0.0, 0.0],\n",
"\t [0.0, 1.0, 0.0, 0.0, 0.0],\n",
"\t [0.0, 0.0, 1.0, 0.0, 0.0],\n",
"\t [0.0, 0.0, 0.0, 1.0, 0.0],\n",
"\t [0.0, 0.0, 0.0, 0.0, 1.0]], dtype=float)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"print(np.eye(5))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Specifying the dimensions of the matrix"
]
},
{
"cell_type": "code",
"execution_count": 21,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-10T17:19:08.627749Z",
"start_time": "2019-10-10T17:19:08.612026Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"array([[1, 0, 0, 0],\n",
"\t [0, 1, 0, 0],\n",
"\t [0, 0, 1, 0],\n",
"\t [0, 0, 0, 1],\n",
"\t [0, 0, 0, 0],\n",
"\t [0, 0, 0, 0]], dtype=int8)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"print(np.eye(4, M=6, dtype=np.int8))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Shifting the diagonal"
]
},
{
"cell_type": "code",
"execution_count": 22,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-10T17:19:31.618997Z",
"start_time": "2019-10-10T17:19:31.600377Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"array([[0, 0, 0, 0],\n",
"\t [1, 0, 0, 0],\n",
"\t [0, 1, 0, 0],\n",
"\t [0, 0, 1, 0],\n",
"\t [0, 0, 0, 1],\n",
"\t [0, 0, 0, 0]], dtype=int16)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"print(np.eye(4, M=6, k=-1, dtype=np.int16))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## inv\n",
"\n",
"A square matrix, provided that it is not singular, can be inverted by calling the `inv` function that takes a single argument. The inversion is based on successive elimination of elements in the lower left triangle, and raises a `ValueError` exception, if the matrix turns out to be singular (i.e., one of the diagonal entries is zero)."
]
},
{
"cell_type": "code",
"execution_count": 538,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-20T06:55:25.025726Z",
"start_time": "2019-10-20T06:55:24.982557Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"array([[-2.166666, 1.499999, -0.8333326, 1.0],\n",
"\t [1.666666, -3.333331, 1.666666, -4.768516e-08],\n",
"\t [0.1666672, 2.166666, -0.8333327, -1.0],\n",
"\t [-0.1666666, -0.3333334, 4.96705e-08, 0.5]], dtype=float)\n",
"\n"
]
}
],
"source": [
"%%micropython -pyboard 1\n",
"\n",
"import ulab as np\n",
"\n",
"m = np.array([[1, 2, 3, 4], [4, 5, 6, 4], [7, 8.6, 9, 4], [3, 4, 5, 6]])\n",
"\n",
"print(np.inv(m))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Computation expenses\n",
"\n",
"Note that the cost of inverting a matrix is approximately twice as many floats (RAM), as the number of entries in the original matrix, and approximately as many operations, as the number of entries. Here are a couple of numbers: "
]
},
{
"cell_type": "code",
"execution_count": 552,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-20T07:10:39.190734Z",
"start_time": "2019-10-20T07:10:39.138872Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"2 by 2 matrix:\n",
"execution time: 65 us\n",
"\n",
"4 by 4 matrix:\n",
"execution time: 105 us\n",
"\n",
"8 by 8 matrix:\n",
"execution time: 299 us\n",
"\n"
]
}
],
"source": [
"%%micropython -pyboard 1\n",
"\n",
"import ulab as np\n",
"\n",
"@timeit\n",
"def invert_matrix(m):\n",
" return np.inv(m)\n",
"\n",
"m = np.array([[1, 2,], [4, 5]])\n",
"print('2 by 2 matrix:')\n",
"invert_matrix(m)\n",
"\n",
"m = np.array([[1, 2, 3, 4], [4, 5, 6, 4], [7, 8.6, 9, 4], [3, 4, 5, 6]])\n",
"print('\\n4 by 4 matrix:')\n",
"invert_matrix(m)\n",
"\n",
"m = np.array([[1, 2, 3, 4, 5, 6, 7, 8], [0, 5, 6, 4, 5, 6, 4, 5], \n",
" [0, 0, 9, 7, 8, 9, 7, 8], [0, 0, 0, 10, 11, 12, 11, 12], \n",
" [0, 0, 0, 0, 4, 6, 7, 8], [0, 0, 0, 0, 0, 5, 6, 7], \n",
" [0, 0, 0, 0, 0, 0, 7, 6], [0, 0, 0, 0, 0, 0, 0, 2]])\n",
"print('\\n8 by 8 matrix:')\n",
"invert_matrix(m)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The above-mentioned scaling is not obeyed strictly. The reason for the discrepancy is that the function call is still the same for all three cases: the input must be inspected, the output array must be created, and so on. "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## dot\n",
"\n",
"numpy: https://docs.scipy.org/doc/numpy/reference/generated/numpy.dot.html\n",
"\n",
"\n",
"**WARNING:** numpy applies upcasting rules for the multiplication of matrices, while `ulab` simply returns a float matrix. \n",
"\n",
"Once you can invert a matrix, you might want to know, whether the inversion is correct. You can simply take the original matrix and its inverse, and multiply them by calling the `dot` function, which takes the two matrices as its arguments. If the matrix dimensions do not match, the function raises a `ValueError`. The result of the multiplication is expected to be the unit matrix, which is demonstrated below."
]
},
{
"cell_type": "code",
"execution_count": 556,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-20T07:13:30.102776Z",
"start_time": "2019-10-20T07:13:30.073704Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"m:\n",
" array([[1, 2, 3],\n",
"\t [4, 5, 6],\n",
"\t [7, 10, 9]], dtype=uint8)\n",
"\n",
"m^-1:\n",
" array([[-1.25, 1.0, -0.25],\n",
"\t [0.5, -1.0, 0.5],\n",
"\t [0.4166667, 0.3333334, -0.25]], dtype=float)\n",
"\n",
"m*m^-1:\n",
" array([[1.0, 2.384186e-07, -1.490116e-07],\n",
"\t [-2.980232e-07, 1.000001, -4.172325e-07],\n",
"\t [-3.278255e-07, 1.311302e-06, 0.9999992]], dtype=float)\n",
"\n"
]
}
],
"source": [
"%%micropython -pyboard 1\n",
"\n",
"import ulab as np\n",
"\n",
"m = np.array([[1, 2, 3], [4, 5, 6], [7, 10, 9]], dtype=np.uint8)\n",
"n = np.inv(m)\n",
"print(\"m:\\n\", m)\n",
"print(\"\\nm^-1:\\n\", n)\n",
"# this should be the unit matrix\n",
"print(\"\\nm*m^-1:\\n\", np.dot(m, n))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Note that for matrix multiplication you don't necessarily need square matrices, it is enough, if their dimensions are compatible (i.e., the the left-hand-side matrix has as many columns, as does the right-hand-side matrix rows):"
]
},
{
"cell_type": "code",
"execution_count": 37,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-10T17:33:17.921324Z",
"start_time": "2019-10-10T17:33:17.900587Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"array([[1, 2, 3, 4],\n",
"\t [5, 6, 7, 8]], dtype=uint8)\n",
"array([[1, 2],\n",
"\t [3, 4],\n",
"\t [5, 6],\n",
"\t [7, 8]], dtype=uint8)\n",
"array([[7.0, 10.0],\n",
"\t [23.0, 34.0]], dtype=float)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"m = np.array([[1, 2, 3, 4], [5, 6, 7, 8]], dtype=np.uint8)\n",
"n = np.array([[1, 2], [3, 4], [5, 6], [7, 8]], dtype=np.uint8)\n",
"print(m)\n",
"print(n)\n",
"print(np.dot(m, n))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## det\n",
"\n",
"numpy: https://docs.scipy.org/doc/numpy/reference/generated/numpy.linalg.det.html\n",
"\n",
"The `det` function takes a square matrix as its single argument, and calculates the determinant. The calculation is based on successive elimination of the matrix elements, and the return value is a float, even if the input array was of integer type."
]
},
{
"cell_type": "code",
"execution_count": 495,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-19T13:27:24.246995Z",
"start_time": "2019-10-19T13:27:24.228698Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"-2.0\n",
"\n"
]
}
],
"source": [
"%%micropython -pyboard 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array([[1, 2], [3, 4]], dtype=np.uint8)\n",
"print(np.det(a))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Benchmark\n",
"\n",
"Since the routine for calculating the determinant is pretty much the same as for finding the [inverse of a matrix](#inv), the execution times are similar:"
]
},
{
"cell_type": "code",
"execution_count": 557,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-20T07:14:59.778987Z",
"start_time": "2019-10-20T07:14:59.740021Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"execution time: 294 us\n",
"\n"
]
}
],
"source": [
"%%micropython -pyboard 1\n",
"\n",
"@timeit\n",
"def matrix_det(m):\n",
" return np.inv(m)\n",
"\n",
"m = np.array([[1, 2, 3, 4, 5, 6, 7, 8], [0, 5, 6, 4, 5, 6, 4, 5], \n",
" [0, 0, 9, 7, 8, 9, 7, 8], [0, 0, 0, 10, 11, 12, 11, 12], \n",
" [0, 0, 0, 0, 4, 6, 7, 8], [0, 0, 0, 0, 0, 5, 6, 7], \n",
" [0, 0, 0, 0, 0, 0, 7, 6], [0, 0, 0, 0, 0, 0, 0, 2]])\n",
"\n",
"matrix_det(m)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## eig\n",
"\n",
"numpy: https://docs.scipy.org/doc/numpy/reference/generated/numpy.linalg.eig.html\n",
"\n",
"The `eig` function calculates the eigenvalues and the eigenvectors of a real, symmetric square matrix. If the matrix is not symmetric, a `ValueError` will be raised. The function takes a single argument, and returns a tuple with the eigenvalues, and eigenvectors. With the help of the eigenvectors, amongst other things, you can implement sophisticated stabilisation routines for robots."
]
},
{
"cell_type": "code",
"execution_count": 496,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-19T13:27:30.977472Z",
"start_time": "2019-10-19T13:27:30.943326Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"eigenvectors of a:\n",
" array([-1.165288, 0.8029362, 5.585626, 13.77673], dtype=float)\n",
"\n",
"eigenvalues of a:\n",
" array([[0.8151754, -0.4499267, -0.1643907, 0.3256237],\n",
"\t [0.2211193, 0.7847154, 0.08373602, 0.5729892],\n",
"\t [-0.1340859, -0.3100657, 0.8742685, 0.3486182],\n",
"\t [-0.5182822, -0.2926556, -0.4490192, 0.6664218]], dtype=float)\n",
"\n"
]
}
],
"source": [
"%%micropython -pyboard 1\n",
"\n",
"import ulab as np\n",
"\n",
"a = np.array([[1, 2, 1, 4], [2, 5, 3, 5], [1, 3, 6, 1], [4, 5, 1, 7]], dtype=np.uint8)\n",
"x, y = np.eig(a)\n",
"print('eigenvectors of a:\\n', x)\n",
"print('\\neigenvalues of a:\\n', y)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The same matrix diagonalised with `numpy` yields:"
]
},
{
"cell_type": "code",
"execution_count": 389,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-19T09:47:46.835349Z",
"start_time": "2019-10-19T09:47:46.785592Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"eigenvectors of a:\n",
" [13.77672606 -1.16528837 0.80293655 5.58562576]\n",
"\n",
"eigenvalues of a:\n",
" [[ 0.32561419 0.815156 0.44994112 -0.16446602]\n",
" [ 0.57300777 0.22113342 -0.78469926 0.08372081]\n",
" [ 0.34861093 -0.13401142 0.31007764 0.87427868]\n",
" [ 0.66641421 -0.51832581 0.29266348 -0.44897499]]\n"
]
}
],
"source": [
"a = array([[1, 2, 1, 4], [2, 5, 3, 5], [1, 3, 6, 1], [4, 5, 1, 7]], dtype=np.uint8)\n",
"x, y = eig(a)\n",
"print('eigenvectors of a:\\n', x)\n",
"print('\\neigenvalues of a:\\n', y)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"When comparing results, we should keep two things in mind: \n",
"\n",
"1. the eigenvalues and eigenvectors are not necessarily sorted in the same way\n",
"2. an eigenvector can be multiplied by an arbitrary non-zero scalar, and it is still an eigenvector with the same eigenvalue. This is why all signs of the eigenvector belonging to 5.58, and 0.80 are flipped in `ulab` with respect to `numpy`. This difference, however, is of absolutely no consequence. "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Computation expenses\n",
"\n",
"Since the function is based on [Givens rotations](https://en.wikipedia.org/wiki/Givens_rotation) and runs till convergence is achieved, or till the maximum number of allowed rotations is exhausted, there is no universal estimate for the time required to find the eigenvalues. However, an order of magnitude can, at least, be guessed based on the measurement below:"
]
},
{
"cell_type": "code",
"execution_count": 559,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-20T07:18:52.520515Z",
"start_time": "2019-10-20T07:18:52.499653Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"execution time: 111 us\n",
"\n"
]
}
],
"source": [
"%%micropython -pyboard 1\n",
"\n",
"import ulab as np\n",
"\n",
"@timeit\n",
"def matrix_eig(a):\n",
" return np.eig(a)\n",
"\n",
"a = np.array([[1, 2, 1, 4], [2, 5, 3, 5], [1, 3, 6, 1], [4, 5, 1, 7]], dtype=np.uint8)\n",
"\n",
"matrix_eig(a)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Polynomials"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## polyval\n",
"\n",
"numpy: https://docs.scipy.org/doc/numpy/reference/generated/numpy.polyval.html\n",
"\n",
"polyval takes two arguments, both arrays or other iterables."
]
},
{
"cell_type": "code",
"execution_count": 187,
"metadata": {
"ExecuteTime": {
"end_time": "2019-11-01T12:53:22.448303Z",
"start_time": "2019-11-01T12:53:22.435176Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"coefficients: [1, 1, 1, 0]\n",
"independent values: [0, 1, 2, 3, 4]\n",
"\n",
"values of p(x): array([0.0, 3.0, 14.0, 39.0, 84.0], dtype=float)\n",
"\n",
"ndarray (a): array([0.0, 1.0, 2.0, 3.0, 4.0], dtype=float)\n",
"value of p(a): array([0.0, 3.0, 14.0, 39.0, 84.0], dtype=float)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"p = [1, 1, 1, 0]\n",
"x = [0, 1, 2, 3, 4]\n",
"print('coefficients: ', p)\n",
"print('independent values: ', x)\n",
"print('\\nvalues of p(x): ', np.polyval(p, x))\n",
"\n",
"# the same works with one-dimensional ndarrays\n",
"a = np.array(x)\n",
"print('\\nndarray (a): ', a)\n",
"print('value of p(a): ', np.polyval(p, a))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## polyfit\n",
"\n",
"numpy: https://docs.scipy.org/doc/numpy/reference/generated/numpy.polyfit.html\n",
"\n",
"polyfit takes two, or three arguments. The last one is the degree of the polynomial that will be fitted, the last but one is an array or iterable with the `y` (dependent) values, and the first one, an array or iterable with the `x` (independent) values, can be dropped. If that is the case, `x` will be generated in the function, assuming uniform sampling. \n",
"\n",
"If the length of `x`, and `y` are not the same, the function raises a `ValueError`."
]
},
{
"cell_type": "code",
"execution_count": 189,
"metadata": {
"ExecuteTime": {
"end_time": "2019-11-01T12:54:08.326802Z",
"start_time": "2019-11-01T12:54:08.311182Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"independent values:\t array([0.0, 1.0, 2.0, 3.0, 4.0, 5.0, 6.0], dtype=float)\n",
"dependent values:\t array([9.0, 4.0, 1.0, 0.0, 1.0, 4.0, 9.0], dtype=float)\n",
"fitted values:\t\t array([1.0, -6.0, 9.000000000000004], dtype=float)\n",
"\n",
"dependent values:\t array([9.0, 4.0, 1.0, 0.0, 1.0, 4.0, 9.0], dtype=float)\n",
"fitted values:\t\t array([1.0, -6.0, 9.000000000000004], dtype=float)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"x = np.array([0, 1, 2, 3, 4, 5, 6])\n",
"y = np.array([9, 4, 1, 0, 1, 4, 9])\n",
"print('independent values:\\t', x)\n",
"print('dependent values:\\t', y)\n",
"print('fitted values:\\t\\t', np.polyfit(x, y, 2))\n",
"\n",
"# the same with missing x\n",
"print('\\ndependent values:\\t', y)\n",
"print('fitted values:\\t\\t', np.polyfit(y, 2))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Execution time\n",
"\n",
"`polyfit` is based on the inversion of a matrix (there is more on the background in https://en.wikipedia.org/wiki/Polynomial_regression), and it requires the intermediate storage of `2*N*(deg+1)` floats, where `N` is the number of entries in the input array, and `deg` is the fit's degree. The additional computation costs of the matrix inversion discussed in [inv](#inv) also apply. The example from above needs around 150 microseconds to return:"
]
},
{
"cell_type": "code",
"execution_count": 560,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-20T07:24:39.002243Z",
"start_time": "2019-10-20T07:24:38.978687Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"execution time: 153 us\n",
"\n"
]
}
],
"source": [
"%%micropython -pyboard 1\n",
"\n",
"import ulab as np\n",
"\n",
"@timeit\n",
"def time_polyfit(x, y, n):\n",
" return np.polyfit(x, y, n)\n",
"\n",
"x = np.array([0, 1, 2, 3, 4, 5, 6])\n",
"y = np.array([9, 4, 1, 0, 1, 4, 9])\n",
"\n",
"time_polyfit(x, y, 2)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Fourier transforms\n",
"\n",
"numpy: https://docs.scipy.org/doc/numpy/reference/generated/numpy.fft.ifft.html\n",
"\n",
"## fft\n",
"\n",
"Since `ulab`'s `ndarray` does not support complex numbers, the invocation of the Fourier transform differs from that in `numpy`. In `numpy`, you can simply pass an array or iterable to the function, and it will be treated as a complex array:"
]
},
{
"cell_type": "code",
"execution_count": 341,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-17T17:33:38.487729Z",
"start_time": "2019-10-17T17:33:38.473515Z"
}
},
"outputs": [
{
"data": {
"text/plain": [
"array([20.+0.j, 0.+0.j, -4.+4.j, 0.+0.j, -4.+0.j, 0.+0.j, -4.-4.j,\n",
" 0.+0.j])"
]
},
"execution_count": 341,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"fft.fft([1, 2, 3, 4, 1, 2, 3, 4])"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"**WARNING:** The array that is returned is also complex, i.e., the real and imaginary components are cast together. In `ulab`, the real and imaginary parts are treated separately: you have to pass two `ndarray`s to the function, although, the second argument is optional, in which case the imaginary part is assumed to be zero.\n",
"\n",
"**WARNING:** The function, as opposed to `numpy`, returns a 2-tuple, whose elements are two `ndarray`s, holding the real and imaginary parts of the transform separately. "
]
},
{
"cell_type": "code",
"execution_count": 114,
"metadata": {
"ExecuteTime": {
"end_time": "2020-02-16T18:38:07.294862Z",
"start_time": "2020-02-16T18:38:07.233842Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"real part:\t array([5119.996, -5.004663, -5.004798, ..., -5.005482, -5.005643, -5.006577], dtype=float)\r\n",
"\r\n",
"imaginary part:\t array([0.0, 1631.333, 815.659, ..., -543.764, -815.6588, -1631.333], dtype=float)\r\n",
"\r\n",
"real part:\t array([5119.996, -5.004663, -5.004798, ..., -5.005482, -5.005643, -5.006577], dtype=float)\r\n",
"\r\n",
"imaginary part:\t array([0.0, 1631.333, 815.659, ..., -543.764, -815.6588, -1631.333], dtype=float)\r\n",
"\n"
]
}
],
"source": [
"%%micropython -pyboard 1\n",
"\n",
"import ulab as np\n",
"from ulab import numerical\n",
"from ulab import vector\n",
"from ulab import fft\n",
"from ulab import linalg\n",
"\n",
"x = numerical.linspace(0, 10, num=1024)\n",
"y = vector.sin(x)\n",
"z = linalg.zeros(len(x))\n",
"\n",
"a, b = fft.fft(x)\n",
"print('real part:\\t', a)\n",
"print('\\nimaginary part:\\t', b)\n",
"\n",
"c, d = fft.fft(x, z)\n",
"print('\\nreal part:\\t', c)\n",
"print('\\nimaginary part:\\t', d)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## ifft\n",
"\n",
"The above-mentioned rules apply to the inverse Fourier transform. The inverse is also normalised by `N`, the number of elements, as is customary in `numpy`. With the normalisation, we can ascertain that the inverse of the transform is equal to the original array."
]
},
{
"cell_type": "code",
"execution_count": 459,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-19T13:08:17.647416Z",
"start_time": "2019-10-19T13:08:17.597456Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"original vector:\t array([0.0, 0.009775016, 0.0195491, ..., -0.5275068, -0.5357859, -0.5440139], dtype=float)\n",
"\n",
"real part of inverse:\t array([-2.980232e-08, 0.0097754, 0.0195494, ..., -0.5275064, -0.5357857, -0.5440133], dtype=float)\n",
"\n",
"imaginary part of inverse:\t array([-2.980232e-08, -1.451171e-07, 3.693752e-08, ..., 6.44871e-08, 9.34986e-08, 2.18336e-07], dtype=float)\n",
"\n"
]
}
],
"source": [
"%%micropython -pyboard 1\n",
"\n",
"import ulab as np\n",
"\n",
"x = np.linspace(0, 10, num=1024)\n",
"y = np.sin(x)\n",
"\n",
"a, b = np.fft(y)\n",
"\n",
"print('original vector:\\t', y)\n",
"\n",
"y, z = np.ifft(a, b)\n",
"# the real part should be equal to y\n",
"print('\\nreal part of inverse:\\t', y)\n",
"# the imaginary part should be equal to zero\n",
"print('\\nimaginary part of inverse:\\t', z)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Note that unlike in `numpy`, the length of the array on which the Fourier transform is carried out must be a power of 2. If this is not the case, the function raises a `ValueError` exception."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## spectrum\n",
"\n",
"In addition to the Fourier transform and its inverse, `ulab` also sports a function called `spectrum`, which returns the absolute value of the Fourier transform. This could be used to find the dominant spectral component in a time series. The arguments are treated in the same way as in `fft`, and `ifft`."
]
},
{
"cell_type": "code",
"execution_count": 460,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-19T13:08:45.027419Z",
"start_time": "2019-10-19T13:08:44.996271Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"original vector:\t array([0.0, 0.009775016, 0.0195491, ..., -0.5275068, -0.5357859, -0.5440139], dtype=float)\n",
"\n",
"spectrum:\t array([187.8641, 315.3125, 347.8804, ..., 84.4587, 347.8803, 315.3124], dtype=float)\n",
"\n"
]
}
],
"source": [
"%%micropython -pyboard 1\n",
"\n",
"import ulab as np\n",
"\n",
"x = np.linspace(0, 10, num=1024)\n",
"y = np.sin(x)\n",
"\n",
"a = np.spectrum(y)\n",
"\n",
"print('original vector:\\t', y)\n",
"print('\\nspectrum:\\t', a)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"As such, `spectrum` is really just a shorthand for `np.sqrt(a*a + b*b)`:"
]
},
{
"cell_type": "code",
"execution_count": 461,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-19T13:08:58.196537Z",
"start_time": "2019-10-19T13:08:58.155036Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"\n",
"spectrum calculated the hard way:\t array([187.8641, 315.3125, 347.8804, ..., 84.4587, 347.8803, 315.3124], dtype=float)\n",
"\n",
"spectrum calculated the lazy way:\t array([187.8641, 315.3125, 347.8804, ..., 84.4587, 347.8803, 315.3124], dtype=float)\n",
"\n"
]
}
],
"source": [
"%%micropython -pyboard 1\n",
"\n",
"import ulab as np\n",
"\n",
"x = np.linspace(0, 10, num=1024)\n",
"y = np.sin(x)\n",
"\n",
"a, b = np.fft(y)\n",
"\n",
"print('\\nspectrum calculated the hard way:\\t', np.sqrt(a*a + b*b))\n",
"\n",
"a = np.spectrum(y)\n",
"\n",
"print('\\nspectrum calculated the lazy way:\\t', a)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Computation and storage costs"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### RAM\n",
"\n",
"The FFT routine of `ulab` calculates the transform in place. This means that beyond reserving space for the two `ndarray`s that will be returned (the computation uses these two as intermediate storage space), only a handful of temporary variables, all floats or 32-bit integers, are required. "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Speed of FFTs\n",
"\n",
"A comment on the speed: a 1024-point transform implemented in python would cost around 90 ms, and 13 ms in assembly, if the code runs on the pyboard, v.1.1. You can gain a factor of four by moving to the D series \n",
"https://github.com/peterhinch/micropython-fourier/blob/master/README.md#8-performance. "
]
},
{
"cell_type": "code",
"execution_count": 494,
"metadata": {
"ExecuteTime": {
"end_time": "2019-10-19T13:25:40.540913Z",
"start_time": "2019-10-19T13:25:40.509598Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"execution time: 1985 us\n",
"\n"
]
}
],
"source": [
"%%micropython -pyboard 1\n",
"\n",
"import ulab as np\n",
"\n",
"x = np.linspace(0, 10, num=1024)\n",
"y = np.sin(x)\n",
"\n",
"np.fft(y)\n",
"\n",
"@timeit\n",
"def np_fft(y):\n",
" return np.fft(y)\n",
"\n",
"a, b = np_fft(y)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The C implementation runs in less than 2 ms on the pyboard (we have just measured that), and has been reported to run in under 0.8 ms on the D series board. That is an improvement of at least a factor of four. "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Calculating FFTs of real signals\n",
"\n",
"Now, if you have real signals, and you are really pressed for time, you can still gain a bit on speed without sacrificing anything at all. \n",
"\n",
"If you take the FFT of a real-valued signal, the real part of the transform will be symmetric, while the imaginary part will be anti-symmetric in frequency. \n",
"\n",
"If, on the other hand, the signal is imaginary-valued, then the real part of the transform will be anti-symmetric, and the imaginary part will be symmetric in frequency. These two statements follow from the definition of the Fourier transform. \n",
"\n",
"By combining the two observations above, if you place the first signal, $y_1(t)$, into the real part, and the second signal, $y_2(t)$, into the imaginary part of your input vector, i.e., $y(t) = y_1(t) + iy_2(t)$, and take the Fourier transform of the combined signal, then the Fourier transforms of the two components can be recovered as \n",
"\n",
"\\begin{eqnarray}\n",
"Y_1(k) &=& \\frac{1}{2}\\left(Y(k) + Y^*(N-k)\\right)\n",
"\\\\\n",
"Y_2(k) &=& -\\frac{i}{2}\\left(Y(k) - Y^*(N-k)\\right)\n",
"\\end{eqnarray}\n",
"where $N$ is the length of $y_1$, and $Y_1, Y_2$, and $Y$, respectively, are the Fourier transforms of $y_1, y_2$, and $y = y_1 + iy_2$."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Filter routines"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"numpy: https://docs.scipy.org/doc/numpy/reference/generated/numpy.convolve.html\n",
"\n",
"## convolve\n",
"Returns the discrete, linear convolution of two one-dimensional sequences.\n",
"\n",
"Only the ``full`` mode is supported, and the ``mode`` named parameter is not accepted. Note that all other modes can be had by slicing a ``full`` result."
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {
"ExecuteTime": {
"end_time": "2020-02-10T18:46:06.538207Z",
"start_time": "2020-02-10T18:46:06.525851Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"array([1.0, 12.0, 123.0, 1230.0, 2300.0, 3000.0], dtype=float)\n",
"\n",
"\n"
]
}
],
"source": [
"%%micropython -unix 1\n",
"\n",
"import ulab as np\n",
"\n",
"x = np.array((1,2,3))\n",
"y = np.array((1,10,100,1000))\n",
"\n",
"print(np.convolve(x, y))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Extending ulab\n",
"\n",
"As mentioned at the beginning, `ulab` is a set of sub-modules, out of which one, `extra` is explicitly reserved for user code. You should implement your functions in this sub-module, or sub-modules thereof.\n",
"\n",
"The new functions can easily be added to `extra` in a couple of simple steps. At the C level, the type definition of an `ndarray` is as follows:\n",
"\n",
"```c\n",
"typedef struct _ndarray_obj_t {\n",
" mp_obj_base_t base;\n",
" size_t m, n;\n",
" mp_obj_array_t *array;\n",
" size_t bytes;\n",
"} ndarray_obj_t;\n",
"```\n",
"\n",
"## Creating a new ndarray\n",
"\n",
"A new `ndarray` can be created by calling \n",
"\n",
"```c\n",
"ndarray_obj_t *new_ndarray = create_new_ndarray(m, n, typecode);\n",
"```\n",
"where `m`, and `n` are the number of rows and columns, respectively, and `typecode` is one of the values from \n",
"\n",
"```c\n",
"enum NDARRAY_TYPE {\n",
" NDARRAY_UINT8 = 'B',\n",
" NDARRAY_INT8 = 'b',\n",
" NDARRAY_UINT16 = 'H', \n",
" NDARRAY_INT16 = 'h',\n",
" NDARRAY_FLOAT = 'f',\n",
"};\n",
"```\n",
"or\n",
"```c\n",
"enum NDARRAY_TYPE {\n",
" NDARRAY_UINT8 = 'B',\n",
" NDARRAY_INT8 = 'b',\n",
" NDARRAY_UINT16 = 'H', \n",
" NDARRAY_INT16 = 'h',\n",
" NDARRAY_FLOAT = 'd',\n",
"};\n",
"```\n",
"The ambiguity is caused by the fact that not all platforms implement `double`, and there one has to take `float`s. But you haven't actually got to be concerned by this, because at the very beginning of `ndarray.h`, this is already taken care of: the pre-processor figures out, what the `float` implementation of the hardware platform is, and defines the `NDARRAY_FLOAT` typecode accordingly. All you have to keep in mind is that wherever you would use `float` or `double`, you have to use `mp_float_t`. That type is defined in `py/mpconfig.h` of the micropython code base.\n",
"\n",
"Therefore, a 4-by-5 matrix of type float can be created as\n",
"\n",
"```c\n",
"ndarray_obj_t *new_ndarray = create_new_ndarray(4, 5, NDARRAY_FLOAT);\n",
"```\n",
"\n",
"This function also fills in the `ndarray` structure's `m`, `n`, and `bytes` members, as well initialises the `array` member with zeros. \n",
"\n",
"Alternatively, a one-to-one copy of an `ndarray` can be gotten by calling \n",
"\n",
"```c\n",
"mp_obj_t copy_of_input_object = ndarray_copy(object_in);\n",
"```\n",
"\n",
"Note, however, that this function takes an input object of type `mp_obj_t`, and returns a copy of type `mp_obj_t`, i.e., something that can be taken from, and can immediately be returned to the interpreter. If you want to work on the data in the copy, you still have to create a pointer to it\n",
"\n",
"```c\n",
"ndarray_obj_t *ndarray = MP_OBJ_TO_PTR(copy_of_input_object);\n",
"```\n",
"\n",
"The values stored in `array` can be modified or retrieved by accessing `array->items`. Note that `array->items` is a void pointer, therefore, it must be cast before trying to access the elements. `array` has at least two useful members. One of them is `len`, which is the number of elements that the array holds, while the second one is the `typecode` that we passed to `create_new_ndarray` earlier. \n",
"\n",
"Thus, staying with our example of a 4-by-5 float matrix, we can loop through all entries as\n",
"\n",
"```c\n",
"mp_float_t *items = (mp_float_t *)new_ndarray->array->items;\n",
"mp_float_t item;\n",
"\n",
"for(size_t i=0; i < new_ndarray->array->len; i++) {\n",
" item = items[i];\n",
" // do something with item...\n",
"}\n",
"```\n",
"or, since the data are stored in the pointer in a C-like fashion, as \n",
"\n",
"```c\n",
"mp_float_t *items = (mp_float_t *)new_ndarray->array->items;\n",
"mp_float_t item;\n",
"\n",
"for(size_t m=0; m < new_ndarray->m; m++) { // loop through the rows\n",
" for(size_t n=0; n < new_ndarray->n; n++) { // loop through the columns\n",
" item = items[m*new_ndarray->n+n]; // get the (m,n) entry\n",
" // do something with item...\n",
" }\n",
"}\n",
"```\n",
"\n",
"\n",
"## Accessing data in the ndarray\n",
"\n",
"We have already seen, how the entries of an array can be accessed. If the object in question comes from the user (i.e., via the micropython interface), we can get a pointer to it by calling \n",
"\n",
"```c\n",
"ndarray_obj_t *ndarray = MP_OBJ_TO_PTR(object_in);\n",
"```\n",
"Once the pointer is at our disposal, we can get a pointer to the underlying numerical array as discussed earlier.\n",
"\n",
"If you need to find out the typecode of the array, you can get it by accessing the `typecode` member of `array`, i.e., \n",
"\n",
"```c\n",
"ndarray->array->typecode\n",
"```\n",
"should be equal to `B`, `b`, `H`, `h`, or `f`. The size of a single item is returned by the function `mp_binary_get_size('@', ndarray->array->typecode, NULL)`. This number is equal to 1, if the typecode is `B`, or `b`, 2, if the typecode is `H`, or `h`, 4, if the typecode is `f`, and 8 for `d`. \n",
"\n",
"Alternatively, the size can be found out by dividing `ndarray->bytes` with the product of `m`, and `n`, i.e., \n",
"\n",
"```c\n",
"ndarray->bytes/(ndarray->m*ndarray*n) \n",
"```\n",
"is equal to `mp_binary_get_size('@', ndarray->array->typecode, NULL)`.\n",
"\n",
"## Making certain that we have an ndarray\n",
"\n",
"A number of operations make sense for `ndarray`s only, therefore, before doing any heavy work on the data, it might be reasonable to check, whether the input argument is of the proper type. This you do by evaluating\n",
"\n",
"```c\n",
"mp_obj_is_type(object_in, &ulab_ndarray_type)\n",
"```\n",
"which should return `true`. "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Boilerplate of sorts\n",
"\n",
"To summarise the contents of the previous three sections, here is a useless function that prints out the size (`m`, and `n`) of an array, creates a copy of the input, and fills up the resulting matrix with 13.\n",
"\n",
"```c\n",
"mp_obj_t useless_function(mp_obj_t object_in) {\n",
" if(!mp_obj_is_type(object_in, &ulab_ndarray_type)) {\n",
" mp_raise_TypeError(\"useless_function takes only ndarray arguments\");\n",
" }\n",
" \n",
" mp_obj_t object_out = ndarray_copy(object_int);\n",
" \n",
" ndarray_obj_t *ndarray = MP_OBJ_TO_PTR(object_out);\n",
" printf(\"\\nsize (m): %ld, size (n): %ld\\n\", ndarray->m, ndarray->n);\n",
" printf(\"\\nlength (len): %ld, typecode: %d\\n\", ndarray->array->len, ndarray->array->typecode);\n",
" if(ndarray->array->typecode == NDARRAY_UINT8) {\n",
" // cast the pointer to the items\n",
" uint8_t *items = (uint8_t *)ndarray->array->items;\n",
" // retrieve the length of the array, and loop over the elements\n",
" for(size_t i=0; i < ndarray->array->len; i++) items[i] = 13;\n",
" } else if(ndarray->array->typecode == NDARRAY_INT8) {\n",
" int8_t *items = (int8_t *)ndarray->array->items;\n",
" for(size_t i=0; i < ndarray->array->len; i++) items[i] = 13;\n",
" } else if(ndarray->array->typecode == NDARRAY_UINT16) {\n",
" uint16_t *items = (uint16_t *)ndarray->array->items;\n",
" for(size_t i=0; i < ndarray->array->len; i++) items[i] = 13;\n",
" } else if(ndarray->array->typecode == NDARRAY_INT16) {\n",
" int16_t *items = (int16_t *)ndarray->array->items;\n",
" for(size_t i=0; i < ndarray->array->len; i++) items[i] = 13;\n",
" } else {\n",
" float *items = (mp_float_t *)ndarray->array->items;\n",
" for(size_t i=0; i < ndarray->array->len; i++) items[i] = 13;\n",
" }\n",
" return object_out;\n",
"}\n",
"```\n",
"If, on the other hand, you want to create an `ndarray` from scratch, and return that, you could work along the following lines:\n",
"\n",
"```c\n",
"mp_obj_t useless_function(mp_obj_t object_in) {\n",
" uint16_t m = mp_obj_get_int(object_in);\n",
" \n",
" ndarray_obj_t *ndarray = create_new_ndarray(1, m, NDARRAY_UINT8);\n",
" \n",
" uint8_t *items = (uint8_t *)ndarray->array->items;\n",
" // do something with the array's entries\n",
" // ...\n",
" \n",
" // and at the very end, return an mp_object_t\n",
" return MP_PTR_TO_OBJ(ndarray);\n",
"}\n",
"```\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Once the function implementation is done, you should add the function object to the globals dictionary of the `extra` sub-module as \n",
"\n",
"```c\n",
"...\n",
" MP_DEFINE_CONST_FUN_OBJ_1(useless_function_obj, userless_function);\n",
"...\n",
" STATIC const mp_map_elem_t extra_globals_table[] = {\n",
"...\n",
" { MP_OBJ_NEW_QSTR(MP_QSTR_useless), (mp_obj_t)&useless_function_obj },\n",
"...\n",
"};\n",
"```\n",
"\n",
"The exact form of the function object definition \n",
"\n",
"```c\n",
" MP_DEFINE_CONST_FUN_OBJ_1(useless_function_obj, userless_function);\n",
"```\n",
"depends naturally on what exactly you implemented, i.e., how many arguments the function takes, whether only positional arguments allowed and so on. For a thorough discussion on how function objects have to be defined, see, e.g., the [user module programming manual](#https://micropython-usermod.readthedocs.io/en/latest/).\n",
"\n",
"And with that, you are done. You can simply call the compiler as \n",
"\n",
"```bash\n",
"make BOARD=PYBV11 USER_C_MODULES=../../../ulab all\n",
"```\n",
"and the rest is taken care of."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"In the boilerplate above, we cast the pointer to `array->items` to the required type. There are certain operations, however, when you do not need the casting. If you do not want to change the array's values, only their position within the array, you can get away with copying the memory content, regardless the type. A good example for such a scenario is the transpose function in https://github.com/v923z/micropython-ulab/blob/master/code/linalg.c."
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.7.4"
},
"toc": {
"base_numbering": 1,
"nav_menu": {},
"number_sections": true,
"sideBar": true,
"skip_h1_title": false,
"title_cell": "Table of Contents",
"title_sidebar": "Contents",
"toc_cell": false,
"toc_position": {
"height": "calc(100% - 180px)",
"left": "10px",
"top": "150px",
"width": "382.797px"
},
"toc_section_display": true,
"toc_window_display": true
},
"varInspector": {
"cols": {
"lenName": 16,
"lenType": 16,
"lenVar": 40
},
"kernels_config": {
"python": {
"delete_cmd_postfix": "",
"delete_cmd_prefix": "del ",
"library": "var_list.py",
"varRefreshCmd": "print(var_dic_list())"
},
"r": {
"delete_cmd_postfix": ") ",
"delete_cmd_prefix": "rm(",
"library": "var_list.r",
"varRefreshCmd": "cat(var_dic_list()) "
}
},
"types_to_exclude": [
"module",
"function",
"builtin_function_or_method",
"instance",
"_Feature"
],
"window_display": false
}
},
"nbformat": 4,
"nbformat_minor": 4
}