(************** Content-type: application/mathematica **************
                     CreatedBy='Mathematica 5.2'

                    Mathematica-Compatible Notebook

This notebook can be used with any Mathematica-compatible
application, such as Mathematica, MathReader or Publicon. The data
for the notebook starts with the line containing stars above.

To get the notebook into a Mathematica-compatible application, do
one of the following:

* Save the data starting with the line of stars above into a file
  with a name ending in .nb, then open the file inside the
  application;

* Copy the data starting with the line of stars above to the
  clipboard, then use the Paste menu command inside the application.

Data for notebooks contains only printable 7-bit ASCII and can be
sent directly in email or through ftp in text mode.  Newlines can be
CR, LF or CRLF (Unix, Macintosh or MS-DOS style).

NOTE: If you modify the data for this notebook not in a Mathematica-
compatible application, you must delete the line below containing
the word CacheID, otherwise Mathematica-compatible applications may
try to use invalid cache data.

For more information on notebooks and Mathematica-compatible 
applications, contact Wolfram Research:
  web: http://www.wolfram.com
  email: info@wolfram.com
  phone: +1-217-398-0700 (U.S.)

Notebook reader applications are available free of charge from 
Wolfram Research.
*******************************************************************)

(*CacheID: 232*)


(*NotebookFileLineBreakTest
NotebookFileLineBreakTest*)
(*NotebookOptionsPosition[     19223,        624]*)
(*NotebookOutlinePosition[     20601,        670]*)
(*  CellTagsIndexPosition[     20528,        664]*)
(*WindowFrame->Normal*)



Notebook[{

Cell[CellGroupData[{
Cell["Denoising", "Title"],

Cell["\<\
Wavelet Workshop
June 6-9, 2007
University of St. Thomas\
\>", "Subtitle"],

Cell[CellGroupData[{

Cell["Objectives", "Section"],

Cell["\<\
The purpose of this notebook is to introduce you to some basic wavelet-based \
denoising methods.
\
\>", "Text"]
}, Open  ]],

Cell[CellGroupData[{

Cell["DiscreteWavelets", "Section"],

Cell[TextData[{
  StyleBox["You should run this cell each time you open this notebook!!",
    FontColor->RGBColor[1, 0, 0]],
  "  It loads the ",
  StyleBox["Mathematica",
    FontSlant->"Italic"],
  " package ",
  StyleBox["DiscreteWavelets",
    FontFamily->"Courier"],
  " for use in subsequent computations."
}], "Text"],

Cell[BoxData[
    \(<< DiscreteWavelets`DiscreteWavelets`\)], "Input",
  FontColor->GrayLevel[0.500008]]
}, Open  ]],

Cell[CellGroupData[{

Cell["Help on DiscreteWavelets", "Section"],

Cell[TextData[{
  "If you ever need help with ",
  StyleBox["DiscreteWavelets",
    FontFamily->"Courier"],
  ", go to ",
  StyleBox["Help",
    FontSlant->"Italic"],
  ", then ",
  StyleBox["Help Browser",
    FontSlant->"Italic"],
  ", and click on ",
  StyleBox["AddOns & Links",
    FontSlant->"Italic"],
  ".  If you scroll down you will find ",
  StyleBox["DiscreteWavelets",
    FontFamily->"Courier"],
  ".  "
}], "Text"]
}, Open  ]],

Cell[CellGroupData[{

Cell["The Shrinkage Function", "Section"],

Cell[TextData[{
  "The shrinkage function is defined in DiscreteWavelets and called ",
  StyleBox["ShrinkageFunction",
    FontFamily->"Courier"],
  ".  Here is the shrinkage function plotted for \[Lambda]=2."
}], "Text"],

Cell[BoxData[
    RowBox[{
      RowBox[{"Plot", "[", 
        RowBox[{
          RowBox[{
            StyleBox["ShrinkageFunction",
              FontColor->GrayLevel[0.500008]], "[", \(v, 2\), "]"}], 
          ",", \({v, \(-5\), 5}\), 
          ",", \(Ticks \[Rule] {{\(-4\), \(-3\), \(-2\), \(-1\), 1, 2, 3, 
                4}, {\(-3\), \(-2\), \(-1\), 1, 2, 3}}\), 
          ",", \(PlotStyle \[Rule] Blue\)}], "]"}], ";"}]], "Input"],

Cell["Let's try it on an example", "Text"],

Cell[BoxData[{\(\(v = {\(-3.2\), \(-2.5\), 1, 1, 1.1, 5, 3.2, 
        1.9, \(-4\), \(-6\)}\)\(\[IndentingNewLine]\)
    \), "\[IndentingNewLine]", \(lambda = 2;\), "\[IndentingNewLine]", 
    RowBox[{"Map", "[", 
      RowBox[{
        RowBox[{
          RowBox[{
            StyleBox["ShrinkageFunction",
              FontColor->GrayLevel[0.500008]], "[", \(#, lambda\), "]"}], 
          "&"}], ",", "v"}], "]"}]}], "Input"]
}, Open  ]],

Cell[CellGroupData[{

Cell["A Sample Signal ", "Section"],

Cell["\<\
For our example, we will use a sample signal suggested by Donoho.  He calls \
this function the Heavisine function.
    \
\>", "Text"],

Cell[BoxData[
    \(\(heavisine[t_] := 
        4*Sin[4*Pi*t] - Sign[t -  .3] - Sign[ .72 - t];\)\)], "Input"],

Cell[TextData[{
  "Let's sample this function to create our true signal ",
  StyleBox["v",
    FontWeight->"Bold"],
  " and then ",
  StyleBox["ListPlot",
    FontWeight->"Bold"],
  " it."
}], "Text"],

Cell[BoxData[{
    \(\(n = 2048;\)\), "\n", 
    \(\(v = 
        Table[heavisine[k/n], {k, 0, n - 1}];\)\), "\[IndentingNewLine]", 
    \(\(ListPlot[v, PlotStyle \[Rule] Maroon];\)\)}], "Input"]
}, Open  ]],

Cell[CellGroupData[{

Cell["Creating Noise", "Section"],

Cell["These packages will be needed to create white noise.", "Text"],

Cell[BoxData[{
    \(<< Statistics`ContinuousDistributions`\), "\n", 
    \(<< Statistics`DescriptiveStatistics`\)}], "Input"],

Cell["Now we need some Gaussian white noise.", "Text"],

Cell[BoxData[{
    \(\(nd = NormalDistribution[0. , 1. ];\)\), "\n", 
    \(\(SeedRandom[];\)\), "\n", 
    \(\(noise = 
        Table[Random[nd], {k, 1, 
            n}];\)\[IndentingNewLine]\), "\[IndentingNewLine]", 
    \(\(ListPlot[noise];\)\)}], "Input"],

Cell["If you want, you can actually \"play\" the noise:", "Text"],

Cell[BoxData[
    \(\(ListPlay[noise, SampleRate \[Rule] 11025];\)\)], "Input"],

Cell["\<\
We create the noisy vector by picking a noise level and simply adding v to a \
scaled version of the noise.\
\>", "Text"],

Cell[BoxData[{
    \(\(sigma =  .5;\)\), "\[IndentingNewLine]", 
    \(\(y = v + sigma*noise;\)\[IndentingNewLine]\), "\[IndentingNewLine]", 
    \(\(plt = 
        ListPlot[y, 
          Ticks \[Rule] {{0,  .2,  .4,  .6,  .8, 1}, {\(-6\), \(-4\), \(-2\), 
                0, 2, 4, 6}}, PlotStyle \[Rule] Brown];\)\)}], "Input"]
}, Open  ]],

Cell[CellGroupData[{

Cell["Estimating the Noise Level", "Section"],

Cell[TextData[{
  "To estimate the noise level, we use the built-in ",
  StyleBox["MedianDeviation",
    FontWeight->"Bold"],
  " function in the statistics package.  We have to first compute the wavelet \
transform.  We'll use the Daubechies 4 filter for this example and compute 2 \
iterations."
}], "Text"],

Cell[BoxData[{
    RowBox[{
      RowBox[{"h", "=", 
        RowBox[{"N", "[", 
          RowBox[{
            StyleBox["Daub",
              FontColor->GrayLevel[0.500008]], "[", "4", "]"}], "]"}]}], 
      ";"}], "\[IndentingNewLine]", 
    RowBox[{\(Print["\<We will use the Daubechies 4-term filter: \>", 
          h, "\<.\>"];\), 
      "\[IndentingNewLine]"}], "\[IndentingNewLine]", \(i = 
        5;\), "\[IndentingNewLine]", 
    RowBox[{
      RowBox[{"wt", "=", 
        RowBox[{
          StyleBox["WT1D",
            FontColor->GrayLevel[0.500008]], "[", \(y, h, 
          NumIterations \[Rule] i\), "]"}]}], ";"}], "\[IndentingNewLine]", 
    RowBox[{
      RowBox[{
        StyleBox["WaveletVectorPlot",
          FontColor->GrayLevel[0.500008]], "[", \(wt, 
        NumIterations \[Rule] i, Frame \[Rule] True, 
        DivideLines \[Rule] True, UseColors \[Rule] True, 
        PlotRange \[Rule] All\), "]"}], ";"}]}], "Input"],

Cell[TextData[{
  "We need to extract the first highpass portion.  There is a command in ",
  StyleBox["DiscreteWavelets",
    FontFamily->"Courier"],
  " to do this:"
}], "Text"],

Cell[BoxData[{
    RowBox[{
      RowBox[{
        RowBox[{"hp1", "=", 
          RowBox[{
            StyleBox["ExtractVectorPart",
              FontColor->GrayLevel[0.500008]], "[", \(wt, 
            NumIterations \[Rule] i, Iteration \[Rule] 1, 
            Region \[Rule] HighPass\), "]"}]}], ";"}], 
      "\[IndentingNewLine]"}], "\[IndentingNewLine]", \(ListPlot[hp1, 
        PlotStyle \[Rule] Olive];\)}], "Input"],

Cell["Now we estimate the noise:", "Text"],

Cell[BoxData[{
    RowBox[{
      RowBox[{
        RowBox[{"sigmahat", "=", 
          RowBox[{
            RowBox[{
              StyleBox["MAD",
                FontColor->GrayLevel[0.500008]], "[", "hp1", "]"}], "/", 
            ".6745"}]}], ";"}], 
      "\[IndentingNewLine]"}], "\[IndentingNewLine]", \(Print["\<Our estimate \
of the noise is \>", sigmahat, "\<.\>"];\)}], "Input"]
}, Open  ]],

Cell[CellGroupData[{

Cell["VISUShrink", "Section"],

Cell[TextData[{
  "To perform VISUShrink, we compute the universal tolerance and use it in \
the shrink functions on the highpass portions.  \n\nWe use the function ",
  StyleBox["WaveletVectorToList",
    FontWeight->"Bold"],
  " to easily extract all the highpass portions."
}], "Text"],

Cell[BoxData[{
    RowBox[{
      RowBox[{"wtlist", "=", 
        RowBox[{
          StyleBox["WaveletVectorToList",
            FontColor->GrayLevel[0.500008]], "[", \(wt, 
          NumIterations \[Rule] i\), "]"}]}], 
      ";"}], "\[IndentingNewLine]", \(lowpass = 
        First[wtlist];\), "\[IndentingNewLine]", 
    RowBox[{\(highpass = Flatten[Drop[wtlist, 1]];\), 
      "\[IndentingNewLine]"}], "\[IndentingNewLine]", \(lambda = 
        sigmahat*Sqrt[2*Log[Length[highpass]]];\), "\[IndentingNewLine]", 
    RowBox[{\(Print["\<The value of the universal threshold is \>", 
          lambda, "\<.\>"];\), "\[IndentingNewLine]", 
      "\[IndentingNewLine]"}], "\[IndentingNewLine]", 
    RowBox[{
      RowBox[{
        RowBox[{"newhighpass", "=", 
          RowBox[{"Map", "[", 
            RowBox[{
              RowBox[{
                RowBox[{
                  StyleBox["ShrinkageFunction",
                    FontColor->GrayLevel[0.500008]], "[", \(#, lambda\), 
                  "]"}], "&"}], ",", "highpass"}], "]"}]}], ";"}], 
      "\[IndentingNewLine]"}], "\[IndentingNewLine]", \(newwt = 
        Join[lowpass, newhighpass];\), "\[IndentingNewLine]", 
    RowBox[{
      RowBox[{
        StyleBox["WaveletVectorPlot",
          FontColor->GrayLevel[0.500008]], "[", \(newwt, 
        NumIterations \[Rule] i, DivideLines \[Rule] True, 
        Frame \[Rule] True, UseColors \[Rule] True, PlotRange \[Rule] All\), 
        "]"}], ";"}]}], "Input"],

Cell["\<\
We apply the inverse wavelet transform to obtain our estimate of v:\
\>", "Text"],

Cell[BoxData[{
    RowBox[{
      RowBox[{"vhat", "=", 
        RowBox[{
          StyleBox["IWT1D",
            FontColor->GrayLevel[0.500008]], "[", \(newwt, h, 
          NumIterations \[Rule] i\), "]"}]}], 
      ";"}], "\[IndentingNewLine]", \(ListPlot[vhat, 
        PlotStyle \[Rule] Brown];\)}], "Input"]
}, Open  ]],

Cell[CellGroupData[{

Cell["SUREShrink", "Section"],

Cell[TextData[{
  "Let's use an audio clip for testing the SUREShrink method.  In the cell \
below are two commands that can be used to access the audio files included in \
the ",
  StyleBox["DiscreteWavelets",
    FontFamily->"Courier"],
  " package.  They are similar to ",
  StyleBox["ImageList",
    FontWeight->"Bold"],
  " and ",
  StyleBox["ImageNames",
    FontWeight->"Bold"],
  "."
}], "Text"],

Cell[BoxData[{
    RowBox[{
      StyleBox["AudioList",
        FontColor->GrayLevel[0.500008]], "[", "]"}], "\[IndentingNewLine]", 
    RowBox[{
      RowBox[{"d", "=", 
        RowBox[{
          StyleBox["AudioNames",
            FontColor->GrayLevel[0.500008]], "[", "]"}]}], ";"}]}], "Input",
  CellTags->"b:0.2.1"],

Cell["We will need the Audio package to continue.", "Text"],

Cell[BoxData[
    \(<< Miscellaneous`Audio`\)], "Input"],

Cell[BoxData[{\(file = d[\([2]\)];\), "\[IndentingNewLine]", \(data\  = \ 
        ReadSoundFile[file, 
          PrintHeader \[Rule] True];\), "\[IndentingNewLine]", 
    RowBox[{
      RowBox[{"data", "=", 
        RowBox[{
          StyleBox["ChopVector",
            FontColor->GrayLevel[0.500008]], "[", \(data, 
          PowersOfTwo \[Rule] 4\), "]"}]}], 
      ";"}], "\[IndentingNewLine]", \(Print["\<The length of the audio clip \
is \>", Length[data], "\<.\>"];\)}], "Input",
  CellTags->"b:0.2.1"],

Cell["The cell below plays the clip.", "Text",
  CellTags->"b:0.2.1"],

Cell[BoxData[{
    \(\(srate = 22050;\)\), "\[IndentingNewLine]", 
    \(\(ListPlay[data, SampleRate \[Rule] srate];\)\)}], "Input",
  CellTags->"b:0.2.1"],

Cell["Now let's create some noise and add it to the audio clip.", "Text"],

Cell[BoxData[{
    \(<< Statistics`ContinuousDistributions`\), "\n", 
    \(<< Statistics`DescriptiveStatistics`\), "\[IndentingNewLine]", 
    \(\(n = Length[data];\)\), "\[IndentingNewLine]", 
    \(\(nd = NormalDistribution[0. , 1. ];\)\), "\[IndentingNewLine]", 
    \(\(SeedRandom[];\)\), "\n", 
    \(\(noise = Table[Random[nd], {k, 1, n}];\)\), "\[IndentingNewLine]", 
    \(\(sigma = 1000. ;\)\), "\[IndentingNewLine]", 
    \(\(v = data + sigma*noise;\)\), "\[IndentingNewLine]", 
    \(\(ListPlay[v, SampleRate \[Rule] srate];\)\)}], "Input",
  CellTags->"b:0.2.1"],

Cell["\<\
We next compute four iterations of the wavelet transform.  We'll use the \
Coiflet C6 filter.\
\>", "Text",
  CellTags->"b:0.2.1"],

Cell[BoxData[{\(i = 4;\), "\[IndentingNewLine]", 
    RowBox[{
      RowBox[{"wt", "=", 
        RowBox[{
          StyleBox["WT1D",
            FontColor->GrayLevel[0.500008]], "[", \(v, Coif[1], 
          NumIterations \[Rule] i\), "]"}]}], ";"}], "\[IndentingNewLine]", 
    RowBox[{
      RowBox[{
        StyleBox["WaveletVectorPlot",
          FontColor->GrayLevel[0.500008]], "[", \(wt, 
        NumIterations \[Rule] i, PlotRange \[Rule] All, 
        PlotLabel -> "\<Wavelet Transform\>", UseColors \[Rule] True\), "]"}],
       ";"}]}], "Input",
  CellTags->"b:0.2.1"],

Cell["\<\
Our first task is to estimate the noise level.Of course,we know the noise \
level but it is good to see how the estimator works.\
\>", "Text"],

Cell[BoxData[{
    RowBox[{
      RowBox[{"wtlist", "=", 
        RowBox[{
          StyleBox["WaveletVectorToList",
            FontColor->GrayLevel[0.500008]], "[", \(wt, 
          NumIterations \[Rule] i\), "]"}]}], 
      ";"}], "\n", \(hp1 = Last[wtlist];\), "\n", 
    RowBox[{
      RowBox[{"sigmahat", "=", 
        RowBox[{
          RowBox[{
            StyleBox["MAD",
              FontColor->GrayLevel[0.500008]], "[", "hp1", "]"}], "/", 
          ".6745"}]}], 
      ";"}], "\n", \(Print["\<The estimate of the noise level is \>", 
        sigmahat, "\<.\>"];\)}], "Input",
  CellTags->"b:0.2.1"],

Cell[TextData[{
  "To apply SUREShrink, we want to compute a different \[Lambda] for each \
highpass portion.  Let's first grab the highpass portions using ",
  StyleBox["WaveletVectorToList",
    FontWeight->"Bold"],
  ".  Now highpass is a list of ",
  StyleBox["i",
    FontWeight->"Bold"],
  " elements.  Note that we divide the highpass elements by the noise \
estimate to obtain a noise level of 1."
}], "Text"],

Cell[BoxData[{
    \(\(lowpass = First[wtlist];\)\), "\[IndentingNewLine]", 
    \(\(highpass = Drop[wtlist, 1]/sigmahat;\)\)}], "Input"],

Cell[TextData[{
  "Remember, for Stein's theorem to work, the variance must be one.    So we \
need to divide all the highpass portions by ",
  StyleBox["sigmahat",
    FontWeight->"Bold"],
  " when computing the tolerances.  We will do all the computations with a ",
  StyleBox["Map",
    FontWeight->"Bold"],
  " command:"
}], "Text"],

Cell[BoxData[{
    RowBox[{
      RowBox[{"lambdaS", "=", 
        RowBox[{"Map", "[", 
          RowBox[{
            StyleBox["DonohoSure",
              FontColor->GrayLevel[0.500008]], ",", "highpass"}], "]"}]}], 
      ";"}], "\n", \(Table[
        Print["\<The SURE threshold for highpass portion \>", 
          i - k + 1, "\< is \>", lambdaS[\([\)\(k\)\(]\)], "\<.\n\>"], {k, 1, 
          i}];\)}], "Input"]
}, Open  ]],

Cell[CellGroupData[{

Cell["Apply the SUREShrink Threshold", "Section"],

Cell["\<\
We will again use a Table to apply the tolerances to each of the highpass \
portions:\
\>", "Text"],

Cell[BoxData[
    RowBox[{
      RowBox[{"newhighpass", "=", 
        RowBox[{"Table", "[", 
          RowBox[{
            RowBox[{"Map", "[", 
              RowBox[{
                RowBox[{
                  RowBox[{
                    StyleBox["ShrinkageFunction",
                      FontColor->GrayLevel[0.500008]], "[", \(#, 
                    lambdaS[\([k]\)]\), "]"}], "&"}], 
                ",", \(highpass[\([k]\)]\)}], "]"}], ",", \({k, 1, i}\)}], 
          "]"}]}], ";"}]], "Input"],

Cell["Now we put the wavelet transform back together and plot it:", "Text"],

Cell[BoxData[{
    RowBox[{\(newwt = Join[lowpass, sigmahat*Flatten[newhighpass]];\), 
      "\[IndentingNewLine]"}], "\n", 
    RowBox[{
      RowBox[{
        RowBox[{
          StyleBox["WaveletVectorPlot",
            FontColor->GrayLevel[0.500008]], "[", \(wt, Frame \[Rule] True, 
          PlotRange \[Rule] All, DivideLines \[Rule] True, 
          UseColors \[Rule] True, NumIterations \[Rule] i\), "]"}], ";"}], 
      "\[IndentingNewLine]"}], "\n", 
    RowBox[{
      RowBox[{
        StyleBox["WaveletVectorPlot",
          FontColor->GrayLevel[0.500008]], "[", \(newwt, Frame \[Rule] True, 
        PlotRange \[Rule] All, DivideLines \[Rule] True, 
        UseColors \[Rule] True, NumIterations \[Rule] i\), "]"}], 
      ";"}]}], "Input"],

Cell["Finally, we invert to obtained the denoised signal.", "Text"],

Cell[BoxData[{
    RowBox[{
      RowBox[{"vhat", "=", 
        RowBox[{
          StyleBox["IWT1D",
            FontColor->GrayLevel[0.500008]], "[", \(newwt, Coif[1], 
          NumIterations \[Rule] i\), "]"}]}], 
      ";"}], "\[IndentingNewLine]", \(ListPlay[vhat, 
        SampleRate \[Rule] srate];\)}], "Input"],

Cell["Here is the original signal.", "Text"],

Cell[BoxData[
    \(\(ListPlay[data, SampleRate \[Rule] srate];\)\)], "Input"],

Cell["Here is the noisy signal.", "Text"],

Cell[BoxData[
    \(\(ListPlay[v, SampleRate \[Rule] srate];\)\)], "Input"]
}, Open  ]],

Cell[CellGroupData[{

Cell["Extra Credit Problem", "Section"],

Cell["Read in the Apollo 11 clip:", "Text"],

Cell[BoxData[{
    RowBox[{
      StyleBox["AudioList",
        FontColor->GrayLevel[0.500008]], "[", "]"}], "\[IndentingNewLine]", 
    RowBox[{
      RowBox[{"d", "=", 
        RowBox[{
          StyleBox["AudioNames",
            FontColor->GrayLevel[0.500008]], "[", "]"}]}], ";"}]}], "Input"],

Cell[BoxData[{
    \(\(data = 
        ReadSoundFile[d[\([1]\)], PrintHeader \[Rule] True];\)\), "\n", 
    \(\(srate = 11025;\)\), "\[IndentingNewLine]", 
    \(\(ListPlay[data, SampleRate \[Rule] srate];\)\)}], "Input"],

Cell["\<\
This is a noisy clip.  Use either the VISUShrink method or the SUREShrink \
method to denoise the clip.\
\>", "Text"],

Cell[BoxData[
    \( (*\ Put\ your\ code\ here\ *) \)], "Input"]
}, Open  ]],

Cell[CellGroupData[{

Cell["Extra Extra Credit Problem", "Section"],

Cell["\<\
Repeat the above exercise with the Laurel and Hardy sound clip.\
\>", "Text"]
}, Open  ]]
}, Open  ]]
},
FrontEndVersion->"5.2 for Microsoft Windows",
ScreenRectangle->{{0, 1024}, {0, 685}},
AutoGeneratedPackage->None,
ScreenStyleEnvironment->"Presentation",
WindowSize->{1016, 651},
WindowMargins->{{0, Automatic}, {Automatic, 0}},
ShowSelection->True,
StyleDefinitions -> "Report.nb"
]

(*******************************************************************
Cached data follows.  If you edit this Notebook file directly, not
using Mathematica, you must remove the line containing CacheID at
the top of  the file.  The cache data will then be recreated when
you save this file from within Mathematica.
*******************************************************************)

(*CellTagsOutline
CellTagsIndex->{
  "b:0.2.1"->{
    Cell[11351, 378, 320, 9, 77, "Input",
      CellTags->"b:0.2.1"],
    Cell[11795, 394, 509, 11, 129, "Input",
      CellTags->"b:0.2.1"],
    Cell[12307, 407, 69, 1, 40, "Text",
      CellTags->"b:0.2.1"],
    Cell[12379, 410, 155, 3, 77, "Input",
      CellTags->"b:0.2.1"],
    Cell[12613, 417, 575, 10, 259, "Input",
      CellTags->"b:0.2.1"],
    Cell[13191, 429, 140, 4, 40, "Text",
      CellTags->"b:0.2.1"],
    Cell[13334, 435, 579, 14, 129, "Input",
      CellTags->"b:0.2.1"],
    Cell[14071, 456, 612, 17, 129, "Input",
      CellTags->"b:0.2.1"]}
  }
*)

(*CellTagsIndex
CellTagsIndex->{
  {"b:0.2.1", 19937, 643}
  }
*)

(*NotebookFileOutline
Notebook[{

Cell[CellGroupData[{
Cell[1776, 53, 26, 0, 117, "Title"],
Cell[1805, 55, 84, 4, 143, "Subtitle"],

Cell[CellGroupData[{
Cell[1914, 63, 29, 0, 96, "Section"],
Cell[1946, 65, 122, 4, 62, "Text"]
}, Open  ]],

Cell[CellGroupData[{
Cell[2105, 74, 35, 0, 96, "Section"],
Cell[2143, 76, 324, 10, 63, "Text"],
Cell[2470, 88, 104, 2, 51, "Input"]
}, Open  ]],

Cell[CellGroupData[{
Cell[2611, 95, 43, 0, 96, "Section"],
Cell[2657, 97, 429, 17, 64, "Text"]
}, Open  ]],

Cell[CellGroupData[{
Cell[3123, 119, 41, 0, 96, "Section"],
Cell[3167, 121, 221, 5, 63, "Text"],
Cell[3391, 128, 441, 10, 77, "Input"],
Cell[3835, 140, 42, 0, 40, "Text"],
Cell[3880, 142, 428, 9, 129, "Input"]
}, Open  ]],

Cell[CellGroupData[{
Cell[4345, 156, 35, 0, 96, "Section"],
Cell[4383, 158, 144, 4, 62, "Text"],
Cell[4530, 164, 110, 2, 51, "Input"],
Cell[4643, 168, 200, 8, 40, "Text"],
Cell[4846, 178, 195, 4, 103, "Input"]
}, Open  ]],

Cell[CellGroupData[{
Cell[5078, 187, 33, 0, 96, "Section"],
Cell[5114, 189, 68, 0, 40, "Text"],
Cell[5185, 191, 126, 2, 77, "Input"],
Cell[5314, 195, 54, 0, 40, "Text"],
Cell[5371, 197, 260, 6, 155, "Input"],
Cell[5634, 205, 65, 0, 40, "Text"],
Cell[5702, 207, 79, 1, 51, "Input"],
Cell[5784, 210, 131, 3, 40, "Text"],
Cell[5918, 215, 328, 6, 155, "Input"]
}, Open  ]],

Cell[CellGroupData[{
Cell[6283, 226, 45, 0, 96, "Section"],
Cell[6331, 228, 309, 7, 62, "Text"],
Cell[6643, 237, 946, 24, 207, "Input"],
Cell[7592, 263, 179, 5, 41, "Text"],
Cell[7774, 270, 425, 10, 103, "Input"],
Cell[8202, 282, 42, 0, 40, "Text"],
Cell[8247, 284, 388, 10, 103, "Input"]
}, Open  ]],

Cell[CellGroupData[{
Cell[8672, 299, 29, 0, 96, "Section"],
Cell[8704, 301, 288, 6, 84, "Text"],
Cell[8995, 309, 1472, 33, 363, "Input"],
Cell[10470, 344, 91, 2, 40, "Text"],
Cell[10564, 348, 312, 8, 77, "Input"]
}, Open  ]],

Cell[CellGroupData[{
Cell[10913, 361, 29, 0, 96, "Section"],
Cell[10945, 363, 403, 13, 63, "Text"],
Cell[11351, 378, 320, 9, 77, "Input",
  CellTags->"b:0.2.1"],
Cell[11674, 389, 59, 0, 40, "Text"],
Cell[11736, 391, 56, 1, 51, "Input"],
Cell[11795, 394, 509, 11, 129, "Input",
  CellTags->"b:0.2.1"],
Cell[12307, 407, 69, 1, 40, "Text",
  CellTags->"b:0.2.1"],
Cell[12379, 410, 155, 3, 77, "Input",
  CellTags->"b:0.2.1"],
Cell[12537, 415, 73, 0, 40, "Text"],
Cell[12613, 417, 575, 10, 259, "Input",
  CellTags->"b:0.2.1"],
Cell[13191, 429, 140, 4, 40, "Text",
  CellTags->"b:0.2.1"],
Cell[13334, 435, 579, 14, 129, "Input",
  CellTags->"b:0.2.1"],
Cell[13916, 451, 152, 3, 62, "Text"],
Cell[14071, 456, 612, 17, 129, "Input",
  CellTags->"b:0.2.1"],
Cell[14686, 475, 417, 10, 84, "Text"],
Cell[15106, 487, 137, 2, 77, "Input"],
Cell[15246, 491, 336, 9, 62, "Text"],
Cell[15585, 502, 416, 10, 103, "Input"]
}, Open  ]],

Cell[CellGroupData[{
Cell[16038, 517, 49, 0, 96, "Section"],
Cell[16090, 519, 109, 3, 40, "Text"],
Cell[16202, 524, 502, 13, 77, "Input"],
Cell[16707, 539, 75, 0, 40, "Text"],
Cell[16785, 541, 753, 17, 207, "Input"],
Cell[17541, 560, 67, 0, 40, "Text"],
Cell[17611, 562, 319, 8, 77, "Input"],
Cell[17933, 572, 44, 0, 40, "Text"],
Cell[17980, 574, 78, 1, 51, "Input"],
Cell[18061, 577, 41, 0, 40, "Text"],
Cell[18105, 579, 75, 1, 51, "Input"]
}, Open  ]],

Cell[CellGroupData[{
Cell[18217, 585, 39, 0, 96, "Section"],
Cell[18259, 587, 43, 0, 40, "Text"],
Cell[18305, 589, 297, 8, 77, "Input"],
Cell[18605, 599, 221, 4, 103, "Input"],
Cell[18829, 605, 127, 3, 40, "Text"],
Cell[18959, 610, 64, 1, 51, "Input"]
}, Open  ]],

Cell[CellGroupData[{
Cell[19060, 616, 45, 0, 96, "Section"],
Cell[19108, 618, 87, 2, 40, "Text"]
}, Open  ]]
}, Open  ]]
}
]
*)



(*******************************************************************
End of Mathematica Notebook file.
*******************************************************************)