Skip to content

M
mooc-rr
Commit 6509bb73 authored by Victor-M-Gomes's avatar Victor-M-Gomes
Browse files
Options

Module2 exo4

parent 33786a46
Show whitespace changes
Inline Side-by-side
Showing with 674 additions and 74 deletions
journal/journal.org
View file @ 6509bb73
# -*- mode: org -*- # -*- mode: org -*-
# -*- coding: utf-8 -*- # -*- coding: utf-8 -*-
#+STARTUP: overview indent inlineimages logdrawer #+STARTUP: overview indent inlineimages logdrawer
#+TITLE: Journal #+TITLE: Journal MOOC
#+AUTHOR: Victor Martins Gomes #+AUTHOR: Victor Martins Gomes
#+LANGUAGE: en #+LANGUAGE: en
#+TAGS: LIG(L) HOME(H) Europe(E) Blog(B) noexport(n) Stats(S) #+TAGS: LIG(L) HOME(H) Europe(E) Blog(B) noexport(n) Stats(S)
......
module2/exo1/toy_document_orgmode_python_en.org
View file @ 6509bb73
...@@ -50,6 +50,8 @@ print(matplot_lib_filename) ...@@ -50,6 +50,8 @@ print(matplot_lib_filename)
#+RESULTS: #+RESULTS:
[[file:figure_pi_mc2.png]] [[file:figure_pi_mc2.png]]
fig, ax = plt.subplots(1)
figure_pi_mc2.png]]
It is then straightforward to obtain a (not really good) approximation to $\pi$ by counting how many times, on average, $X^2 + Y^2$ is smaller than 1: It is then straightforward to obtain a (not really good) approximation to $\pi$ by counting how many times, on average, $X^2 + Y^2$ is smaller than 1:
......
module2/exo4/exercice_python_en.html 0 → 100644
View file @ 6509bb73
<?xml version="1.0" encoding="utf-8"?>
<?xml version="1.0" encoding="utf-8"?>
<!DOCTYPE html PUBLIC "-//W3C//DTD XHTML 1.0 Strict//EN"
"http://www.w3.org/TR/xhtml1/DTD/xhtml1-strict.dtd">
<html xmlns="http://www.w3.org/1999/xhtml" lang="en" xml:lang="en">
<head>
<!-- 2020-04-29 mer. 15:24 -->
<meta http-equiv="Content-Type" content="text/html;charset=utf-8" />
<meta name="viewport" content="width=device-width, initial-scale=1" />
<title>Reading and Analyzing some experimental hdf5 data</title>
<meta name="generator" content="Org mode" />
<meta name="author" content="Victor Martins Gomes" />
<style type="text/css">
<!--/*--><![CDATA[/*><!--*/
.title { text-align: center;
margin-bottom: .2em; }
.subtitle { text-align: center;
font-size: medium;
font-weight: bold;
margin-top:0; }
.todo { font-family: monospace; color: red; }
.done { font-family: monospace; color: green; }
.priority { font-family: monospace; color: orange; }
.tag { background-color: #eee; font-family: monospace;
padding: 2px; font-size: 80%; font-weight: normal; }
.timestamp { color: #bebebe; }
.timestamp-kwd { color: #5f9ea0; }
.org-right { margin-left: auto; margin-right: 0px; text-align: right; }
.org-left { margin-left: 0px; margin-right: auto; text-align: left; }
.org-center { margin-left: auto; margin-right: auto; text-align: center; }
.underline { text-decoration: underline; }
#postamble p, #preamble p { font-size: 90%; margin: .2em; }
p.verse { margin-left: 3%; }
pre {
border: 1px solid #ccc;
box-shadow: 3px 3px 3px #eee;
padding: 8pt;
font-family: monospace;
overflow: auto;
margin: 1.2em;
}
pre.src {
position: relative;
overflow: visible;
padding-top: 1.2em;
}
pre.src:before {
display: none;
position: absolute;
background-color: white;
top: -10px;
right: 10px;
padding: 3px;
border: 1px solid black;
}
pre.src:hover:before { display: inline;}
/* Languages per Org manual */
pre.src-asymptote:before { content: 'Asymptote'; }
pre.src-awk:before { content: 'Awk'; }
pre.src-C:before { content: 'C'; }
/* pre.src-C++ doesn't work in CSS */
pre.src-clojure:before { content: 'Clojure'; }
pre.src-css:before { content: 'CSS'; }
pre.src-D:before { content: 'D'; }
pre.src-ditaa:before { content: 'ditaa'; }
pre.src-dot:before { content: 'Graphviz'; }
pre.src-calc:before { content: 'Emacs Calc'; }
pre.src-emacs-lisp:before { content: 'Emacs Lisp'; }
pre.src-fortran:before { content: 'Fortran'; }
pre.src-gnuplot:before { content: 'gnuplot'; }
pre.src-haskell:before { content: 'Haskell'; }
pre.src-hledger:before { content: 'hledger'; }
pre.src-java:before { content: 'Java'; }
pre.src-js:before { content: 'Javascript'; }
pre.src-latex:before { content: 'LaTeX'; }
pre.src-ledger:before { content: 'Ledger'; }
pre.src-lisp:before { content: 'Lisp'; }
pre.src-lilypond:before { content: 'Lilypond'; }
pre.src-lua:before { content: 'Lua'; }
pre.src-matlab:before { content: 'MATLAB'; }
pre.src-mscgen:before { content: 'Mscgen'; }
pre.src-ocaml:before { content: 'Objective Caml'; }
pre.src-octave:before { content: 'Octave'; }
pre.src-org:before { content: 'Org mode'; }
pre.src-oz:before { content: 'OZ'; }
pre.src-plantuml:before { content: 'Plantuml'; }
pre.src-processing:before { content: 'Processing.js'; }
pre.src-python:before { content: 'Python'; }
pre.src-R:before { content: 'R'; }
pre.src-ruby:before { content: 'Ruby'; }
pre.src-sass:before { content: 'Sass'; }
pre.src-scheme:before { content: 'Scheme'; }
pre.src-screen:before { content: 'Gnu Screen'; }
pre.src-sed:before { content: 'Sed'; }
pre.src-sh:before { content: 'shell'; }
pre.src-sql:before { content: 'SQL'; }
pre.src-sqlite:before { content: 'SQLite'; }
/* additional languages in org.el's org-babel-load-languages alist */
pre.src-forth:before { content: 'Forth'; }
pre.src-io:before { content: 'IO'; }
pre.src-J:before { content: 'J'; }
pre.src-makefile:before { content: 'Makefile'; }
pre.src-maxima:before { content: 'Maxima'; }
pre.src-perl:before { content: 'Perl'; }
pre.src-picolisp:before { content: 'Pico Lisp'; }
pre.src-scala:before { content: 'Scala'; }
pre.src-shell:before { content: 'Shell Script'; }
pre.src-ebnf2ps:before { content: 'ebfn2ps'; }
/* additional language identifiers per "defun org-babel-execute"
in ob-*.el */
pre.src-cpp:before { content: 'C++'; }
pre.src-abc:before { content: 'ABC'; }
pre.src-coq:before { content: 'Coq'; }
pre.src-groovy:before { content: 'Groovy'; }
/* additional language identifiers from org-babel-shell-names in
ob-shell.el: ob-shell is the only babel language using a lambda to put
the execution function name together. */
pre.src-bash:before { content: 'bash'; }
pre.src-csh:before { content: 'csh'; }
pre.src-ash:before { content: 'ash'; }
pre.src-dash:before { content: 'dash'; }
pre.src-ksh:before { content: 'ksh'; }
pre.src-mksh:before { content: 'mksh'; }
pre.src-posh:before { content: 'posh'; }
/* Additional Emacs modes also supported by the LaTeX listings package */
pre.src-ada:before { content: 'Ada'; }
pre.src-asm:before { content: 'Assembler'; }
pre.src-caml:before { content: 'Caml'; }
pre.src-delphi:before { content: 'Delphi'; }
pre.src-html:before { content: 'HTML'; }
pre.src-idl:before { content: 'IDL'; }
pre.src-mercury:before { content: 'Mercury'; }
pre.src-metapost:before { content: 'MetaPost'; }
pre.src-modula-2:before { content: 'Modula-2'; }
pre.src-pascal:before { content: 'Pascal'; }
pre.src-ps:before { content: 'PostScript'; }
pre.src-prolog:before { content: 'Prolog'; }
pre.src-simula:before { content: 'Simula'; }
pre.src-tcl:before { content: 'tcl'; }
pre.src-tex:before { content: 'TeX'; }
pre.src-plain-tex:before { content: 'Plain TeX'; }
pre.src-verilog:before { content: 'Verilog'; }
pre.src-vhdl:before { content: 'VHDL'; }
pre.src-xml:before { content: 'XML'; }
pre.src-nxml:before { content: 'XML'; }
/* add a generic configuration mode; LaTeX export needs an additional
(add-to-list 'org-latex-listings-langs '(conf " ")) in .emacs */
pre.src-conf:before { content: 'Configuration File'; }
table { border-collapse:collapse; }
caption.t-above { caption-side: top; }
caption.t-bottom { caption-side: bottom; }
td, th { vertical-align:top; }
th.org-right { text-align: center; }
th.org-left { text-align: center; }
th.org-center { text-align: center; }
td.org-right { text-align: right; }
td.org-left { text-align: left; }
td.org-center { text-align: center; }
dt { font-weight: bold; }
.footpara { display: inline; }
.footdef { margin-bottom: 1em; }
.figure { padding: 1em; }
.figure p { text-align: center; }
.equation-container {
display: table;
text-align: center;
width: 100%;
}
.equation {
vertical-align: middle;
}
.equation-label {
display: table-cell;
text-align: right;
vertical-align: middle;
}
.inlinetask {
padding: 10px;
border: 2px solid gray;
margin: 10px;
background: #ffffcc;
}
#org-div-home-and-up
{ text-align: right; font-size: 70%; white-space: nowrap; }
textarea { overflow-x: auto; }
.linenr { font-size: smaller }
.code-highlighted { background-color: #ffff00; }
.org-info-js_info-navigation { border-style: none; }
#org-info-js_console-label
{ font-size: 10px; font-weight: bold; white-space: nowrap; }
.org-info-js_search-highlight
{ background-color: #ffff00; color: #000000; font-weight: bold; }
.org-svg { width: 90%; }
/*]]>*/-->
</style>
<script type="text/javascript">
// @license magnet:?xt=urn:btih:1f739d935676111cfff4b4693e3816e664797050&amp;dn=gpl-3.0.txt GPL-v3-or-Later
<!--/*--><![CDATA[/*><!--*/
function CodeHighlightOn(elem, id)
{
var target = document.getElementById(id);
if(null != target) {
elem.cacheClassElem = elem.className;
elem.cacheClassTarget = target.className;
target.className = "code-highlighted";
elem.className = "code-highlighted";
}
}
function CodeHighlightOff(elem, id)
{
var target = document.getElementById(id);
if(elem.cacheClassElem)
elem.className = elem.cacheClassElem;
if(elem.cacheClassTarget)
target.className = elem.cacheClassTarget;
}
/*]]>*///-->
// @license-end
</script>
<script type="text/x-mathjax-config">
MathJax.Hub.Config({
displayAlign: "center",
displayIndent: "0em",
"HTML-CSS": { scale: 100,
linebreaks: { automatic: "false" },
webFont: "TeX"
},
SVG: {scale: 100,
linebreaks: { automatic: "false" },
font: "TeX"},
NativeMML: {scale: 100},
TeX: { equationNumbers: {autoNumber: "AMS"},
MultLineWidth: "85%",
TagSide: "right",
TagIndent: ".8em"
}
});
</script>
<script type="text/javascript"
src="https://cdnjs.cloudflare.com/ajax/libs/mathjax/2.7.0/MathJax.js?config=TeX-AMS_HTML"></script>
</head>
<body>
<div id="content">
<h1 class="title">Reading and Analyzing some experimental hdf5 data</h1>
<div id="table-of-contents">
<h2>Table of Contents</h2>
<div id="text-table-of-contents">
<ul>
<li><a href="#orgd274c28">1. Data reading using h5py package</a></li>
</ul>
</div>
</div>
<p>
I will read and get the statistics of some data I usually acquire in the lab. I
will use this log as a future reference to when I need to read hdf5 files
contained my data.
</p>
<div id="outline-container-orgd274c28" class="outline-2">
<h2 id="orgd274c28"><span class="section-number-2">1</span> Data reading using h5py package</h2>
<div class="outline-text-2" id="text-1">
<p>
To read the data we make use of the h5py package which contains most of what
might be necessary to read hdf5 files.
</p>
<div class="org-src-container">
<pre class="src src-python"><span style="color: #51afef;">import</span> h5py <span style="color: #51afef;">as</span> h5
<span style="color: #dcaeea;">f</span> = h5.File(<span style="color: #98be65;">'testdata.h5'</span>,<span style="color: #98be65;">'r'</span>)
<span style="color: #51afef;">print</span>(f.keys())
</pre>
</div>
<pre class="example">
&lt;KeysViewHDF5 ['FileType', 'Frame', 'Waveforms']&gt;
</pre>
<p>
Now that we know what is inside, let&rsquo;s try to see which waveforms we have.
</p>
<div class="org-src-container">
<pre class="src src-python"><span style="color: #51afef;">print</span>(f[<span style="color: #98be65;">'Waveforms'</span>].keys())
</pre>
</div>
<pre class="example">
&lt;KeysViewHDF5 ['Channel 1', 'Channel 2']&gt;
</pre>
<p>
We have two channels, thus another subgroup. Gotta investigate deeper:
</p>
<div class="org-src-container">
<pre class="src src-python"><span style="color: #51afef;">print</span>(f[<span style="color: #98be65;">'Waveforms/Channel 1'</span>].keys())
<span style="color: #51afef;">print</span>(<span style="color: #98be65;">'---'</span>)
<span style="color: #51afef;">print</span>(f[<span style="color: #98be65;">'Waveforms/Channel 2'</span>].keys())
</pre>
</div>
<pre class="example">
&lt;KeysViewHDF5 ['Channel 1Data']&gt;
---
&lt;KeysViewHDF5 ['Channel 2Data']&gt;
</pre>
<p>
Now we can get the attributes of the recorded channels. They are going to be
necessary to properly map the data to a <code>numpy</code> array.
</p>
<div class="org-src-container">
<pre class="src src-python"><span style="color: #51afef;">print</span>(f[<span style="color: #98be65;">'Waveforms/Channel 1'</span>].attrs.keys())
</pre>
</div>
<pre class="example">
&lt;KeysViewHDF5 ['Count', 'DispInterpFactor', 'InterpSetting', 'MaxBandwidth', 'MinBandwidth', 'NumPoints', 'NumSegments', 'SavedInterpFactor', 'Start', 'WavAttr', 'WaveformType', 'XDispOrigin', 'XDispRange', 'XInc', 'XOrg', 'XUnits', 'YDispOrigin', 'YDispRange', 'YInc', 'YOrg', 'YUnits']&gt;
</pre>
<p>
Besides the channels we&rsquo;ve read, a time series should be created. Each of
its elements should correspond to a data sample. To do that we will use the
information stored in one of the channels (supposing they were acquired using
the same oscilloscope config).
</p>
<div class="org-src-container">
<pre class="src src-python"><span style="color: #51afef;">import</span> numpy <span style="color: #51afef;">as</span> np
<span style="color: #51afef;">import</span> matplotlib.pyplot <span style="color: #51afef;">as</span> plt
<span style="color: #dcaeea;">idata1</span> = f[<span style="color: #98be65;">'Waveforms/Channel 1'</span>].attrs.get(<span style="color: #98be65;">'XDispRange'</span>)
<span style="color: #5B6268;"># </span><span style="color: #5B6268;">XDispRange and YDispRange are the total range covered in the oscilloscope display.</span>
<span style="color: #dcaeea;">idata2</span> = f[<span style="color: #98be65;">'Waveforms/Channel 1'</span>].attrs.get(<span style="color: #98be65;">'XOrg'</span>)
<span style="color: #5B6268;"># </span><span style="color: #5B6268;">XOrg is the x-origin of the data and will be the same as XDispOrigin if the saved</span>
<span style="color: #5B6268;"># </span><span style="color: #5B6268;">data is constrained to the displayed area.</span>
<span style="color: #dcaeea;">idata3</span> = f[<span style="color: #98be65;">'Waveforms/Channel 1'</span>].attrs.get(<span style="color: #98be65;">'XInc'</span>)
<span style="color: #5B6268;"># </span><span style="color: #5B6268;">Xinc is the incremental value of X, or the sampling period (i.e. XDispRange/NumPoints).</span>
<span style="color: #dcaeea;">idata4</span> = f[<span style="color: #98be65;">'Waveforms/Channel 1'</span>].attrs.get(<span style="color: #98be65;">'NumPoints'</span>)
<span style="color: #dcaeea;">temp1</span> = idata3*(np.arange(<span style="color: #da8548; font-weight: bold;">1</span>,idata4+<span style="color: #da8548; font-weight: bold;">1</span>))+idata2
plt.figure()
plt.plot(temp1,<span style="color: #98be65;">'b'</span>)
plt.grid(linestyle=<span style="color: #98be65;">'dotted'</span>)
plt.axis(<span style="color: #98be65;">'tight'</span>)
plt.title(<span style="color: #98be65;">'Times series'</span>)
plt.savefig(filename)
<span style="color: #51afef;">print</span>(filename)
</pre>
</div>
<div class="figure">
<p><img src="time_plot.png" alt="time_plot.png" />
</p>
</div>
<p>
Once the time array has been created we can move to reading the data.
</p>
<div class="org-src-container">
<pre class="src src-python"><span style="color: #dcaeea;">yinc</span> = f[<span style="color: #98be65;">'Waveforms/Channel 1'</span>].attrs.get(<span style="color: #98be65;">'YInc'</span>) <span style="color: #5B6268;"># </span><span style="color: #5B6268;">Vertical sampling period</span>
<span style="color: #dcaeea;">yorg</span> = f[<span style="color: #98be65;">'Waveforms/Channel 1'</span>].attrs.get(<span style="color: #98be65;">'YOrg'</span>) <span style="color: #5B6268;"># </span><span style="color: #5B6268;">Closest value to origin (y=0)</span>
<span style="color: #dcaeea;">data1</span> = yorg + yinc*f[<span style="color: #98be65;">'Waveforms/Channel 1/Channel 1Data'</span>]
<span style="color: #dcaeea;">yinc</span> = f[<span style="color: #98be65;">'Waveforms/Channel 2'</span>].attrs.get(<span style="color: #98be65;">'YInc'</span>)
<span style="color: #dcaeea;">yorg</span> = f[<span style="color: #98be65;">'Waveforms/Channel 2'</span>].attrs.get(<span style="color: #98be65;">'YOrg'</span>)
<span style="color: #dcaeea;">data2</span> = yorg + yinc*f[<span style="color: #98be65;">'Waveforms/Channel 2/Channel 2Data'</span>]
<span style="color: #51afef;">print</span>(<span style="color: #98be65;">'Data type data:'</span>, data1.dtype)
</pre>
</div>
<pre class="example">
Data type data: float64
</pre>
<p>
Now we can plot the data to see it.
</p>
<div class="org-src-container">
<pre class="src src-python">plt.close(<span style="color: #98be65;">'all'</span>)
plt.figure()
plt.plot(temp1*1e6,data1,<span style="color: #98be65;">'b'</span>)
plt.plot(temp1*1e6,data2,<span style="color: #98be65;">'r'</span>)
plt.grid(linestyle=<span style="color: #98be65;">'dotted'</span>)
plt.axis(<span style="color: #98be65;">'tight'</span>)
plt.xlabel(<span style="color: #98be65;">'Time ($\mu$s)'</span>)
plt.ylabel(<span style="color: #98be65;">'Amplitude (Volts)'</span>)
plt.title(<span style="color: #98be65;">'Trigger (blue) and Data (red)'</span>)
plt.savefig(filename)
<span style="color: #51afef;">print</span>(filename)
</pre>
</div>
<div class="figure">
<p><img src="data_plot.png" alt="data_plot.png" />
</p>
</div>
<p>
Data looks nice and amplitudes agree with what was observed in the oscilloscope
display. Therefore out reading method works.
</p>
<p>
Some statistics from the data can now be computed, once we can trust the data
was correctly read.
</p>
<div class="org-src-container">
<pre class="src src-python"><span style="color: #dcaeea;">avg_data</span> = np.mean(data2)
<span style="color: #dcaeea;">med_data</span> = np.median(data2)
<span style="color: #dcaeea;">std_data</span> = np.std(data2,ddof=<span style="color: #da8548; font-weight: bold;">1</span>)
<span style="color: #dcaeea;">min_data</span> = np.amin(data2)
<span style="color: #dcaeea;">max_data</span> = np.amax(data2)
<span style="color: #5B6268;">#</span><span style="color: #5B6268;">------------------------</span>
<span style="color: #51afef;">print</span>(<span style="color: #98be65;">'Data average (Volts): '</span>,avg_data)
<span style="color: #51afef;">print</span>(<span style="color: #98be65;">'Data median (Volts): '</span>,med_data)
<span style="color: #51afef;">print</span>(<span style="color: #98be65;">'Data standard deviation (Volts): '</span>,std_data)
<span style="color: #51afef;">print</span>(<span style="color: #98be65;">'Data minimum (Volts): '</span>,min_data)
<span style="color: #51afef;">print</span>(<span style="color: #98be65;">'Data maximum (Volts): '</span>,max_data)
</pre>
</div>
<pre class="example">
Data average (Volts): 0.0019089614972975556
Data median (Volts): 0.004743201433087485
Data standard deviation (Volts): 0.039513536115368106
Data minimum (Volts): -0.1612688487249741
Data maximum (Volts): 0.15652564729188667
</pre>
<p>
The data was recorded when the input channel was turned off, thus from that we
can estimated the signature &ldquo;noise&rdquo; of the system when the input will be on.
For that we will look at the histogram and see if looks like a Normal/Gaussian
continuous probability distribution.
</p>
<div class="org-src-container">
<pre class="src src-python"><span style="color: #51afef;">from</span> scipy.stats <span style="color: #51afef;">import</span> norm
plt.figure()
plt.hist(data2,color=<span style="color: #98be65;">'b'</span>,histtype=<span style="color: #98be65;">'bar'</span>,ec=<span style="color: #98be65;">'black'</span>,density=<span style="color: #a9a1e1;">True</span>)
<span style="color: #dcaeea;">x</span> = np.linspace(min_data,max_data,idata4)
<span style="color: #dcaeea;">p</span> = norm.pdf(x,avg_data,std_data)
plt.plot(x,p,<span style="color: #98be65;">'k'</span>,linewidth=<span style="color: #da8548; font-weight: bold;">2</span>)
plt.grid(linestyle=<span style="color: #98be65;">'dotted'</span>)
plt.axis(<span style="color: #98be65;">'tight'</span>)
<span style="color: #dcaeea;">title</span> = <span style="color: #98be65;">"Fitting Gaussian to histogram: $\mu$ = %.5e, $\sigma$ = %.5e"</span> % (avg_data,std_data)
plt.title(title)
plt.savefig(filename)
<span style="color: #51afef;">print</span>(filename)
</pre>
</div>
<div class="figure">
<p><img src="data_hist.png" alt="data_hist.png" />
</p>
</div>
<p>
At last we can use the D&rsquo;Agostino&rsquo;s K\(^2\) test to infer goodness-of-fit or if
the assumption of a Normal distribution really holds.
</p>
<div class="org-src-container">
<pre class="src src-python"><span style="color: #51afef;">from</span> scipy.stats <span style="color: #51afef;">import</span> normaltest
<span style="color: #dcaeea;">alpha</span> = 1e-<span style="color: #da8548; font-weight: bold;">3</span>
<span style="color: #dcaeea;">k2</span>, <span style="color: #dcaeea;">p</span> = normaltest(data2)
<span style="color: #51afef;">print</span>(<span style="color: #98be65;">"p = {:.10e}"</span>.<span style="color: #c678dd;">format</span>(p))
</pre>
</div>
<p>
<a href="p = 0.0000000000e+00">p = 0.0000000000e+00</a>
</p>
<p>
Since the p-value is (way) smaller than the alpha value we can conclude the null
hypothesis can be rejected, that is, a Gaussian distribution does not represent
the data.
</p>
</div>
</div>
</div>
<div id="postamble" class="status">
<p class="date">Date: 2020-04-28 Tuesday</p>
<p class="author">Author: Victor Martins Gomes</p>
<p class="date">Created: 2020-04-29 mer. 15:24</p>
</div>
</body>
</html>
module2/exo4/exercice_python_en.org
View file @ 6509bb73
#+TITLE: Your title #+TITLE: Reading and Analyzing some experimental hdf5 data
#+AUTHOR: Your name #+AUTHOR: Victor Martins Gomes
#+DATE: Today's date #+DATE: 2020-04-28 Tuesday
#+LANGUAGE: en #+LANGUAGE: en
# #+PROPERTY: header-args :eval never-export #+PROPERTY: header-args :session :exports both :eval never-export
#+HTML_HEAD: <link rel="stylesheet" type="text/css" href="http://www.pirilampo.org/styles/readtheorg/css/htmlize.css"/> I will read and get the statistics of some data I usually acquire in the lab. I
#+HTML_HEAD: <link rel="stylesheet" type="text/css" href="http://www.pirilampo.org/styles/readtheorg/css/readtheorg.css"/> will use this log as a future reference to when I need to read hdf5 files
#+HTML_HEAD: <script src="https://ajax.googleapis.com/ajax/libs/jquery/2.1.3/jquery.min.js"></script> contained my data.
#+HTML_HEAD: <script src="https://maxcdn.bootstrapcdn.com/bootstrap/3.3.4/js/bootstrap.min.js"></script>
#+HTML_HEAD: <script type="text/javascript" src="http://www.pirilampo.org/styles/lib/js/jquery.stickytableheaders.js"></script>
#+HTML_HEAD: <script type="text/javascript" src="http://www.pirilampo.org/styles/readtheorg/js/readtheorg.js"></script>
* Some explanations * Data reading using h5py package
To read the data we make use of the h5py package which contains most of what
might be necessary to read hdf5 files.
This is an org-mode document with code examples in R. Once opened in #+begin_src python :results output :session *python* :exports both
Emacs, this document can easily be exported to HTML, PDF, and Office import h5py as h5
formats. For more information on org-mode, see f = h5.File('testdata.h5','r')
https://orgmode.org/guide/. print(f.keys())
#+end_src
When you type the shortcut =C-c C-e h o=, this document will be #+RESULTS:
exported as HTML. All the code in it will be re-executed, and the : <KeysViewHDF5 ['FileType', 'Frame', 'Waveforms']>
results will be retrieved and included into the exported document. If
you do not want to re-execute all code each time, you can delete the #
and the space before ~#+PROPERTY:~ in the header of this document.
Like we showed in the video, Python code is included as follows (and Now that we know what is inside, let's try to see which waveforms we have.
is exxecuted by typing ~C-c C-c~):
#+begin_src python :results output :exports both #+begin_src python :results output :session *python* :exports both
print("Hello world!") print(f['Waveforms'].keys())
#+end_src #+end_src
#+RESULTS: #+RESULTS:
: Hello world! : <KeysViewHDF5 ['Channel 1', 'Channel 2']>
And now the same but in an Python session. With a session, Python's We have two channels, thus another subgroup. Gotta investigate deeper:
state, i.e. the values of all the variables, remains persistent from #+begin_src python :results output :session *python* :exports both
one code block to the next. The code is still executed using ~C-c print(f['Waveforms/Channel 1'].keys())
C-c~. print('---')
print(f['Waveforms/Channel 2'].keys())
#+end_src
#+begin_src python :results output :session :exports both #+RESULTS:
import numpy : <KeysViewHDF5 ['Channel 1Data']>
x=numpy.linspace(-15,15) : ---
print(x) : <KeysViewHDF5 ['Channel 2Data']>
Now we can get the attributes of the recorded channels. They are going to be
necessary to properly map the data to a =numpy= array.
#+begin_src python :results output :session *python* :exports both
print(f['Waveforms/Channel 1'].attrs.keys())
#+end_src #+end_src
#+RESULTS: #+RESULTS:
#+begin_example : <KeysViewHDF5 ['Count', 'DispInterpFactor', 'InterpSetting', 'MaxBandwidth', 'MinBandwidth', 'NumPoints', 'NumSegments', 'SavedInterpFactor', 'Start', 'WavAttr', 'WaveformType', 'XDispOrigin', 'XDispRange', 'XInc', 'XOrg', 'XUnits', 'YDispOrigin', 'YDispRange', 'YInc', 'YOrg', 'YUnits']>
[-15. -14.3877551 -13.7755102 -13.16326531 -12.55102041
-11.93877551 -11.32653061 -10.71428571 -10.10204082 -9.48979592 Besides the channels we've read, a time series should be created. Each of
-8.87755102 -8.26530612 -7.65306122 -7.04081633 -6.42857143 its elements should correspond to a data sample. To do that we will use the
-5.81632653 -5.20408163 -4.59183673 -3.97959184 -3.36734694 information stored in one of the channels (supposing they were acquired using
-2.75510204 -2.14285714 -1.53061224 -0.91836735 -0.30612245 the same oscilloscope config).
0.30612245 0.91836735 1.53061224 2.14285714 2.75510204 #+begin_src python :results output file :var filename="time_plot.png" :exports both :session *python*
3.36734694 3.97959184 4.59183673 5.20408163 5.81632653 import numpy as np
6.42857143 7.04081633 7.65306122 8.26530612 8.87755102
9.48979592 10.10204082 10.71428571 11.32653061 11.93877551
12.55102041 13.16326531 13.7755102 14.3877551 15. ]
#+end_example
Finally, an example for graphical output:
#+begin_src python :results output file :session :var matplot_lib_filename="./cosxsx.png" :exports results
import matplotlib.pyplot as plt import matplotlib.pyplot as plt
plt.figure(figsize=(10,5)) idata1 = f['Waveforms/Channel 1'].attrs.get('XDispRange')
plt.plot(x,numpy.cos(x)/x) # XDispRange and YDispRange are the total range covered in the oscilloscope display.
plt.tight_layout() idata2 = f['Waveforms/Channel 1'].attrs.get('XOrg')
# XOrg is the x-origin of the data and will be the same as XDispOrigin if the saved
# data is constrained to the displayed area.
idata3 = f['Waveforms/Channel 1'].attrs.get('XInc')
# Xinc is the incremental value of X, or the sampling period (i.e. XDispRange/NumPoints).
idata4 = f['Waveforms/Channel 1'].attrs.get('NumPoints')
temp1 = idata3*(np.arange(1,idata4+1))+idata2
plt.figure()
plt.plot(temp1,'b')
plt.grid(linestyle='dotted')
plt.axis('tight')
plt.title('Times series')
plt.savefig(filename)
print(filename)
#+end_src
#+RESULTS:
[[file:time_plot.png]]
Once the time array has been created we can move to reading the data.
#+begin_src python :results output :session *python* :exports both
yinc = f['Waveforms/Channel 1'].attrs.get('YInc') # Vertical sampling period
yorg = f['Waveforms/Channel 1'].attrs.get('YOrg') # Closest value to origin (y=0)
data1 = yorg + yinc*f['Waveforms/Channel 1/Channel 1Data']
yinc = f['Waveforms/Channel 2'].attrs.get('YInc')
yorg = f['Waveforms/Channel 2'].attrs.get('YOrg')
data2 = yorg + yinc*f['Waveforms/Channel 2/Channel 2Data']
print('Data type data:', data1.dtype)
#+end_src
#+RESULTS:
: Data type data: float64
Now we can plot the data to see it.
#+begin_src python :results output file :var filename="data_plot.png" :exports both :session *python*
plt.close('all')
plt.figure()
plt.plot(temp1*1e6,data1,'b')
plt.plot(temp1*1e6,data2,'r')
plt.grid(linestyle='dotted')
plt.axis('tight')
plt.xlabel('Time ($\mu$s)')
plt.ylabel('Amplitude (Volts)')
plt.title('Trigger (blue) and Data (red)')
plt.savefig(filename)
print(filename)
#+end_src
#+RESULTS:
[[file:data_plot.png]]
Data looks nice and amplitudes agree with what was observed in the oscilloscope
display. Therefore out reading method works.
Some statistics from the data can now be computed, once we can trust the data
was correctly read.
#+begin_src python :results output :session *python* :exports both
avg_data = np.mean(data2)
med_data = np.median(data2)
std_data = np.std(data2,ddof=1)
min_data = np.amin(data2)
max_data = np.amax(data2)
#------------------------
print('Data average (Volts): ',avg_data)
print('Data median (Volts): ',med_data)
print('Data standard deviation (Volts): ',std_data)
print('Data minimum (Volts): ',min_data)
print('Data maximum (Volts): ',max_data)
#+end_src
#+RESULTS:
: Data average (Volts): 0.0019089614972975556
: Data median (Volts): 0.004743201433087485
: Data standard deviation (Volts): 0.039513536115368106
: Data minimum (Volts): -0.1612688487249741
: Data maximum (Volts): 0.15652564729188667
The data was recorded when the input channel was turned off, thus from that we
can estimated the signature "noise" of the system when the input will be on.
For that we will look at the histogram and see if looks like a Normal/Gaussian
continuous probability distribution.
#+begin_src python :results output file :var filename="data_hist.png" :exports both :session *python*
from scipy.stats import norm
plt.figure()
plt.hist(data2,color='b',histtype='bar',ec='black',density=True)
x = np.linspace(min_data,max_data,idata4)
p = norm.pdf(x,avg_data,std_data)
plt.plot(x,p,'k',linewidth=2)
plt.savefig(matplot_lib_filename) plt.grid(linestyle='dotted')
print(matplot_lib_filename) plt.axis('tight')
title = "Fitting Gaussian to histogram: $\mu$ = %.5e, $\sigma$ = %.5e" % (avg_data,std_data)
plt.title(title)
plt.savefig(filename)
print(filename)
#+end_src
#+RESULTS:
[[file:data_hist.png]]
At last we can use the D'Agostino's K$^2$ test to infer goodness-of-fit or if
the assumption of a Normal distribution really holds.
#+begin_src python :results output file :var filename="data_plot.png" :exports both :session *python*
from scipy.stats import normaltest
alpha = 1e-3
k2, p = normaltest(data2)
print("p = {:.10e}".format(p))
#+end_src #+end_src
#+RESULTS: #+RESULTS:
[[file:./cosxsx.png]] [[file:p = 0.0000000000e+00]]
Note the parameter ~:exports results~, which indicates that the code Since the p-value is (way) smaller than the alpha value we can conclude the null
will not appear in the exported document. We recommend that in the hypothesis can be rejected, that is, a Gaussian distribution does not represent
context of this MOOC, you always leave this parameter setting as the data.
~:exports both~, because we want your analyses to be perfectly
transparent and reproducible.
Watch out: the figure generated by the code block is /not/ stored in
the org document. It's a plain file, here named ~cosxsx.png~. You have
to commit it explicitly if you want your analysis to be legible and
understandable on GitLab.
Finally, don't forget that we provide in the resource section of this
MOOC a configuration with a few keyboard shortcuts that allow you to
quickly create code blocks in Python by typing ~<p~, ~<P~ or ~<PP~
followed by ~Tab~.
Now it's your turn! You can delete all this information and replace it
by your computational document.
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment