From 59b88de195388e6ee1643624f5aacb2e5181d74a Mon Sep 17 00:00:00 2001 From: t-pollington Date: Wed, 7 Aug 2019 16:53:28 +0100 Subject: [PATCH] Cleanup of dev --- .gitignore | 3 +- README.md | 30 ++- anRpackage/DESCRIPTION | 12 - anRpackage/NAMESPACE | 3 - anRpackage/R/RcppExports.R | 11 - anRpackage/Read-and-delete-me | 9 - anRpackage/cleanup | 4 - anRpackage/configure | 11 - anRpackage/man/anRpackage-package.Rd | 34 --- anRpackage/man/rcpp_hello_world.Rd | 17 -- anRpackage/src/RcppExports.cpp | 45 ---- anRpackage/src/RcppExports.o | Bin 113888 -> 0 bytes anRpackage/src/anRpackage.so | Bin 75480 -> 0 bytes anRpackage/src/get_tau.cpp | 77 ------- anRpackage/src/get_tau.o | Bin 28144 -> 0 bytes anRpackage/src/rcpp_hello_world.cpp | 13 -- anRpackage/src/rcpp_hello_world.o | Bin 28504 -> 0 bytes get_tau.cpp | 92 ++++---- run_get_tau.R | 26 ++- run_get_taupersonal.R | 42 ---- tau-statistic-speedup.Rproj | 5 + testloadbinary.cpp | 328 --------------------------- 22 files changed, 88 insertions(+), 674 deletions(-) delete mode 100644 anRpackage/DESCRIPTION delete mode 100644 anRpackage/NAMESPACE delete mode 100644 anRpackage/R/RcppExports.R delete mode 100644 anRpackage/Read-and-delete-me delete mode 100755 anRpackage/cleanup delete mode 100755 anRpackage/configure delete mode 100644 anRpackage/man/anRpackage-package.Rd delete mode 100644 anRpackage/man/rcpp_hello_world.Rd delete mode 100644 anRpackage/src/RcppExports.cpp delete mode 100644 anRpackage/src/RcppExports.o delete mode 100755 anRpackage/src/anRpackage.so delete mode 100644 anRpackage/src/get_tau.cpp delete mode 100644 anRpackage/src/get_tau.o delete mode 100644 anRpackage/src/rcpp_hello_world.cpp delete mode 100644 anRpackage/src/rcpp_hello_world.o delete mode 100644 run_get_taupersonal.R delete mode 100644 testloadbinary.cpp diff --git a/.gitignore b/.gitignore index b66fee6..73978b9 100644 --- a/.gitignore +++ b/.gitignore @@ -3,4 +3,5 @@ *.RData *.Ruserdata *.bin -*.out \ No newline at end of file +*.out +.Rproj.user \ No newline at end of file diff --git a/README.md b/README.md index 0fbe69a..ec18b65 100644 --- a/README.md +++ b/README.md @@ -9,7 +9,7 @@ An optimised implementation of the tau statistic (relative prevalence ratio form * [References and credits](#References-and-credits) ## The statistic -I was evaluating the ['elevated prevalence' form](https://journals.plos.org/plosone/article/file?id=10.1371/journal.pone.0155249.s003&type=supplementary#page=6 "Lessler et al. Appendix 5, p6") of the tau statistic [[Lessler et al]](#References-and-credits) as we had data on the underlying population (i.e. non-cases as well as cases) containing months of disease onset *ti* and UTM coordinates of their household (*xi*,*yi*). I optimised the implementation of the tau statistic from the development repo of the `IDSpatialStats::get.tau()` function, leading to a ~**52x speedup**. +I was evaluating the ['elevated prevalence' form](https://journals.plos.org/plosone/article/file?id=10.1371/journal.pone.0155249.s003&type=supplementary#page=6 "Lessler et al. Appendix 5, p6") of the tau statistic [[Lessler et al]](#References-and-credits) as we had data on the underlying population (i.e. non-cases as well as cases) containing months of disease onset *ti* and UTM coordinates of their household (*xi*,*yi*). I optimised the implementation of the tau statistic from the development repo of the `IDSpatialStats::get.tau()` function, leading to a ~**76x speedup**. @@ -23,11 +23,11 @@ which is the proportion of related case pairs within a specified distance versus The relatedness of a case pair *zij*, is determined here using temporal information, if then *z_ij*=1 else 0. ## How the speedup was done -Rather than running `IDSpatialStats::get.tau()` function in an R script as described by `?get.tau()`, I isolated the C function responsible (`get_tau` on [line 403](https://github.com/HopkinsIDD/IDSpatialStats/blob/master/src/spatialfuncs.c#L403) from `spatialfuncs.c` source file for the [2782d6d](https://github.com/HopkinsIDD/IDSpatialStats/commit/2782d6dcc9ee4be9855b5e468ce789425b81d49a "Commit 2782d6d on 17 Dec 2018") commit), coded the four features summarised below, then ran it in an R script by sourcing it with `Rcpp::sourceCpp()` and then calling `getTau()` without needing `library(IDSpatialStats)` on [lines 6-10](https://github.com/t-pollington/tau-statistic-speedup/blob/master/run_get_tau.R#L6). I have provided both the R script file `run_get_tau.R` and the C file `get_tau.cpp` containing useful comments. +Rather than running `IDSpatialStats::get.tau()` function in an R script as described by `?get.tau()`, I isolated the C function responsible (`get_tau` on [line 403](https://github.com/HopkinsIDD/IDSpatialStats/blob/master/src/spatialfuncs.c#L403) from `spatialfuncs.c` source file for the [2782d6d](https://github.com/HopkinsIDD/IDSpatialStats/commit/2782d6dcc9ee4be9855b5e468ce789425b81d49a "Commit 2782d6d on 17 Dec 2018") commit), coded the four features summarised below, then ran it in an R script by sourcing it with `Rcpp::sourceCpp()` and then calling `getTau()` without needing `library(IDSpatialStats)` on [lines 6-10](https://github.com/t-pollington/tau-statistic-speedup/blob/master/run_get_tau.R#L6). I have provided commented R script file `run_get_tau.R` (in the root) and the C file `get_tau.cpp` (in the /src of the package) containing useful comments. 1. **Stop calls to R within C** (~26x speedup) -*Currently*: Previously the R function `get.tau()` would call the `get_tau()` C function on lines [403-449](https://github.com/HopkinsIDD/IDSpatialStats/blob/master/src/spatialfuncs.c#L403 +*Previously*: Previously the R function `get.tau()` would call the `get_tau()` C function on lines [403-449](https://github.com/HopkinsIDD/IDSpatialStats/blob/master/src/spatialfuncs.c#L403 ) (and internally `get_pi()` on [line 427](https://github.com/HopkinsIDD/IDSpatialStats/blob/master/src/spatialfuncs.c#L427 )). My `getTau()` function skips that step for easier reading here, not contributing to the speedup; so in essence the heart of the code isolated was described by `get_pi()` on [lines 64-148](https://github.com/HopkinsIDD/IDSpatialStats/blob/master/src/spatialfuncs.c#L64 ). `get_pi()` has 3 nested loops: over distance windows, then a double loop for paired links between people. The slowdown occurs at [lines 126-129](https://github.com/HopkinsIDD/IDSpatialStats/blob/master/src/spatialfuncs.c#L126 @@ -37,19 +37,25 @@ Rather than running `IDSpatialStats::get.tau()` function in an R script as descr 2. **Stop repeat evaluations of undirected edges** (~2x speedup) -*Currently*: In the sum over all people the same link will be visited twice ie *i*->*j* and *j*->*i* but this isn't necessary as is symmetric to *i*,*j* switching due to the modulus function. In most temporally-related or serotype-shared scenarios the pairs are undirected and so the summation should really be 'upper triangular' style i.e. . +*Previously*: In the sum over all people the same link will be visited twice ie *i*->*j* and *j*->*i* but this isn't necessary as is symmetric to *i*,*j* switching due to the modulus function. In most temporally-related or serotype-shared scenarios the pairs are undirected and so the summation should really be 'upper triangular' style i.e. . *Change*: Understandably by halving the loop we get a ~2x speedup. -3. **Split the `posmat` data matrix into multiple vectors** (~20% speedup) +3. **Explicitly force the compiler to use AVX2 CPU vectorisation** (~47% speedup) -*Currently*: `posmat` is the matrix argument that gets fed into `get.tau()`. As the loops sequentially go through `i`, `j` & `k`, if posmat is large columnwise then the next row's value in memory may not be that close to the current location and this extra memory access time will add delays. +*Previously*: The standard compiler flags are '-O2' and it is unclear if the loops were fully vectorised given the `if(...) continue;` statements within them. + +*Change*: Since the first version I have built the .cpp file into a minimal package using `Rcpp::Rcpp.package.skeleton()` with empty description and help files. The purpose is that the _configure_ file can set the 'g++' compiler flags to utilise '-O3' optimisation and '-mavx2' vectorisation which is common in CPU architectures; to tweak to your architecture please read [How to change and set Rcpp compile arguments](https://stackoverflow.com/questions/32586455/how-to-change-and-set-rcpp-compile-arguments). All `if()` statements in the inner loops have been replaced with boolean arithmetic. + +4. **Split the `posmat` data matrix into multiple vectors** (~20% speedup) + +*Previously*: `posmat` is the matrix argument that gets fed into `get.tau()`. As the loops sequentially go through `i`, `j` & `k`, if posmat is large columnwise then the next row's value in memory may not be that close to the current location and this extra memory access time will add delays. *Change*: Using vectors for each variable guarantees that the next observation for a variable will be next in memory. -4. **Work with squared distances to avoid `sqrt()`** (negligible speedup) +5. **Work with squared distances to avoid `sqrt()`** (negligible speedup) -*Currently*: To calculate the distance separation *d_ij* a Euclidean distance was calculated. +*Previously*: To calculate the distance separation *d_ij* a Euclidean distance was calculated. *Change*: Working with squared distance ranges can make `sqrt()` redundant. @@ -59,13 +65,13 @@ Although I have found enormous speedup improvements (constructing 500 bootstrapp ## Replication Unfortunately I can't share the dataset but can describe what you need: -* R v3.5.1 +* R v3.6.0 * library `Rcpp` for `sourceCpp()`. `IDSpatialStats` isn't required. * data = R `matrix`-type object with columns named: "ORIG_ID"; "x"; "y"; "onset" and no missing data. For non-cases, the "onset" column should be numerically coded as -999. ## Features not implemented -* parallel computations across the `for(i){}` loop for *i* in `get_tau.cpp`. I tried using parallel packages in R and C's `#pragma omp parallel for` with `#include ` but to no avail. -* GPU computations. A good starting place for rapid code development would be MATLAB's `gpuArray` class. +* parallel computations across the `for(i){}` loop for *i* in `get_tau.cpp`. I tried to no avail using 'RcppArmadillo' to enable `#pragma omp parallel for` within `get_tau.cpp` as there are [concerns that `Rcpp` alone can't deal with OpenMP](https://stackoverflow.com/questions/48069990/multithreaded-simd-vectorized-mandelbrot-in-r-using-rcpp-openmp). Applying `#pragma omp simd` to the loops after the compiler optimisations above gave no speed improvement; there were difficulties applying it since although the for loop over *i* goes from 1 to N the nested *j* loop is of different length which depends on *i*. +* GPU computations. A good starting place for rapid code development would be MATLAB's `gpuArray` class or ArrayFire. My only reservation is the code would need adapting so the distance matrix is calculated on-the-fly only once rather than over the k-loop for the size of the radius intervals required. This would help reduce the amount of GPU device memory accesses from the GPU cores. * separating distance calculations from the Rfun code didn't lead to a speedup. Even though distances are needless calculated multiple times, the slowdown probably comes from storing these distances in slower-to-access caches, when a just-in-time method works faster. * **Have you found an even faster way to do this? I'm open in principle to pull requests to this repo but message me to check.** * **Found a bug or even a typo? I'd love to know! We're not perfect and we should all be open to criticism so we can do the best science for infectious disease modelling.** @@ -73,4 +79,4 @@ Unfortunately I can't share the dataset but can describe what you need: ## References and credits 1. Lessler J, Salje H, Grabowski MK, Cummings DAT. *Measuring Spatial Dependence for Infectious Disease Epidemiology*. PLoS One 2016; 11: 5–1–13. doi: [10.1371/journal.pone.0155249](https://journals.plos.org/plosone/article?id=10.1371/journal.pone.0155249). 2. `HopkinsIDD/IDSpatialStats` (development repo for R's IDSpatialStats package) maintained by @jlessler and @gilesjohnr. Code from the [CRAN read-only mirror](https://github.com/cran/IDSpatialStats) hasn't been used as it is several months behind the development repo. -* Thanks to [CodeCogs](https://www.codecogs.com/latex/eqneditor.php "CodeCogs LaTeX equation editor, just copy+paste the HTML they provide") for renderring the mathematical formulae, let's hope they don't close down their site—Github **still** doesn't consider renderring LaTeX in README.md is a core functionality of code development! +* Thanks to [CodeCogs](https://www.codecogs.com/latex/eqneditor.php "CodeCogs LaTeX equation editor, just copy+paste the HTML they provide") for renderring the mathematical formulae, let's hope they don't close down their site—Github **still** doesn't consider renderring LaTeX in README.md is a core functionality of code development! \ No newline at end of file diff --git a/anRpackage/DESCRIPTION b/anRpackage/DESCRIPTION deleted file mode 100644 index 88c7e8b..0000000 --- a/anRpackage/DESCRIPTION +++ /dev/null @@ -1,12 +0,0 @@ -Package: anRpackage -Type: Package -Title: What the Package Does in One 'Title Case' Line -Version: 1.0 -Date: 2019-08-06 -Author: Your Name -Maintainer: Your Name -Description: One paragraph description of what the package does as one or more full - sentences. -License: GPL (>= 2) -Imports: Rcpp (>= 1.0.1) -LinkingTo: Rcpp diff --git a/anRpackage/NAMESPACE b/anRpackage/NAMESPACE deleted file mode 100644 index 0ca1d67..0000000 --- a/anRpackage/NAMESPACE +++ /dev/null @@ -1,3 +0,0 @@ -useDynLib(anRpackage, .registration=TRUE) -exportPattern("^[[:alpha:]]+") -importFrom(Rcpp, evalCpp) diff --git a/anRpackage/R/RcppExports.R b/anRpackage/R/RcppExports.R deleted file mode 100644 index 1bbcd19..0000000 --- a/anRpackage/R/RcppExports.R +++ /dev/null @@ -1,11 +0,0 @@ -# Generated by using Rcpp::compileAttributes() -> do not edit by hand -# Generator token: 10BE3573-1514-4C36-9D1C-5A225CD40393 - -getTau <- function(ORIG_ID, x, y, onset, r, r_low, index) { - .Call(`_anRpackage_getTau`, ORIG_ID, x, y, onset, r, r_low, index) -} - -rcpp_hello_world <- function() { - .Call(`_anRpackage_rcpp_hello_world`) -} - diff --git a/anRpackage/Read-and-delete-me b/anRpackage/Read-and-delete-me deleted file mode 100644 index 6d6236d..0000000 --- a/anRpackage/Read-and-delete-me +++ /dev/null @@ -1,9 +0,0 @@ -* Edit the help file skeletons in 'man', possibly combining help files for multiple - functions. -* Edit the exports in 'NAMESPACE', and add necessary imports. -* Put any C/C++/Fortran code in 'src'. -* If you have compiled code, add a useDynLib() directive to 'NAMESPACE'. -* Run R CMD build to build the package tarball. -* Run R CMD check to check the package tarball. - -Read "Writing R Extensions" for more information. diff --git a/anRpackage/cleanup b/anRpackage/cleanup deleted file mode 100755 index 9a8dcf9..0000000 --- a/anRpackage/cleanup +++ /dev/null @@ -1,4 +0,0 @@ -#!/bin/bash -if [ -e "~/.R/Makevars.bak_CustomConfig" ]; then - mv -f ~/.R/Makevars.bak_CustomConfig ~/.R/Makevars -fi \ No newline at end of file diff --git a/anRpackage/configure b/anRpackage/configure deleted file mode 100755 index 1839ac8..0000000 --- a/anRpackage/configure +++ /dev/null @@ -1,11 +0,0 @@ -#!/bin/bash -if [ ! -d "~/.R" ]; then - mkdir ~/.R; touch ~/.R/Makevars - echo "CXXFLAGS= -O3 -std=c++11 -Wall -mavx2" > ~/.R/Makevars -elif [ ! -e "~/.R/Makevars" ]; then - touch ~/.R/Makevars - echo "CXXFLAGS= -O3 -std=c++11 -Wall -mavx2" > ~/.R/Makevars -else - mv ~/.R/Makevars ~/.R/Makevars.bak_CustomConfig - echo "CXXFLAGS= -O3 -std=c++11 -Wall -mavx2" > ~/.R/Makevars -fi \ No newline at end of file diff --git a/anRpackage/man/anRpackage-package.Rd b/anRpackage/man/anRpackage-package.Rd deleted file mode 100644 index a23ac20..0000000 --- a/anRpackage/man/anRpackage-package.Rd +++ /dev/null @@ -1,34 +0,0 @@ -\name{anRpackage-package} -\alias{anRpackage-package} -\alias{anRpackage} -\docType{package} -\title{ - A short title line describing what the package does -} -\description{ - A more detailed description of what the package does. A length - of about one to five lines is recommended. -} -\details{ - This section should provide a more detailed overview of how to use the - package, including the most important functions. -} -\author{ -Your Name, email optional. - -Maintainer: Your Name -} -\references{ - This optional section can contain literature or other references for - background information. -} -\keyword{ package } -\seealso{ - Optional links to other man pages -} -\examples{ - \dontrun{ - ## Optional simple examples of the most important functions - ## These can be in \dontrun{} and \donttest{} blocks. - } -} diff --git a/anRpackage/man/rcpp_hello_world.Rd b/anRpackage/man/rcpp_hello_world.Rd deleted file mode 100644 index e4f90bf..0000000 --- a/anRpackage/man/rcpp_hello_world.Rd +++ /dev/null @@ -1,17 +0,0 @@ -\name{rcpp_hello_world} -\alias{rcpp_hello_world} -\docType{package} -\title{ -Simple function using Rcpp -} -\description{ -Simple function using Rcpp -} -\usage{ -rcpp_hello_world() -} -\examples{ -\dontrun{ -rcpp_hello_world() -} -} diff --git a/anRpackage/src/RcppExports.cpp b/anRpackage/src/RcppExports.cpp deleted file mode 100644 index c909457..0000000 --- a/anRpackage/src/RcppExports.cpp +++ /dev/null @@ -1,45 +0,0 @@ -// Generated by using Rcpp::compileAttributes() -> do not edit by hand -// Generator token: 10BE3573-1514-4C36-9D1C-5A225CD40393 - -#include - -using namespace Rcpp; - -// getTau -NumericVector getTau(const NumericVector ORIG_ID, const NumericVector x, const NumericVector y, const NumericVector onset, const NumericVector r, const NumericVector r_low, SEXP index); -RcppExport SEXP _anRpackage_getTau(SEXP ORIG_IDSEXP, SEXP xSEXP, SEXP ySEXP, SEXP onsetSEXP, SEXP rSEXP, SEXP r_lowSEXP, SEXP indexSEXP) { -BEGIN_RCPP - Rcpp::RObject rcpp_result_gen; - Rcpp::RNGScope rcpp_rngScope_gen; - Rcpp::traits::input_parameter< const NumericVector >::type ORIG_ID(ORIG_IDSEXP); - Rcpp::traits::input_parameter< const NumericVector >::type x(xSEXP); - Rcpp::traits::input_parameter< const NumericVector >::type y(ySEXP); - Rcpp::traits::input_parameter< const NumericVector >::type onset(onsetSEXP); - Rcpp::traits::input_parameter< const NumericVector >::type r(rSEXP); - Rcpp::traits::input_parameter< const NumericVector >::type r_low(r_lowSEXP); - Rcpp::traits::input_parameter< SEXP >::type index(indexSEXP); - rcpp_result_gen = Rcpp::wrap(getTau(ORIG_ID, x, y, onset, r, r_low, index)); - return rcpp_result_gen; -END_RCPP -} -// rcpp_hello_world -List rcpp_hello_world(); -RcppExport SEXP _anRpackage_rcpp_hello_world() { -BEGIN_RCPP - Rcpp::RObject rcpp_result_gen; - Rcpp::RNGScope rcpp_rngScope_gen; - rcpp_result_gen = Rcpp::wrap(rcpp_hello_world()); - return rcpp_result_gen; -END_RCPP -} - -static const R_CallMethodDef CallEntries[] = { - {"_anRpackage_getTau", (DL_FUNC) &_anRpackage_getTau, 7}, - {"_anRpackage_rcpp_hello_world", (DL_FUNC) &_anRpackage_rcpp_hello_world, 0}, - {NULL, NULL, 0} -}; - -RcppExport void R_init_anRpackage(DllInfo *dll) { - R_registerRoutines(dll, NULL, CallEntries, NULL, NULL); - R_useDynamicSymbols(dll, FALSE); -} diff --git a/anRpackage/src/RcppExports.o b/anRpackage/src/RcppExports.o deleted file mode 100644 index e349a8972e3d0a2d4b8e236280bb06e8c58d3ae9..0000000000000000000000000000000000000000 GIT binary patch literal 0 HcmV?d00001 literal 113888 zcmeFa34B!5`96M=$v`BKiHaDlOT^IPQk$?xToOohqJu^Ph-(NTfhZ&~nFv(TXcA=_ zqG+{s7Zt5-ZC$E}d)&aixYW9XOD$tus9L2i<^Q~Ax%14Oxrkl9pa1Xo``%A(?!C`> z&wJkUp7)&dp5@*-yCQN}W=2N9F-L}Th?8L-=vvU8KuZ-~>g<`ZM~0xs;*{%5b3M*n zk2lv7%=JWEPr`W)&T|cRGOjf^m*70hV6$+o#d$W)4xGY!HqJVnb8yz<6s!T~T%7Z8 zHsTcQ6r87;`}w#oF!#;4wwU|VaBVgBZMa6weLJqFoBM^hE;9Fvas8>eKLgh@&3z}X zUFN@m;1xL$zsLi4;F*NbpoY@RQ{^-`R_Fwd9cdIioaabAVy0>XGS9!l^=6#E!})uI{RY?H;=IK?--_#PIBz%4ci{R5oOhb% zyK!BO^B(hjAFlV~e84IK9B1Q=Ke)< zeF@h;nfsT`_0PDzV(wox*Vk};-Q2%nu5aS{mbqVVt{ZTD+uXlnuJ7Xdp1J>vxxSC< zU(NkST;t~c16=#f{fD@IWbQx4^%HZy3D>`w`%iKG4Cm)K|BmwuoL}OU^4pB_A2`3l z`8CdOaDIzZ!sMS^r9#Oy(_91Qnq{s*a~)u=1I;ztT(>jVLFSrcu7k~Wh`HwC8Z!5J zxNdLmcffT=b3YWCv(C*StCY5;&80IFuS-f z-2FwUORPT}`#K!^vNHC$sM85|ZVrWGkB4Jfg$pO(A*21Ive2rm!lje*z@RYPvtu|m zIeSf3rW5V~pE=>tnrr*P;m*hMPdKUU#Ilnf4aeRuejyxNuQ*l;w}6Tq1C!VXk=Xl@ z*e5Ei|F_}pzUWQ;ah-PW!6JgXuP*)qf>l7y-kk(2eo~Xf)}K)Pg%r9I+OyPYI3TpA z)NJH;M%i?e+tQ@0PAKjRtqKZBS6}7dZA1HmQSJ606RLFrh8 zBtct#Oe}oP0sn2}GXtc$fOO%8V~VDx}`m~dVCNl9@c8c=2{thOh5N1LxmwH1&4c#y{rxQj8Y*FWrB-|a*`7-R+&<9 zA*6|`Yl#**EXFZ$iDRroQN+iKN>j*9wB8d7^Q@jw;ey`?tz~JX?>Q4RosZ^0{BX(o zB;KQhw+)wipBx@G^%NofA6vemRV^dEnu))StW|b?lo9H#RpM$9R0h<);w|ydCHGdg zp8TZj$A6$RiGqs{kXzKV%Geq_Mh7M8PbBtP{8hz8N;tkpLY;|(LVI>V=c+xZjD1@+ z4Gr6LQ^TtxJ?&1U#LMcqHL@!fw=!=q||?vPOkqFdAK&Y z|C%0-{iW>K%Gj2&8I`dw%BD|=Jzh3#O6-~H@Hroj6dUNNK;L<3pU58TU9B$Nx8Xva zUBwey?dx965;nVw{e^o1;aK+*h)#C=j{$B6)ce3XghgWi+<2!nD7_DCP>)}1yip!| z?owj(KClVs#)~ziMCyIOK^bp6TSv+hGl+e>u|+=?sK@s=*5WbJvv+uH#Zr+a((^kV z^4CbLzQe5<)sQwaaX%>nt@=o8K?g~jGZ5jiPi(E*FC42_8eZGI2#IK}H|pd+iNt=_ zB0>o&RqhVOc<)5F5Q*I=1jPFjXac22ycJ)tpO;`SIfVdVKa1J#orb*-5c_rLM?LI? zFWAq-#Yiu@)$FU7{Tz${BNLOQ7Xo5`dn$Y33-;jzdtIO+eJQgan})p*5c{dA$`QU` zFV)^iugyo;7cu*WH0*_d*#Cj0hmi`?aBoG&mb-3)M#8;ImRE&)dnL=^-b+QaaPM8h zBGTL2qAuN{ZC#~%>X$=TGmlS{b-3pNv3l1S7O1MVBfO#L?zX@XNY%->V43|6-&2f8DRLPzPBn#n^wQe%NmUMW42MSV&9>`0B zIcYF4D=$B>1urHAC|HV{`g6E6red%ehI^(w;X0|2*oIJdqck_Y^{<>DLxQrC%W6(S z4<&YQln))Y?-TA(tx%+7L#XR~mE^imSDV}zmn(iU9uUtI9Mm1^5=b$LAA>t|QRR}7 za;tRdSU(V^Wcu28;A>%1Dout-Khz~-G1y0_>pBEl!Qz+Tp<%gwlomqRZGr{iJyl#;(7|B{|4;%rWM0UmVt_BPz%@Q>a2x7zS) zID8o%bo^G8KKR$F^o=5pcH=*1<3GUhe=G61{H4r}`xZ8|ZR3ywK$;T&v(4Ka{|S($ z#NQ1MuAUC?KQUlpB@Fm?TNpD>}nCw zrl~Gu;zZNT0#N-21FG|2V3xG%45&VX0CX7iqL(abodg4_hhRW;4=$h^ z;|tLk>x}h(tgrDi0ii)re(%fdT~CR9I3@OX)u~RL80z}eZP7O-n7nP7Aj5`lx4uG8 zMhl^?uidDB!<~1W9O+%M!WO16r$}te6gal8e<)3yF&sP2?^&{hg7sOw&E(#3Nd?zolEQfodso<@2eFeVr3y1(ZHe3U^wd{{rzA^v9J%?vA@?FnND=fl2$~mc)*;g^9*uZel+1cUV$@vMp z!npqgHP!h)M<4E7i&?OA zXSNY796Pqpif(LMsaQ=u>$~1F!kS`oSSvQ+UsZ>DPy#bIT4otvCz6;PwrKAeNO4EQ z!wdy_zcUrlNwWU_&ys({0sG4GThdYOwm>gpjz0Ej*i55@W3QYbwvJT?tsjJTK9>Ec z_xZcxcl>lqOJz?{bZ+;Pi}nlGJ&E;%?kA)99fvsWgZj5`Sl)5aFekcQxbu;WaLH5c ze{)v^u&fcCsH`ORCDMc;M+@zk4R??Q(E1a?u}4%CcYez-KsgNz2NrGYSy<5dTn5#E z=?|c_*rFCwXvRgYq?CLdEk?B1XQfM{g)pMhS>cJwCfW}S_f9JusGJB{$WSAc*qf%A zji}WVJ*ywv*uJ-3KfxHOLlhnE1+RPIia1-=$udl+yubfx6;nbT=l^f=uzhc(J^QZq z;=9_z|62V2PW?|7w_$mry(w*1VlZELKnX^XZtHqAZoz?vyFbkGo4e=c(g@p>rSnSw$|3wV-E!6t7)r`BKw58qNe$?;ezP)6S%-g;p#kRhCq!iE=$6 zqllMEI(AmKvSdNhlzj4vbEsSP0q9oKPBA^|w&!XI+V$^yh8D`QbvPzMqp!kRmR%c@ zHSOx`_-F{J)-y4A&@D@`$Xt961$N*G?Qvp$N#+-hMbNgDd3=A6ri{W$)$c~i!HLBOq`yX))@o~kZe?2(%zh+Rd(C8`{ z=04YC?D0sJ^IXP&o-F6G%;2Vgb8aBGHp_V=5WGFh=?@6<&j75J4yrqHH#D~{sEr<2FuAr#sy|wU z)9YKK1r4pu3kqsm=e94XZ;BQ~7q`?GM4Jm5o1*n|>st#Dsi3{BzM!@Z?L}*2?d-x9K&ZeczRh-+uHlzByOLlM-gX1T4srKPzwTCbXBjz6Nw8)@tG zlT`t`%7I$ynp;s~(fNxDqV=r{8k=g7O2O>K(B_+vbHwZh3K@QwsF3e!?u@d zf6UCJ0uC4Vb(D+>7V6qNfJNvv-5c;wP}wu4Fx-1wMn&;Q;ocz<9_~48V+^NGCJr2PH~)}PctZp@S2m|D7|m<5OlW6p5pI&T5mu--}QDvGxp zwB5q82&6CWo074iKgZp_5qr8ftE0Dk(z?)t)fwU5BQtt`x~{Zi%XX*#C45<#T1u&m zmF9(3Au#rI$5;F69X1_b9RO9K|IsAHNaW6g4p=xMd=Idyj54?14To-cxGh%+wF=4t z(di#{e6=Wa{zIZot<=z}7of9?l+HxKWeutiE2MOnxFg>OR`rwqHV^3mlg>*y-n0H`S36D!TM z$_}mCF@6(rQdyTVv8!-75)AeJ45YnPD5k8^ve@>Oy;-pA3E}Q7p{^`Rx2H4@`?5k^ zA6}{y&mfN2vT*NAyh|z3YABbt->4+Gg~8MirM+f zh&9#UK~?;M+yD0ty*S+Ya7OR6Di~x&?;>rB@o*jXr9B;g4f44Kgz{9SGU9fiLu5mr zdW_enC6f8Z&L^@o$)vO-!l(0zOii*U9>h9e8mO1bYBFIp@yC&QSA_nX6SsXHYB9F> zFJ*$ti$*M8V46f zu(J%aC#6MZyoz-n2>W1xlv-@P6uJr@To`YWOh~T{p27;<+%(~rW%-F01w>+CpCG8h zJw<3GiXwU>jI}n!0ilNdS|icW!ZHkSC7MZXXkq*f$Xuy58&-+h+%`)LUAK0s%_{qH zy|KZbs{D?3j_!P_AXb$hKUH$v`DiBFyI5U%d;c6EuDw5vhb`-rao4Sj-z9U4J+~BU z`}`=>H3jC@u|(|i8mOpjO_!?D8S&9#p#PS3#;a1=HtpmHE!C~NA+%J7XpOB`YA=l4 zxz$r4`AZaW1&qa%!6#eR#m7m*w`y``xVJhp9DA{IV-Xs7)6NFe)9IfUzYza5*4U86 zZ76#fR3%bVy&)mUj$t!2=*8O^W5qvK6YoDBOlAa!R_g|R# zx^RzpW!WdlUX~6K5VUsYTX3|II%Zqrf|mL96|Jq!txmLc@sI1(9fqiL8f7HnxZ@9z zL!=1z9I6G>JM?7cg;%Y4{vi^>yp5T^m>9sqR$_4*%V}o6aS3+#J1W9HD(2lVM(g)8t7>Ahc#bOan} z_G`uZF%1`573k`V?t9RP_8){}uZeN)^Qkh(0SQ(29MXn-+Vf&OF3C!9EDayTr%qK9i4Xr|fmMzsI z&2W#l5cjQ}Q1|wdnR6Z!L!xTeu)MPO>=xZjXrqR+TLc;hujQ=mDjX@dIDi(*tzZ)j zCikWGUR9N@LV4T5HYf-V+hE*R+54VSw1JB@hGRxnXtT(~iDbaz#h4U}&aF2;bERT-(_j)r&>2o7RtyeoRhgP%!Z7z0Xq*tq7-Dg6puQ~ zawGQ?y|VZ4S(0E8!W1hR$w|-vcf*Z zn9_~;sAHhl2vi-rL}A#B1G9%F-n)J3$))FE-6% ze`)XF0Uag#wdbtKfOsYp)jC|V4!(qmhWfq{u@qw!T{7eD-BKM!xW~-z-LEVFl9Y<@ zTmd>C&Tj~ZnjRJ@<>m=7nTF8(hs9<$2z!!OF%=_MS_CH%zzAW(U|J?osX|howz`+1 zI$e^bc)GT)f<+$1C4MkQc4~hunr2BuYy$CDAer{#2{`4uqG?tfKaM1b*%zA6#|tnv5w(lq1HNtSE%b(0&R;` z{6Yv-(y{w>DzFMA53E8`kwa24R7a0(Zc ztge;uUto?wuTQ;z$->UBp!@5Q>Cmb!CBRbbJ23jRwa!)@QOlaFot5I1K(fVdxsX30 zT0jjj{3aNF>i}zhiXaWYBFnE2r15hh+VXHLeT~(n(z^&&r%NW@sP(ZG`l7I*SM2u& zd-O+W%Qhi3vB<9L&!qtEh6EPv|1hsy-Ylsq|Ns9ScN!XpQmPWEn&1$TRcOt~Z zi}d_fl0hxNyNpPtOeTDyXG?>3ocKnnKV8%LHfS=IDvv4t>=$+-Zm^$6fXLa z=uDe$xc4gEy!4)m(pVpk-J%}Zs@I3_u?sX}oEs!E>_{|nZ&OqrvcvbtlAA19g=1kn znT4wcJmrP&(W_Lot$4~`TUsD%a1F7l!ZoEoY!gdBt+>Q+)h(H<9UB3xi}+}`<;4TAg41|>3|Hh%wuJiu~QGE+K0TLC&BK{sY!+_g)~@tLA& zp-QyIXlF!wZX4?g>lQ#-859?Bq1^*kuZwz#6mNev-;Qu znh3n%B9k%USJ}BHFH-VlsC#cD5|@AdZy-3n3+~YKz_~Vi6;Oah;(1PRZKFpzm}OZT zz$N|_j6Fq(7-RkPgNzt@iDd(VSMTiHH8A*PzVqb3l27uTr8&V*hhgcZ;oY5`mv@Yz z^?7}#;D()?%{%Rh?XALeWv76jcxI|F!Y=iongcCN?|g30yy z!Kd?`2i5(W{NP9V&hK^+!aH_R~L~-sl;>pY!2Hh@UhPdj)`O`Z}?2 zRyc;$lHN14KONEQ+{ItOmAefzXqT0w;x^wb&o~fty9xrf+F%miXYQx;G-OZdIV(@~ zoxq>4G_X2py@ht2!79;#ZH&a`K#UA51L#SHh<@cALK){-e7t?S3Of*4uFG+)?5%+3 zTwXW}H>FbIiZiE&_ z8BJbgBKprx$cHL{ow;; z^oWgGEmR{Xt;w1IDP%|%6;~E-fT}U3)Ozt_K0wTH$(fK=>Px6=hX9_%iJ5{{q1_3` zrb|BZLR~=#5I*2qpnY((*dwiul&c~k;-@?fG%o+#{;M}qSSlIoorZB?3`1C~$drb< z&H&TMJ{$Dr9W2;j2pl^9m#_uxyw~H-W~)jPmn?;mr2`>Kw6IiV$L3^%M6)>o^|Zop zsH;x0mo7UqB*0{+RI*bg9Z=YlozjGJYpqz>7Ygt;kb3N_vUiHdR~DnuzA`K`!OK=6 z*`O=pBj4IKWjV~4KkN~?&6;It6`eL~ZGN$h(mhGV#Q=GzT8s^K%kE7Kkfnap02%(1 zHVHALuJ^>fkT1p1{ys6yE2T+{t?Qp4{$Cr4$-L1t2raAiQ_wd^ z9~rQdeGydUJSm3Kc$xB8%0=VYcLE;yTfX;C5doqGI40isQ8GJOyS)Uoa zEyuY#yY}3{&VxA+{+hwVaH|eA&8Id{sU0ds)VH9GP+O2xlFB^b2 zwxq?kxC}`nI;(x%E!+=FZxd z7VO_~4lHmEaGb-`jvd*})X<1cJwF0z{ecJClt%6&{;eEk#~Mw$AA&)j44m}VRrm^k z;y#s;M4d9e7%D~O!~rQ2=7uX?!m-Up&q%|-9q%C?3LrJ=#F3${u58p%nmE^;M(j14 zKH?@%;4^J3u0NMU#Bjzap(NtDjK6u|>*H}#FM`B<1@2%UJ%83W?%s%RJ zrQJy{b5wdXP+2kn@ApfC?RpSu*CsTq%S>xSdIvKdY$i^yc;8d3NK(p@6_|*wcQm4^ z2(LU9vvc!+P`8Y`Yz6AvjJFA8t=9t3*oEI@DUF3o;Le1)_2)74q>@aoT`&Q2E$V#% z8^=|?UGPD~yioT?Xn-w=iWu5TM3GSotVrL0waY^r7BK{W@FFto+v!<*7?JOC; zsUB)NC{aN~wM7Z4hejf8qntAfi{u+tR`r&(IBc^}Q^gytQj=}Ufy1=6bl7kU7a|u{ z8`1^C%a@5*S#rU@aOnIs;%KR<*bk8dcym*nr4cwo@Cd`=cfo&nPDp1)dgm9m=udSi z2OV6MOQyVYaL(3!TP7T^Jn_0Q2r#ls$1jgE`0qfL!uM97NpG)>Hilz);d`-@Qnl&f zSQQ+DX`n;*sbVZr$*qRb%W-yiKN<6jtmjRJ;hPrZUTr9*hw8n=O!-Py@vuVDL%Pf5 z+Jn`lD)`mUhq~U8Y-3tadYYt z?P*SvMJzqN=Fa0$decpMinvFo*QZ#S^pMnQd4@)8Dv6AjVC05CG)KY)qUGy;S(~9z zb5s`3W}sYWK#hr<$}le+I|`VjdK7vdOti#~!aXu&{A81&&bZiya@y_p8l*fFtx}dJ z=*+^F%Ff3zEDs5tkUFqgR*F)Pr z#J`dk#KCfcRD`gQVL92hrD9ado~oo$%#ug!J99f+ZEgUoq>%ICxo`=kN+BObsJh85 z1*)+Jn=K-2qr@6A{ zkoYSo19YVE(ejcxuK&cdG5f^(3us{DBT+qZMRObK+7I`Mf-`eX$>blP&6as12#wzM zPu2|0)RCF4Gvkx{RTY)@s1q+yOG&3#owcWM=i^xBT_5Vc2dQXjW9y`ygbAQgh<__t zSLycO%?&on1#qZKK9EN8OV(fw@j`GzcswHg4W5Jg8REyM>fFn1yP=jDe`YQcXCs=D zCw?75oh?`1oRM)&#(qC^ocRU69O~SVG2)gToIhs#@bw*>Z!+dyIMn$lGkEP#=h;9I zp??q5-MNFaK5K2}=Q}!&4h-VqH`(J5{7LrC4-9p#+HT16L!GsQf_DrJ+?f;Huw!85 z;NY)@27W&z`0fsYzvc$d$#cHS9SQp1hJxVmWGD#wb7kAiuXYF`-IsR|!Y8HE`(#In z|IARKI9q?C)tOU|Ee3Pv*Xy^3=6c>9n!MkB6>3My0^@><2k+oLBlU_A%**aikxrQgQ*RTYI)=_aOF~o^ zb_u!#>=<5ZYB}mVIxgu%=%papOQa@mR1=^&!nOL1BGqU~qKTORy>YZO&+cjM&Qu=P ziG|rydfW3NG3gX!NgycwQIc>NAS&T5a)k>ege&dBvKT}wTaz^aHSk-s>XDuiP>F*{ zj-hz6JJ?^Jf-#sRll>rKI}hdIpR%FKiI_1+9wxv=Fy{Vr)nlrSh&7;K=&g4fa$o9@oPIWva=2E_w%8 zfVa<$6fS~WJ6KIWDYtgHusr86sek0w+WDi{oBbVFA1TK?-q5nx`uKRXaEdse#lth=*vqN2f0FC7lzX!=#|9$K7;!VV4d8N0i zidR^9b>+u)9o+(;H<6|7uVn#{h5)1=8QvZO0V{i_=fxX9XXKEP!zXx*#5VOWP0(BY z$NJ-0NOQSlR}@gDJ^^|k7&*ls#kZhvYr^#Dbx}r;TxTQ+yt(*Ebh(L%$NnxiV}@VsmHxg|tSt6P|F6;LR!-clkUhF^Y-9gL zs4e}=rc4~8TIDwr6}|uY{`dOdktS0bDRT0p$e&W5&=t3$1dfA#<3+uIfmfr{r$_pG zFhka|B4fZcxz4kh!Ec8+Z)FCz3~@GQHk^~|d=eOOd9HJL*0iUGIG<-lHx6;G9T5Cv zu=CJ>0l=RhxCg!*o*leB*SSACcvG&kF?+y0xlYG+!FPu^S8R9U%E8XLIl&b}oJVqk z>u^6fcqi@$qkH_1^w;yG<4fuil}nEPx=BJ*e?2YHC8`!db&8(;dL3xKSAT6QXO-=r z`;=cN8M93RQiuOwyXXsCd1OD3@gF>0^jBzl{*UXT|AdAjVTrnnei)eTJ`V(2x`4{@< z1-8FaYV!8a@_~L^p|+ua4nt8X{c|Vrh~JOoY=5?Od42llu$5Pqa=YcVKh@vfl-%DA z2bDCws}uU6{&qY|3&fpAKukyn-a5qj*>=HAInE=4f*-0L z_^QFu19uE@`f~)oKsu_nXmg8$^{@I?_q!rb3X0O$4r0A7`UpAH!E z&_L(PfrX9}bh4HNmSir;SfYN>1$EkU$Xex78NqvoI3H#PUmWbb9tidgcGhk;`?4X< zbAyk-50DHQhFi$Gd&W*b7>E5m8N0+Y@R5k%HyMHVGRCN1g;@foU^ITO;MGBApYh9s zPV0b6gU&S><1s$#%sdppiY;GTl1a%L3@D1&vW)btVeV9K&;Pl_5hz72JZo$FAm-V94;La1kGJTf}r`+kRU?; zJ|qarOLK#uxhppans;)8?_+ZkUZL6{=bcdShFx;5&%@5AoOkkq@8f>^;Jv%#yt94q z!d-K&-ywMJ4|1;Aamd|eUd%$BRU zYAlD#jGU}8=mj+1&E&{f&knM&5oj}7b2+dHC_@DS%a$7)vO_7)+}nL1VVS1}o^NTe zH>g1RYG4HhGrO5qYN58Wp-FE%lD$mpkepmM;EJgt%wXy*N4PJQdzFM)-7#B2Wai!s zHlRXCpm#QH#;r_Zn4L{|l0OkH+hkN_036k-fZ5#SAj}=NsP&Bwi4|Lq0pcRbl+0j5 zCbfv7cX5M*qA@aImiRJ5hgE=Og-hkZtbN?0I3U7u1?tULIyspXTbtNBg)-B7r#67n z-8;1k@n!Fn=t%ZXy&_lHJM|`Rlnbx+PGQ6)4^S5-*v;mt#c=Xz^ArXNcx%dRp85&y z>@K5HHGUJ-$X25XX7`kEM8=Ek_9^|slx&~+8gXS$1*XW6qW+FX+WkhuWDFoq4RGN zGgRZBDDG_5W~t;QTQl7qRG13Ih8eLV$;`(R&t?YX(+X#Xx=xcErrSX^M*?h)szj(+ zlB1cP991PIJ8jibwv2KVYxCIj*0xaTijZJuE3dnYO53tBD7BGFEI$&{wJ#F8KkGKp zb#mR2>ZC1DZ>N%1h1_M@P`5bin9P#qh9H&nUhC6gq}#A!;i zpXyuXQ+Fsm!EY`|-cN-`Oh=n}E$@D+0zI!aN^vNc++f<@INv%OC^R)Q1{1^1tJDou zm>EN@i0=WBW!_E7IL&^lS$b2{%*3WBSgzerRTUrVrk}c>s#HgvkTz;jTt03J{UzC{Q+Eqo?VrIXQy4rDs79e`H zDuMoA@DK^JImx@Kq|_ZpHl)`zFKEFMwf>?LKIqh1e_DHen|u;TJ_2=M!3lCXctqR8 zf+)U|R3F9dNw!B+za6E%npQyFASSoWej)rqO!P=Qs_fjHzc8%Ub7K#OUwr4wNANRc z8PRbM=0XrlquO;C{unipi(j7Ve5^n$UOrJKc7_S7H&pf#dj=lT@y%8)*`aUOL+y?W znIg6It#=9RR95@3P=-4_6p@Ap4ZX44okTO3IxuUo8kL1KqiZ`!7kf#F%+y?E=VP)O z*Z!U&7Y=gR@802UOFN&f22|$fxlH;dT@zwPM@kAE!ULFI)wRmId*QX}R+6##i0OW} zvu|X$_pH)2I^0}B`>02A;i*wqHA(>eTA*g6U#d=y8QKSv&hB?)HC?a#LYYyULQIZK zsLql8^(UKirE;RBJT`=yB?%q*zn;nD+3k+~2II}s?=mOmSgo=XR9aVd z0#WM$R`nf(QtYY6QDfw*j3w?T94k^uk`PT_TqIi$Yt-hHt;aC*WANcPd`i-`8~$-L z{RFj%v8i_csG>@{ow0sSh5gQ?e{{#nHx*xStFK)!ySm`6c?TCuZ_DK4QM(uq`@C%HaP#{EuiLyp5^YolN$+jnPub{?qCGr}?D_ z>;e9#)3eDBv^W2)>*13pA6PJQ`t0_mXnVo<1I8Rs^yA`o^;CRT@q`14#vD+*kH(b* zLya66Kb`7iEXvFH!Qg?}%i+t*co?r_sRJV?bqH2u?t|tk7~CFpFs24KKpC?ie*t@S zq*KlrLXXdQjY}`uNRf6G27q(glrG$w*FPSjmA%|QysyoLSAy<{oOqGt% z%Wu2N!SwkC^U@2=7!x)_A}L;(mBLKIsomYd%hD9~M}Fw045#yJ$tG2nK(0XfH$3ou~3AVe_^U*0q(e6(($sFjhLtJlX2{ z{;h<)V!}?Ac)E;^vGG3NN?4vT_7fB7&9L$I_X(RINtm!ci8m=XuUB-sEKrtPa9&nf z?#Rwyd2VUuyg|7mft2OugJM!{c3F;F$1;9e6;U;C;FZGR2euv-;Gs|u9L+ehfZvZ(jwv`3G+|AHh!S@L zPB}h8{qU6wQW&|)bC+iXuE=}d4!oA?|G(a}=5QQOerOQNW_ZWu{(q;<|NN%SDGM(Linerli`5MJ7z-*7n(~!e-%7Z9BvD?!hT?RJdh{M9&+bTE_h1|Jfb_Eh8vzU;+f3**^FjAH2Z_pXY<0;)7d{GC4TlG$jM( zDWH1L)S!TKS}F)6;I#SR?LPQIAAGS7eufXe#0SSF@H9ujIVTM&@^|{+-9GsFDLBj@ z>wgJ5?<0zUvy`3^w*g1C>ZBe4rzaHz5^$FJ;JrThg+BOlAN*n;{O3OSB|i8seDGiT z;MR*u4h}fd3#J|cNBMp!0RiXg6ik4CbBzyvtq*>k4}OCWzS0N3DFr_QNQU!I7Z1C4 zBImE@DRCQcew}!ryMXf>dP>{|oLdqPbQf@bM^A~{fOBi&f$jp%?evtm4LI@xA*n~e zxib|65^(PF!B_d)8y`H;FEBwssz90};AHsV z@*@tZN5BdA;C5y&IZnU{`q0ZS5u_dgM_zzQJpzv1A(R{zaJKWIALN7Q_~3(m@F709 z{61Ie5peQ+@Ev^cp+5L9AAA=deAg6Qeudd_QvJ;56#j6E@OmzDX$tA}x7xcF0ibSr!nL+ZZA&`hn|e{uT?I-+dRNg8NC5T@or2>Gv`_AJXN9xZrh|$U?v=)I zj#jwW8XQO7cygSbUB324zG3X7Xj6v<9J?B1sV(lQEG9Sr&oSMBekvUdG|LC~wqlwr zAU|!D=0HD{1{HrYRc}Iwe#(bFRd1SJX4F!TKuT+z}YDk1QKw(t(E1IYF82h`_TlR{{hF_`e=fHv#W>F z;y>`gck{t__rblbhULGf4}F0TzLyU^+y~#=2QT!&NBH1B^ufKYo+b;ZmCV#A0cW37 z5JH0JDk>8Y)`yF;R!Fz_PAdw{5Y47#q-fs z4|k4F!AlgL@an8*Q0s%o6n>)1$F>wV0Z-T4-4A>(=R=p?`k}A*(69HwKLbA8NljP2 z#WLJE+2wEL`4RBF+}$&d+masaLqAR7vr_mhP<%dcaVvkf5B)U?pY77y9`pg=>Gbn~ z5B{|epIzZX4R_|a{70C(5`{OoxJ`Ga557?0^IUq%=NAe;#l>fsyE}mI<$UAf*0P@! zJacmj{vq&m`P~7oUOIdP@ZrvUmyb=i4ESEo7cOqoJyGz?PgC$S6rW}nx1PlHz|*Ds zTSb4GOK)kO^`U=5(YK}0|5M@ZE^gBuk)1x>9|IrmEOhCAX6~AN=+E}S&-1}w^1(mx z!S{npJ=|HGBF{0vtxnuP+{{+=XSk4+f2qQkxVY_!uT}WjDfnuIpOb?BN#UI?ZY};F z3hz$AcONA3cRMX6#I^_}3hzw8rz?EEGt1E1mSwTRPf5YA2A)n2Z~NdIeefWLBEvId zuKdR?`#QdAZ7brY2|7P(^fv1!6*FN}aAN&~~{4F2+03ZA# z;KQ9>H@~vaOOAnqC0~hoQ7ix63Qx?7TKr((>GBo!!K)RY#Js5Gvqa&Ec~Ohs=)>m@ z;OXRlM&ZkKX7Int+-*N3y_^e%NjvtOOK;ofeSPT1`QS$aAMRY@@|j}p@bWkC%UxXd zt;zA651&tb@PQcO4tK8Ak@4U1{|WGP`k&~7AL)ai=!4h!;C~vLKHUYNAMRY|$~nv2 zUEss#Vjuh}#b>30G zKJ+L1;H|)i(|Y$ohR>}&^!+~cBVi=NX&qel0?IK5csl!8=!4(wgMZ_LAA~8Dbo{6I z;P@fXwDfCy@PU|`8BXi<OE-w28<#+~o zI{E+Pga5?`-{gaD@xik(wUjR1?S1fJKKLF!_;4S5jZeMa--mv*555TaUe5DciTH2Z z)!S5i`)3!od0g(p=l4GNW594RsW9MjS7;DZ7`YzymIUl>E9p>(Og>Op1`xO4Hi=!D-$ICwar=wib$-hPMS?BUWwNl5< zyNEn5x;V(xFFi#2lsR5o&{su+=j!SR(oTKn}UEHR7mBRnx;y*BVKa=|6 zyzk;Rt^0h^{f7@;1co-Q%M3T?`{3{T;6K@ijGdixw3Z7mwnfDygiijusc!&8uy#iMG}@M%yauDynPouVV73 zimD^)rXE>6y5@k!w$|DM#!djENOAFGd}g_#Dcagt-v;3D8OM*Bz+Z>*ehR8$f_zr0 zrfu;8P=6m<@9WC7ZIwq=*NmG{Ul(m|t&B7lkEy_O@z|=?`nLMk)9b4NWZ0FQ!6G*9G$h|3CJWqpMVcoHPyE^*40F7=T_*vjVsYvDNf3Y{?Ju&2EL2p z@W|1V$|5y2ZSAvbDr@8g@QOt(&8<=Vz)VXE5~+?xCe+k4Ha9jkM)9aK87Z!*NIZ<4 zTG!kjO$P8fc!hI7WOA|m8$DkBk5ael-vsr)ME;La6h(47R{o7sp!`#jid9IlB2x6l zqb00ZF)JP;_ljAuid?K>7wecBSCK|3E~6BeQHsGR6?K$~HA=-ArG-@Vqcy!sWRwgd=yFQQGchEO^Q_1AaB!VnD_EMLQ4IZGrvleeRVXU?Q)hVB_$RkC5}xh5s;{a39OuX zuwG5dmJ0_Z_%Q0>lwxfgMV(HTY|+I#O51=CKP~8>96)bNE`v{Q4$H zijo*tQ)9l~QPbFjUOSq$e5QhOv@SyXtq(K^MKu-G<7%qJA(=1vu5K==nPS-><%$Xj zGEJHcG_&7f0b}dt*VndIAe~V)l`sRG(b}l5I#FH|noqB9ZJ6J@u%gjH?;gdsFl*}O zom$fX?QEqKlA@WyrKW@{LQ>0+6yRz3O!b*|dQnZ&s4+FEwL$DwGw`R8s>X@1-ES>y zTX+kNsW-kdN*bn>piXfceN!p;I(&l$CN_1!0%t+}g1VN)+QLd|YNGR6n-|v1Z=MSm zPklxQX6hDEqK_I?zo;>4KXi(kRetONt@RBD)YQyfw5X<~zO}6xAIZT_0oI%@ehBaO2XDN#W^Q|J>zta}y3^VlTkBD`53i3-7XP{y9{_`yB}iK9 zVX$qs3PU-wtOQOyAs8XVBv@K$q zNTW~^@)@OSX?PB6ZC;?8SZTzj&Yp^hbydh+QB9)htw76%#s}UWIfBKbQ5Tdah^DYn z`+T~g)zx5j?VOsr+O}xIM^4m$B8;41Ge@*9XsK|EO^Ug(X&weh(KcA*RBd=i&pt)U zeQFIl9`p|F_1fy3X~$Jn9DSI?K^s3eI?uH~HXgyRx5To8!&eKCXh&n9}*fR6>NU<>{<+8(UCUP$cMHh zOY>{n+Nu{XnBA;MVb5hz@R*GtQ7H{>fwb1oZN&HZTBl-&*4Tt0@YI@ixDbn*Y8N!t zX(lo(X1h%IZph}=J$>u({cF)-+5c+7rYTNir`+L-ai%j&WtG)Y*P}v{Blb3J8t21w z$<4aBwehFx{`K^O$;)^2<71S~SHa*DN_HJ&J*bB3Ms38uB!uy(C%gUCcStq9uDNA# z!-6QDO+C)2m^|(1sp=6fikN|E!S$4+c1^^1$lp2L>yr3|f48oA z#vw(`^sRjFY$#Ng#*-YK^wIr_gq&WN(~q!toHrp=RT|rlYM1`6#?FzNZd`kov(c3` zpw9OLgA6njjv2j97PvhjuQCy^lrQ%WPDQK!>aL>eM>zeCZ>*^CWEHiea zWQHB8bQ=WEkWP9#*u{ABa|^CB<{GSD2$>p;S^X@-@#FNJww-O6FGq&?&rKadp>cP6 z^)Z$3h7z4%GOoLrnmKLFHS=nl=D5Qb*bIhQl656iZjXH6gr{ zX-X8gr$`@~-or^*x=G4avZAdj(P*=)tClcD~pS3Y8K2!SBBwgRHm)|b>1ay^O_gp zcU2e5Y~uIGDrO7Ou3*CQ`)Q~Iewdj!`(DY7ZC8_R-z!c@%WO=`$&8##(#h~cFJGiF z<}^Xd9%b=#j9Ed=!u3$>Zos$a)o~u-M-zTM;qq+@OMe^T<%HjDaM{*+Ax=wwKhYmU z^v@9eVxqS@2t>|vh(3gSIfTE2S^oA7buYfBp%>a2#ODBm?+d(!@DhUypW6sOobcNT z7hTCA{O=&Vmhf4GHxm8_!p|c7PQouDypHf22yY?$HwG8^rx7mSDUn0)>4dK_xX90Z zJ|~>z+(8)N@Y4U-2iNa~>vWeRo-MC1(f^d>(eH()mhTKhFBeH$PC4{@;re+d@i~>~ zSr09gF582CpIfJU5%^j;&m%q;6RzLq*6AKg>0V6qobKg>v;MCkd075G5Ix)BYNBWP z*Ao5DNzUg8XZc?wJ}m!hM9=BIO*q^C`^1Ok-$XdK4_^`==ARjW&^V-Ak0X8B??`&> z?I5DRi|8u}|2^Ty5YF}SIKo#E{VaovxNL9o8=P_oCf^jacG5yP^IvFiX=kn?`g4hX z3E|5K=lXt?!L2+t-IYYodb@@2Cm=|2v7k&d^J~9w9tR_@jiMPyAUw*AmX@ z-bFa4+eiF4-RBACbl)XDTu<@q_Ue#)ak_gG&h&@);L`|aKBo}Qd@l6CZ!ox5pEncz z4ASS_M1LvrV9W7A!nvP%l=yHzW%mbr(|w-!us&ZVKHN{eOY~eW8;PFf|B~og{(lnA z@(0X5VJR<`e~`hw@^4Q#%Rh|xu>8Y`p5^~B(Q`dLm~gh|BMIknuO^(^!{ZI^m47zT zbGoMy&iv2z!Oth0)BW1ulCJ2~yiVf`FdMhMc#gq^mdkNhgL~Vh{fN(PWVa=Rv;T7s z@#lK<9O2v!yh?nyygneD%Y6sD-zkT;zW<1D=0Dcp-f}sNaOQs+@nQMTBb@nPLVQ^M z8wqFr>xd8Y{|n*F|1;vl{IdpWJxIQo|4$6=)lWI$%zrBJVgB~LBQO68h!4wu1L4g7 zb|3za5YF;Dcpp>_k^gbBw;c)R{&7#jS)TEPvpkiA*Af5Y3FrFQ=5aH^SKt#}dx+mlMwP zClSu{QNo$NmvE+ENjTF#KseLCL^#ubL^#uD;XPM5yydk!;Y@#+58mK|pXY<$_m2bdKC2u;&+{R(X;%A6Ftv|97Q;=}U)f#{jfYNBWP*AhL;{~X~g|BJ+j<$sOnIo-DjXZhbJ zJ}m!dM9+LS6Ftj6(7Y@q_RsR?63+7PXz&#KCwfkIFTz>=9}*vye>BlEp9w_I@*ht0 zEdNo2v;5V>hxLCv(Q~@BgtPn&#E0c?BYNhui0E1V^NF72zkqO-|L4Sq<-daHIo<0B zXZdd;J}mzqh@SbZCVG~CEzz_5&k@e@zes#o{?~|}(|wz8mj8X?!}4z;JZV0b_;CA{ ziTD2H5PRlv+%Uqqznerj&qp0iIM2r(OE~+tHH5FCcD}*jUjL+t=y^W&bi$ec5`*uH zeDQqjQljVi*oz6TBYAEoob~XC5B{PL{+CWA|P@qF@?M9=#C1L2&nkA3j331>SA?5NbE?3U*ZLWHw@%J&=O5IK2XdLN=+ zMCD#gINQ}>g!4S}k%Y6H^9X191%xwwJK>zKa|!2sT|hX`OaFp!mj8EzGoL#MXL;@> z{BF|E8p4^+bA&UWmkH;2>30ln>$@GNd_eTvE`3Tk>t~1=$VfZJcD0Mag+KeHdk{VI z-<$B&l;3@c54X=_iJti!K{(4Z!-r2T(KDX~KJ<%x=z9t0bT9Yevy$jJ-BpA$ANf9q z9HLLw!?Q%seBL3P`HbC3u}R@mYH)9V^p+-9_5EBb7qhIS)*W^uynyJr-t9y99z@T4 z_9Fb>;Zyh@;`1ZoW7MhiAm4ts^Xe==)7$S3h(4KqFVX|+$B)ndWV%10eC=W6lzyD^ z^%Ihx^EJ-kVxL9eW92S2xbP_^{-+Yo>7GM4rz<|59A151rD0{yobF8q_oh45;KFAV z$+OxA=W-cC^#6{2#)8lH(GREl@92lq{ZHvf^eBgv`!8_XdSUY?;p0hfX4y@ZW63t) z+rs}K;&T+zlSAV1{;MK|YWxtQuOqya@Y4t{Bm6wVClPM>i@l8_yu#3X^%EgJtRLpj z;|S(I@q6-T{rK@`JKvUecqGZs?Rhcbe(l?~w9o$z|841yzAt|>9j^L2u8(_}{CeBT zv7`rXCl4~X*aO$c!-x;JlgAPtZm-&ip5?iKaBe^C7Z#;nusk;zdapdU5g(T4VS@`F zmS+~>T<%Z#;H*z>|Nk9*a{K@9=#%C7cl62f{I~UK>#vmKG-@Brbh~PwrxT9P0I9=E zZ}l(bJA>%=Cwk^Hn(&KoZ`=O^3BQ=|Lk%wcc^<}=msg%C#D~-6ew^oFP9S>L1Gn?+ zCtE*6qpMh?M`Ci>BYR}n68<&b=Fy7P#h?Pnp;UxqkV4`&h1cG5|F z*iL>%^qlTRgmZp>L44RwZXkLtubYXU<-eQgFDLmQB%I}cl=!gxPZ2$*`#j-n4=)oR zmVX1`>?iz{_*_VGen|Ak;NI&0OTzacd<)@RFZAzaX*&@cuzdCtL%<>S!|neVgUf~U zdl1oceLRBbxqPSk(3|B2m2Lycv$qI?L*$u5xYR2-q=Zx}Ouy@_$2oSpIC56@)YY>xmD``Fo;YL3~yb z&U_y9;qw&HUq^giAe{NU=EDb3)giPp7nCnII3)jaSbRKmA%|h(;xcxXW3axs_+Y~) z&x21k^mabM@;}_*cAm}RM;e^{TaiB_8|@gJ0&s zHyHdX5B{FPS9tJ^24CsHKQ#Ew9^Afnbc+Z7+|b|V!8aTHP7nT#!B=~5#|&ie_uyFu z@AKf<27k

tG{9^Ag~^o$1|YUtN_@O*>6^56{ypX0%AG5lrB#r$tGxcK&r-)Zm` z55C&qQ4fB!5=ethX;Sc;9Va48G|qN;Oh*&%!9vV@Z}!-6@y>m!QU|W zWgdKk!LRb*?-_iB2j6J$l^*;egWv4IHyQjE5B|BqZ}Z@r4SuHw|Hj~}J-Brw@Au$Y zS-QOYJb1RjAM@aY4gQ1&&mE-sJmbNK8v1n}Tz2CDk+Z*~1J@^9ujqq}!j}p%EFuk3hmw6AK zPv&{-g~Vq#@jspL6MgWV2wzO}yBXY@?*2qSi|DPKQVciXwECPt^gku~$%OO1#A?Es zKg-W@&iA2jGq}iqCZ&58(Jvv~t`mtpaJoGIzY<|KUork@)+vVL+H;ZK4%(SpsIP3FN;@?AjT8W;^YZcM6JWmlWK9!YclIfpCelFh_;iW`>HQ`LZ(%|0w{?-Ry zLpaO-sSp2u5IyTBm--*Jx19_w@^Cp8_~4TXXM5oH(|H~L6q+CICH-_0E_KuD`DVha z2!EPzrvI35re9(F3n?$g@Abj|O!$0C*S_~E`WZ@i<{$t#gg?^{Gq{)j1fpm8-yodl zY5ztz+s{7;=k_hfY=DsJ$Metjy;Q;Fw0i!5p%?i%-7$o7x(5-?=^p8mZZ+YY?g_+) z)15=~7m}QfM1KtFpWk~eCwv*vv;3D3egVk3`S)%Z@Kae(o=-1;8N}Za-TH=XU-Q!kPYd!kPXY>Q5QJ z+y}pfaBhF@CY;-uhY06(<7vXVe!XpQuf2Ui^sKjk63+bX_+GNj<-V)w$E6M6`n#vW zg`Vs0em?X?KJS8emJKx6CO{zU&%;y;G)pA$aO;37YduWkDy`sDH8G@|Em zjrW^52Nl$}jJT$szR3#AhktJgyhnll)xWWg2&(yS|e&14XmWSsnnEy`Xw=nc&2xmS=_~6!06FJ*SpUp(i_~`~0JK=u8 zuH$+6T<^o@Uc$M)*!3%~K412s|G)>|TNvPwV&-;rKZDDK+l>Q>{&xtoc6BJxv)%Id zA$eZc&i6|+t{3+Ekly;mekQkv{C!B~!{32q|L1Cw^J*%`8wtOL@D~U_l<)%MaEP3R zgoh0-(b(_3+2Eq*-{Z7;b+xqR*ZPm$+3O7|L~=X7r&oZEq? z2sog*6#mD>5ALfJ4^ug`EM=$>~eCRs~Kc4in8_mD*cxrEh zOS#BiAZxcHiT)zuW51i{m1ipPxt;hNOMI>;J{^R!eO~E<|K0~*P52(f|0%+6AiUq; zl3$i{6VbDr0|tQ*hc~}rg9{%X-%TNU9^V~LIFExOio52|?N67%h5wBt|Mi65MEG6A zhxPfG!S}^=C!&uNJ*WFa(8wY5HAH`u!9@?hB3x+X5S-g*$)g;C-%RvU2jvj_*Mu)I zxa9W^oYoG{C3@y_G2zTdd~!L25BGP{=F8!g=l6s&pSy?;``h;tJdn@mtd8ABp~K z!e1u(vj}HDnd#Y2W;@CC@#_i*XZ;r$T$1PaamomvP4Y}3KAi4}gmZfm^}(Md{@fnE zK=j-m_7l$iuiYOccFXny_geYZ~Hb1d=s z8PRjP7ZJ|(^9$lLlla^~^xS`0KUuQP^4ob`(ZB4Gu<}1hd|3WRDP5NTDWd0e`McuW zZ@x@?SpIj358KH`qG$QPBzpPH7c0M=H9{Sd03w(5TE0y+|MDL+mlNP=X7r%{%mi*C3@!1-%DrxJmf?F z8sRLD{cgGBm+?IE{~7;p`)O}L4mm`gLrI^z;9d@I`&MRfq33>iD&gE-ok%$Mn~j8X zf7fDgu{ZX|WDZ0QNtgL_6V7~=5g+!)WK1E4m(O*CGapGy4w2_ZoVI>RTP=r|&jW-r zpXUkZ@_mX z8=JJH#H6rY7MPhP((T#-T$BYa&vF;rm5%bJ?-2#-}}G+`|tnmyZ64C_a?j^<8>OZ$6inN zj4z7d%L(Urt|Xk}=~f&Sug7}Ho{uB1(>QN-pAu-i9=k5$=lz7UpRJ0c;`P`^$e!c5 zi*WXHkK+D%>|wI!_0`wNp4VefMeKh@IQNqm(Yl54;|b^dpRPFE@%Yw6IM1gomvArzQK#NI}KxF}wfTAFs|Kd(O`a#gQK# zy*JE#w_<0ZE9RGI0IsVgWzUJdn7n$D|l0C1lE+m}i_r>Ig=l5k~&$@3VocqJ& z^YzBiP(RFaE|9o5qxh1FH8GT5a<0mN8`PHpAOtUs{wv3UZenjfb7QrN&&u;@FoCK zP_h4J!r8yg*VxVeTgaaMFNSj|7@qz?bZtGeEJ59RZ!_U`zX8V8M#8_L_VDv=!k-{N zHxYiUY@@)>nS{Sqaroi;qs}7x8_0fu>{<7NWPbtKqpzYs_Zq_YDh}Ou5zg}z_s<{e z`~^Q856?$z&+`$tm*=B8(!B(MNb&Q0wc@Dek^L^hZz23gg!6oqp>crcqcaF+d%iE4 z?d^IHGK`-^=L_SD2_%9`##~;4#yOHc|eDL#jGx_1X-AXw7|0Lmj9Cs*= z(}~B`e%lobMOV3u4nMCZd>`Qt5&mPv;is2yd#)LH zAK?edo{#%W#0CZSuOfTg6OID>V}#>93i<=KAJcL6l>nsx`wH35PEePhuM+;3ii5t3 z@YfQ~{?Ak#*w@IuUUB5@WWpCHZgJZ03@jp?buU)j_uoeL*!PHNpm0Crd1RH_o#bbM z*joE@gu{(bDP@Rd1QXD215Cvt15Z3 z_#m?a$A7S9=q z1LOXGwrAfX`%$n~L4iH@&kM-DLAI@bE+HJ8W_^!25yjdy*S$?XLuTB@@gC29vA%zs zhnpVmDPlh>^!+ze)wI*(XMaEwzcbN^27aVKiSul z{jbO#8WztVJ>24XIT|X(+WUywouD}4;dOF0Vt(B$ zRDK>mOn!FP=;sNsx45}o-y%QO20*bqSl@=N`GftpYxMISvY&~cJrO^4t|3kml0MJQ}hwmf2fp8mp z0W}l;L$Yrn{CvXO2>%h`9fUta_-ev`On9E~pAggmM|L4xLtz{tj{~^!kQXV+ZwpK!%oKGHqnz2&~Tri%1ienvadFSJW>k-4* zKkycHVEf@8ehg>(|HMu$Y>)M=m8@<6vNau`X2RJ|3*ncLA3k2~gmXNc2lmgnx=-?k zb)dz~`N4YFaP|+};80rg5Z;#mX{AC=N#@0W7Lw>SKc(K&_ z(OrTx5w0IKC5Yt*_Iis!g4=o@xa9;phLXK*eXF;oCddthAEn!)oP1suZLUzI7m4B7;eA&f!)Y+UjC!-{ed_xd4N(@=;!LNSCbg*loWr1^;g(`vJ0VB;4-ngZ+HMcay!H``B&2zYY6y$^K*=EP&fH_1K*! z9Q(!rO1YMB`)zaJTM2&?04=u@zQ{nSUm*Mf!sqB91Pd9T^@FIuWbz%`wR8u;%GTj-%}_TGWBalMluOdpSxn^@~-lN z?((LgVtaRMb5AL^v9G5%*xF47O`Y9CL+2GoyNg5RzS3aR#@<4?b%Srv*ihV1=+700 zhfBjTlcvE^IoDkp7%G(e*7O%AnI%fieS_uV@L-|8@!XY#fnrZbf1$hBR_gC54lkJp z-Aa29enM&V9LMacjfh!rE8;etlG$>)xa46QE@wstJLa2D`(t zH%@df;n!r$jqqY=Y$*2)ZdzL!9w?O0?J1TEef{UXC3$!8@VdqZx1)VvsJ~q@(Av78 zv#ZqAl#_pb9m~UgbYiHI!sgP3;_%x3(#F<4Puz9ti3DleGz9A1ksdLjcgJk}b17oy zcbgPuI8&UN>s?8=X_mH6>}=HmU|w;wuiUwEX;*h?sMxw;ihdU5Wa`&-cjw&Tm6L(7 zzgX$iQ*gY9s+nCp4?TX(6p>7_og#vO?c~9&X{1~lYVBOUva7r>m+KxKZES3ux27=C z*PW9D_6@FU@9ru$Hg@+GhI8fNLSK1AM)IzLytKAXp*=IV)dWp>qMTPQwGWo7rwL9D zmpdGGo0o7zYIYND+99}^vMKW#-4NPUE|iOJ87>XDF`!jW_s%t)(wy!NIWP^m>e$fQ zH9yx^ZDXZt2O^u1bAPx!eE#*w{5kTH`k)j#R<)(;G|b}Z~_UER^y zx_FB7xH>unnO^suxN@04SnjQgV05=`(0HCZXLIF}oUEJNa_P{{?O^PO(SPt*9ajGn z1^^n}7VzkXMOJix(}kC1H#CayT#eHRJMEkArx6S_}w0w!}>J6f?J3F*>UL~-ltnL2-6v7GBKjK~#B zp}(INAJ*F%gBISz<{JewG})*U#{9<0g%`$*?$TgSAFjf3<)SP;Wm0N)mq68}L~C{A z!66B+Rc7us*pS2FmJjRn{<_N4^2)fLH2e@FoYJEnmt^H?GqAtF^iJ5|i6-5%YhLaU zPd8!O=~n3{t{2dx<;Ku%>A&qO`Wl-=s*)X9p;Q+rJm5{7#Tpx?_QCBODHn#z>xUHT z8|*9R3d6&NO#;~IvZ`xRPgV2eQd17Q*ve9H)#M%>@WPievdm37P1uib>Br-LcrcAO z@A3%NzQ@W@nAaeyU-?92{ZP5s(`p}}Of8Ut?CVXw3Z7abY(I@xxp=F|tJ?csz01a_ zT?UM&@;7aW>@N~uq(z0L=J*^ zpJ*RX@?sr>=wm1O(4w^=*HsFipsT`!u(3wNhY;Sw#(_`}cTJ(h8yW~DY@B@5Dnr`R z%dc37LpM|&t_;Gh^VhB)G=ph&3%^k=L<`+-!xesaQD%?6iTtT|3%|LiI8Yc|*Izue zi;YC>%?pxGA1Xor8GK5l^BRWb^XA-$tiHI(cF8ue-L)83ePd-uzBdGYPH7~*wDVRz zt@b{wn&`6@_5Al_G9v zJ?Ed?lNh{Ko7ywhK);5QuTqR#t^Y`lJvqd-b>z(9TklFIRt$s~esbl=M03Y|~rR8@7 z__+TrEq`@@kNfM=^78@y{eq?C_bT7+LBRbDcvc`S|Jnfmfg1Q*0({)VkyihX0RO=n z`1b|)U#NlqsPbQ-HS&l1jnc+%_oXl&_j9G?KNHZ$eQs&__I)tx<6f4u{1*cH=hwiu z?`>Hh_d}-D$GMKed^~HBmY)sqaX(;MzTMl(`nYd1E#K}NVm|KMOUqvph#&X0rR8@h zf0j#QGJbZ}z|RNxUldrNe|Ho@Mes*s(>*JZAwE2HN5I>GdTE0Ce!us#8fj?XC17ZFL zYT(;_bIixRdTHag=WCe%=^FTU9}D+CJl~pD-|p{Y{ugTCU#a=y_1B{{@OzZc`nXRu zZT$B8i>#05bJFrh1Nz%);NK8vKhCMN`Zovk@tj6lzCF*%?Z@*#Y5Cg%`nZ2ME#IDN zWqmwfl$L)_K>x`a_~U{6f2Riiu7EzCuSgsJlL7rbHSp~|SML9~KR>Ph{(wH7gG$Rk z5NQA28u%{+^zr;uTKyyRepqfl?rTrWpB3Q0rUw400sh%3GfRw@OK9IYii)z-+kcr-&q6S?z3n9T{ZBZ z48*^x2L9dv|FRnR&jk28YT!Q`;D4qD{`1P`@#BFS_}0;wzq1Cu{T+1Xf4&C3{oWDt z@2`PBS5Fk?;)B3Yrm%dIg=GQ{w|Afxt#>v`-)_67~X}@P@ zeBZC%{>})0i|Xs+8WmNaKN8_@OUeIWguf#tAD=0P+CQ#*?Bn9ce?WGI_`8*V4wGfu z=l@)=r2lyN-=aUT(x!a$pDShC?>}<{59!;(UeoCxFXWK^o`C+<0eyTX7}9Uhhs~zb z$Fs5_{p>GY!@2$WfIjBrkbYkE?L2edCg=aJBl^9n-=IX~e@#Gtb3}hk_0KS2t^MPF zeMEn2p#41o{m(}9ccg5;&8aP3|JxU6|Js25-iZG0l|mrT||Er1B(>a@0a(! z{-vrvFR`w^T&oLrfBeBpfb&{uf8|>wfK6qe@2pH_Eo(eq5)8+CQcpZ#w<=MfA6+exCF< z1@ylh(VvBfkfltgKOWJ~{KoO~u1fmZ$>@%dR{Aja2bWp+mRmn#1_=ce*%Jgf0{W3cVd-$TmBoWaF!|1To#&;H)Axc%=I z$l|wo;S`x5q5onT$GQAK^?g0qa#^W-KYzbc{%hG-wvBH!o{Pt}-~Ox;S9@{!@!;>} zND4Ck`f2E(6Sm(sK>xc%Aiw>$McO|)*YW-JyT$GH%JbZj_;)D3!wj|dCVa zJig!l&j~+t{KmEYJb&IQzD(cx?{6aQe^B*(J&W7II$xgu54C@{^7}o4=^Ol^2!GL; zPT^P+(%N`d|3UeF{yLQZ3lmhn-=61%AIe|z|2RSic1l72`LNXH&hpnO0Alg`PfG(> z@|zvr`KRvOPI{`tQ{^|L0dwewodW9BHg?dNZ|@^O97#gG5LBk^z5__t~QL;mg% zM}GXPRNwEvUAl16X`~W@7f$B>Nc>~UUsPeb4L|;SgdaNo&Br;Vyz&wMCuQ4@{{_|e z^FJQZH@>xFh5Seqw*CCukJa4kO#puUONAebf2+oSK*m=T#Q!O&P2a|k_o}|1|9jt@ z&@zALs>-XyZNSr!_*)h`izQ~LwIBc9Nc@{mc1qm;?-Cz={AWpiLiyXG`c2w?3&Yo+ zEBsLYGN(BGyj-84Ab)pD?dxw=eLsJXv?sKTP6O_ZygILyx0!>H_zx%_V?Gx@{@+C6 z&!6U$4k#J%e^$2r_)mGIi^-4w-OC&?Hvd^nZOZri|354LDl^pD+HdtW5&l8t`{%#K zZ~XTw-#>q{%bj5NQ4<63^LMlGL&v}4RZd}#l%;I)B=3zyE9yescWGKD&4L4%1JpsxMf^_r|y5`p!uFW2%2<;-&hMt!^e;pNqu5UHN`{P2b#nIuieCjlV$+ z(SM$jZR1<~Kaa#euKM~lqRCDA{^>+HK|=j!U&{O+EBsLYHwW_n-9Y>otG?fVGCFZj z7k_^w{;cvXPnH17pS@ofiGOz>{_h9k|4=0U7S%t=gtYei!xtm*uTB~N=OgjY(I2gD zQ(g3*eX?!sx4ir`68})j_>Vb0+lyT;Gw?@t2p-x|^1t@?E)q_waAzKH%F)z2#r$A7t`~x%OpLWvB+P`Fm{Jt6Tx6Y9NNA$M_^zoID-~PQ3 z{g#yaPet^11oUSI^j{|TnF~DW__w9hKSKDS{O=0r&gSL5rX;PI2)f8k6fypG=O-kyEh1h|qnL*?U>Apbu7BR@UTJ}!ut G%l`qtz`HO2 diff --git a/anRpackage/src/anRpackage.so b/anRpackage/src/anRpackage.so deleted file mode 100755 index 298ce6156a083310db27707f111dca3eed01d83e..0000000000000000000000000000000000000000 GIT binary patch literal 0 HcmV?d00001 literal 75480 zcmd@7e_Yg6_6Lp+3Pq(2m9%VGO@(DY*1$wXqCg;DQ$vxcS*9bxAS438P}oljD934v z+J0DDt+mx!ZJV`aX=V5e?bk;3W9=u>kL@*5X;x~v#^?Dw_r8V~Bh$9;@Ase2<1@v) z?m73|bI(2Z+;i@)*LYol=j=goaf2QDGuUyaBMxcTxeAe^@zakXR(*0D;~o1uavVoH zjuP4gD`SA?W1AJydV;ec$ZbFTZPxjlH!HmL+;V`BIkp_Co>sX*09h{a6NOLw#3c#q zxm}BCJrfm~BOXuwZPfYE&s0rkJxR+?(U<-*rc*iu+~@hV=4U;vau=W+<@~LmBuzM3 z*B47KBllN?)>Bt6Z5xIxcI2e%J#r^_{n@-u(LhaXF5Q5B%iB5KHFauLghA z3sszy8i*f#eNLj|@cdJ5$e*5e?8)h#*@?Fgil24PX-#cs=eb6<4|TuzM(Yqq{@HgJ z`#DOEZu>fK@S^kMp3S@a&|7Ph$*UHB_4sSRA3vAluMvM&;qN;9*+16}&{8^tz zkywhq$ME+zjeH8vr}6g>{5_Aq7x1?Ne=pU1CGcF2Kbpp?`1|C-J&SfO zem`)>5vQcR^Xa0W?p@mW%8ZGZ9Qu#PzPshrVsG1H$tRv){ZH@urLXL_)xG8*d{yIaQnz2eB%er#X9boB)rQ{H*0=9}}pZ7V;&{py+j9v)iw zTK7LDe44WT^3fl@_q&Ti3+r1wM}B!{-LT;Aw*E3L<&v8>Up9UC>6h)w9lh;@=AYsk z?}o}_el9uC;fS~n8&xp|$RRe^N>(2}F#WMZ2gcum0Xq=?v%q{H{3{0yjGu==I}rVZ z81me5@WAx%!Y>Tet_5iCK>RO=A%9K`{&Ql`FMys0s&_IJdz__Q+j9dP5uF1_w3}pRLftK=`H3f$@*R2@XVmcntcFV#v8IMt_B4@b|^gTWJjaZ;jE9Eiu}g zhd?<{zx*vmyJ}WV>sQw;iTQ1U?Xq{ZNKTa5nwCWbzTV-hk@KRytny#+Dk zJQw;ONY2M%jMK|v@Hr_)yY7l1|E?Hzo)$x%`7y>(NsM}5Paar5voL`eNN*Ry&Ii)N zpcr<(C`P+h#lUZZKOBhvf*AIHPmKQkXAC`@6r>Dd>HnbdyD|SHKIwLau+}=CV!RW-q*CEiH2qr5XL}T%U9c?%dLNYmGV59M)0@xNB&JlIjYUJ)Ll=_jHcq)++KC76*K|G4JAJ&OOk7>`o#!-{Z! zj0b*J;z@i_t%9dRZhp>&{Se=Dt%A=`4@U#+jrh)a3Sr5Aj^^*WU*WSf|6J&g^hq@e zU#9t7CH4kA{6irQ)OsEf!ye$yF|LO?a#m>u`teE(d$=oxJ=~}Jy*5hEA8PvO{@SbU zA?F!YVG8K@p*e|uw7pq&4sytOfFtp2QuTdSqXF!%qK{pIRNoGVPwS^i>t~iOlB)H* zMC+M+`SCz+9#!w@+McIJ*~6h)4@;u(ZK8*R9i8hGUy#Y?dTl>NTK*J`9|46VIZ_-ieCE;^ z^|)U5<1*ckmfb$8^(Gr*Ahr5st=3Pi?#H7vw_)(B((kt`L@Df*pR=J4!AH?^6O?#| zt~W!GJ08{JtL-C~VAQw6u~)aNJ!)KR)b`V+?Ps=>LR`R8>a9=+%MS0?^(JX~ta0(I zreC7zEkDqS1o<@0S9~lv$BRCH?|7ldr84yQI^8cWZHFgne2&I1)BOcdK1-ob+Tk+Y zUTge*rR!~vs&`imKVJfFl(SaL54Q68gV+!3M9Xue#?OI0O1pvzk+1PfwVd0voRY2gbYfmDz60L#pj({?9V9mFR1mo zOAE%A)p)O{Eb~M_EhrxAEjqXObZ-^gU0jpt zEwtF6=7ItUG6ilDS^uqzGfJy`C3OX;XSCOie#W1_#2-_g7*6ABF89@ySJhlm0Ow1W zTsqHNS~1UC4(HvJO4Yi~Zxg2uzLX>FJ{*j;gYH?m)ItTtZjHjxU8TgNaR^Ki@Im zH@~!YfwHhnuh(BuS9688s%9=+o&+bDDf*0|Z*D_LU75F}^sW;B%ePoucWF< z#Z4PYoeu_54}-AbO6JtmA^uLSzM`_a%sbgv-!Pv9a~m7IwZ6Lg8b($B0`KK%5qs!! z`+1`|80MhL^D3*Rl~gsL_Hu7+U5(#Y>SwZ{+5)TrYb%_t3R^j2>ispf1(VN(o&T?z zE6^=Mo8jmEK5rzvL+B~)=xAk;i%NCT@;YDtV(0y;!^rg4xU2nr?!s~s0~8vau6%hh z$Ku&_HS?9%;kceWXEF+w7NN~)-ac!Kc-u2PZ4rIDgKi(sy3?VmQnO zpE3-`l#7ZA&Oe)F;P&SFD;s zkkmJ&6}m0|W%e?L9aK;RgYQ$ao)?xDr+di+(<577{GWu-)5N@p3;!=tjVrCGT~I#X z53f}`TF2=B)pybB!z?T3R+Gz^XfFR-I|c|oefs*ZMyVrbKikmrqH{3sv#a+vhBSKk z{MyOqstJsEC2{j^^1ECz&zHRP)93p97`*jZU1224sO3;j$FQ4-kj2$hJqD))7%cV5 z&Tgni*3{}bH4W8ek92(9ldI4O?6o%7qwo$4-;90n|(uHpND3z!1a}TC2CFs&zd@d^XKJZ z2_t?(g|d#|Xe zt70?L$3`Q?J}T?aYoOEj79qU(>Mr-`*rjGSbFhg}4P)|`B38+K2mPgCw6dleOB04A z$b}A>kKF^mPh#hQ=8t~H2>(A00+fsz1QB~kPg8S4MsByKGJOoHh(u^i1K1&@HNLu1 zUtb4tV`DUzxsX0$cu})ySX5s$^I~jGDQN)tSd~%}qOppls`k8MZ+eEM$YPiRq99~} z==T33qJbCOK_m{0al<|Lxy}A~<&>lgpyNjF_Rv|Y9gKa}$=a0zSj0*W$AX{imY5Mq{7+>vYv*9{&z%o$hHy@+R#7&av}}EVEU{!b>qfmia-TP)U{b+kuHs=Fj;hKz zlSkFp@LBIKD?Q~D0HekdoLh=5OhV~|N)_JgD=YDrAa2cp_d>^Ug^uv$n3R|2J$+P~ zW0J=`G0!`C)R;&ry)QX>RE9!NpAOC%dmIV)j}X)RkQT=o`V-$Tj}Yl_EF&M5OM-*2@4KlYqRwMOB~vFq|VPl zYll@<;shwsQJ`^m+;$hbV1LK?8g4kE5^7FxTpTH%;P?e+bJjU*$?m(Mv;!RFnm%yd zNr3lrG-&vpr{ z>oBA6O@)endlY`l`3k=~3g0R5 z>pU+h3g4#rjEurt_i56i@YcEM_$a(}epwWa*XO3wqwv;wYEu;6y05S-3cpOBo4yc* zx9+KQM&Ye<=*>}h>wbenpL<$*v+f6^Md7XU=J8Q@>;BE+DExE{FNwl8X?SB4-ntjn z8I9NH)7zu)*7@=7D7v_;{q`z6bw@YX$#7ozai zeVz6w{OEPcZZ}5Zt@{d_qwocqeoGYIx-XK_*yL%q;q}mx+)f+*S{r?v9v2AG4m*B{t>5kV-8Mc~+W2g@@v-Aw+CQ^hciQN4 zZ1i@#VWWS-MnBtzx5g^7x5U6#*zl*>_?R*9wKlvC;ZkaQ41A*vZ{4k7f!#6i3vKvI zZSo}A?7@z2v(b;U(I?sH?fAtu`Uy7r6dS!Azr;qb=Dn(r%SLa(gYT|C@!xb7>6x?oVAp#qq_}K;&cNj|p|;lW4qMDYNKzm(J!;%XWH;D*zlz`{6-u85gWeKhSyUK$=z(j>nXC}x7hId3bx=) z8(vQh1>a_C*MovlU(0NG>n+mE0m59t+0EXSxl4a0Kgc%(mg3HhhH*f2a*#Ys2GeVC2(i z!{Z8LQzIp%z7o4#k|6#hNvW@-ZkHs8KUq@hE2Y~d3F1$Yl=^aYyCgyUsghD(BfDLa zApQ@MQeSD^E=drdCMoqbzS|`U;?pIizH+)Lv6T7;6n*B#D;SOK9Dd&YdB5dL4+An!!Ch;UJRI1>Tr_4I|wth zh8+TLC43Cw?R(k&&j`B+ZxQ%I!p9Qs6nGQi;|RA4{2Jlo2`>|P4Pl1R@DhPn5@zTO zw+Z|VVTR0blfaJ?K9O*(zz-2-=nT&m_+G*cnc*UVZzs%98O{-S5n+bNaGJn35N2o$ zy9BUXhQx4_z$X(PN4Qqt;|Px@JX_!+2xk&568KQU6A0%B zd>~xd^cb&Ny142?;y;O8FmP~mGD`Fx9^euC!9lgi@+ZeW@rp| z3cQIhLt?mH;MWM}5nd+n8o~^L;UxmEB+SqkZWH(!!e~BY|$$(UI^eAI_8HiN{v}L2823 z{D}O6370@>2=FFQWDTO*$tZZ#xyS%?2OvD|mq&c65h_YSEji%w5a~jBNk$;~d?o`4 zMR=H>OIs82g>OcN5M^Hrzkr6Fe;4SQ+fgJTm9&#!D5!il(B_e7W4W`2RU|x*hN6+r zLJ$0)k$19WAmK*EKsWNw%T2>*Tv@Ds)K*Hg8mpE#> zc9D`SEM}H*E^^H;#d5NLfnSG|<5(b+yb%QE>ndkBt?2ytyD(y z1RoLmTF3m3gjIY%k9^{>0S{>B7Q)4TqBQu z_H_g4lzpAXWFR3853^Zy=@-H`BUy;DkCt9Y-bSNSy%0~@uc)?A_Gaq(#a(^&bs-uG z`=V|0gGM$eqbcVj(#Y$W-;uC@56H+T9yj0tVVcSGE-EpYaF%Lh1yH(?e?StnW7Xp< zJZ$!52q^aTqq47^w6DjM#y+I+hmwz11G7c-+?%9REqhIpw68?fyo*F2nvIu56;(cv z=2tr*9_*`2^>wFIM#8augfFt{Rhn=)2>a~oRhkR!>lkTf_I4#%K+G)T3FMlsVma9l zq&TR>yOn)afWVxma;}C+m}ilFD8Uts*jFYf`9Yy0l+hf-M@03h%_(UH4?CUeqDf_yW$v{FI9_I0? zOIw6*MuQM#ull_Mw8^*j4a+jp-UORM<>RUA89(;f*RDn2w;|yN;Q2u#-(!~t5<2*Z zG_s5N9SI-s0U7zk<1;)!`Z-9=u~cF(p^Nq!%H9Z+Ze$2a(B7pUbc8nhS}LH}*Hk70 z3DfW}-AZFu3bW+3%4Snj9F&jtsx9YaNz%TSO1EXN5P@hmy*Edw9tP69=?91h`?^l` zb(&O0!uR_Kk7U)eG~qN5D*HlnE`tP9VHTB?By+bB_|pG-l$}iZu#^}&vv*7RDCaRP zRH0+$TCO7cV3bdn7G#_(Y)<+G<&D+93^oGWjn(E^MqG#S_AmYvFg!^!972Y%_%&sR zVcm{~(~Rcq(;z)0gFHhJweJPqJPLxF2Bt9`39A9eQJf}-<@kQBY7fb0hRZ>r;M}R6 zP)$Yr-RW#uB=8<1u*)6z!EF2<>TKTaGy?05 zKz!;I;{n7q%*=Hzk565Q?hL|ygNGS`yu^<9L9lV~DKk!2TpO~n?n;?4Gxw6*nQM)} zztcCdWA=lAy2)*@6b{fPwt52p_5`+>kE2g|o}i{5?FpNw!hE6A(4`KC=|WC%`c5OT zrvP$>k~IK##-8)E%T7moBSL{T;MBnZyD)^ry5x^9sWc% zZ3dL(K~5e-K^v|!f>U8d@A!{l<}crI+$hAsuYeQA=MIeuo{?%~{>z_Wt^fyAI6+tV z45)gO$FRaotCG9eBSBB<)`KZ^^ST_R?%Y{x%iCD%3(#9@JL_$5c``ru zU+>O**V$49o1mcQfLC}L0Eyp{HB)Du#w-tQCt~kV6n`4PETh<=?D@3eKteShW-d_S z4Va9u2avY_13>x@3T>GHE*Q_lDj{Jg`v^SDI5r}bJ(L=`Eh1PfgK`wK? zkTy75uIFSYl#vb@khQ4-d{!t?j!=A#B=vIgnFdnxYn8bdg9B3-ci=tqO5}O6eCRH= z-h+I1^WH(u7XD(5^Ohof10L*b;V+gT74K~QtBi&&bN2S$-c`VY(+2ZI|bI7 zn^6z@t(`1IXoJR)b&Uho@h22Fcfvb{$5EvqKS{{C%-i&H1*;BSl9~iPwa%rKU1kDj zCwRVzg(<~)j#UJ&22JzYB#3Whc7gYG?`8f(GOccW?U=5OMhttBFEmw5z8%?2QS1peI6RpSW_X@+W%+U~|dH2E2menXmQ|Z2l z?ON*rO*tf|RsS7XeM8^HYL+LYqh5N7beIa-)tfvL2v$*;`Ga4NdLKvy6+q4KvKtB# zDV;e~Gk&BW;|<@?#%9tDy$c&5^I;bA*J(6x4DLWdqVhls(1019fGIu~PC&-BC-8ZA zqV9pE=$O{_hE#Xvr_Pphf%IfGIGin+Ordd^C?7rn0LQV)BE+ntne}@Bed7}yH(mJo zhkG)~rdZSNKoi4f3mdEpJel7+TgEd>h%AP0lISoYQt=dn`a4X=%pZiOv*ktU%Z&5k z|B>|>+^zX4${_t3$ly0#Cv^uDB->4G$atEVXLb3_s=QN`e~9G+s{A5VK1k{tS^i07 zHg!NCHn{_zn&$&8u0q?|Q(r6j={k0qhky?Q*6&djs+}O^)|p4Jvx3*aE@Zs+0PJ(z zpO7PL=}RaQ!>;Fj)YNmIs;=h_0)OrwzZk!BFLVcbbEmljJ9DQN2G-|JDGa<(Y+V2K zNZLTK0I~DVcFGyjv&Cq8M^OV+& zC=OG56~T|D_S1e0Kx;X7487cmJWt>k^I1BS(92r{+-W`{tDomYjL^&5fi-VeC=n_2 zGWQlAjuHuLXWS#*9AqOA!2WZXDy9F2CH;p$z(~o|2o*H-KDP`SF+$h06&axrn{0&cq-aLy zIkNDCLbdW}rMAs8rb16aUmJ8a{bJSEkD)n6@MT)P8R}8$oG(M`FwBPhM{&FsKqky_ z2@h3q=PwFl1g~l1v@Ovb1e4)3IQTNFRk&}3HE|-5(0v%eiLip}uy7)zk-5fy5F_YQ zg2o*AWj4Xc+<AoYbD)H$hUb0mjt=9a=l!$ z(C*D7Xpr;N!x~O9g5rfdnOmGKH%XP9&X#&6pq9np8=eIK<#fkE-Od)G%&VAr0a6H2 z`K&41Q3j{~>{lp@{lpt+JPH5A!@Lcs@L+6iNSuQ=X+w)RmJ-qVXoR!nL1a>}6Up-q z0Of52D>NN1YQPmaTiWQkz@V2Lxe?SiUnV`r3;8Bl^si|8w@ClYKJOo0{wXW}0iC~``8NYl zd zf05M7Sw4ySLHP%Xd03b4R`U1GR_$l`1*-g&s{B=?zMSQkYxzrvDbwZms`9U>@{?IU zQ7@|QZ3 ze|Hs|D67C{fl;eK^24+zya~-!B3dCND+_M+2`Q@Lham1%2nim72q*-@n<0osg^-vb zh`EA5z!2ynmFv!`_@YC0w}96@O5F}M+vG~SX&X)OP@fK-I6^lLBM5@~?trK%~ zwoE4Bqa-vG;if*qr-V?=x7et65xUGPN&6a7O8A~n)ia!hrfbwg8nr(->g3SACx)hi$ZKz7kC6z8qT;*yL3E-uoND z>k^GnUW%G^eat2Yx@Lqwhx(;*hwd9ZAhw~zE?o+r;{sTz5kNkV;fel^NPyV^8mc&Rc4=5s=mQuC+JnuCpJD|;X zuqPkbVPM-N5gIfCn`dxdk5h(3^FrtYxy`tMCGey^nB|v7mEUYD-*jbftvi_JpWFIs z<0(ezt2jSsebt}Rbf%->z@9_O+nUZe+~GgKXkHU%WWLt$UHIxxv|d|x`m;o#c$TPKD+gCK#b9sr$>;F%K@XpTEH22z+90 z5otVfTB3IG17^eV>bwNHYjWUYgurVhrXCL7VA~z%%jLmBg+~3hsPjdZ2w1{q5KzO zBk+T-1JS?ygTGe(K3)D{t9-D3dE0pTfabr?czK=JsPQt6jS3}fu0oc06&as4uVViA zFWubXtgfw3V2*!qj5YpQS-i7LU}@uRX7VI7WMJAs>5t_M$l2Gb>WlnCO~bd|zwsv(0P|b{vjH z{A$5DZ10F;kkaFXek2i~`66V;=Jx(GdPg__V{*6TJO>`xkK$p5zy#?}B`zd%;bAU9 zQ+nP8@>M~uH7`TB)sYyIc&0jcg zB%1s!TsJpaZ*Q477Y)RTYx1K9N}_)oX4IE4ol^+uMwco@clRo z+F*v+NI4TOP7I$8AtlaX@z(kb7@%?CqXpR{{BT|lk<56o=*j%b|Eie_Ry}JV zP>$5J-b@0&o-T0g`4HH1fCUTlXO!z%Zs`;I`+HwWgC-bK63uL68aOHvkI$td)-4P6 zvPsE@?xaXvMlfS3MvGZ0Rw}}#h}_t}0ty!wGB%&9*x_{8W}_n+pFd~gc6SLUvx&Oq z$s*1G*={-GyGk?512e>L%l_rQ!1R6$jr8`W%n{DkKP5uf2$zBP+=2hN1G@@g+`%H# z*pUCbB?ChJHsI z2%c?v;6G;zxrNzR)Vvn==boXgvd()jC*xQ)GA)O>0h3VgP6Vv4Jw ztZAjL=x25MGtG|uF>!J9jQ2xYD=c?j3BXn+D=*rHVrrrZC z|I%2TD_hI%Ku(f#IWhyUH~oCFx^d9-^C+ka{zQ?aqmpOO7t$mv4g3ZWEnR!MQJkiGyl3Quq&P=63Aw z2&bYYvc6CqfE;vR_+Y_8?^14Oe*EtfovlA%AfZg0Tqrjp8ftnY2l{E=@?7)AL?jAY zdT|*9x$e^S=pke9MEH$~NuI&wW-c6%>WbW^mr}dgj$GVENr4L$>u5>c0#q7jZU*Ki zW7Lg@P#)aoVjK%##x&i!@&I4uT$LAh2K06e%7%(MspF=f4{+YNA6j|FADk^2@hB2_ z9WzO^U2A+we4=ytVc|#75_f4_R!eFdDs+YxfizTvZi>&z4IJVQ#ly148?C+0mUz}1 z%z>|m(-{lE6Bv@(+vy1oNo{S<4Qy~OoWq?X8niOu!cuttwQ-(6W0LAa);EIQDtN+FTYc`#m<&uQuPC%PH*0tJrTW>0ELoZQ`qZ}bE({WZ{)<~yC_geBsTR1fFG z=?#*~EzreIm;52A@YG6steuI{1SZXd704%0!8U|t#{*Br>h_tq`^*n^; zaPWJraV2*H5ai8&lCps9V;3JHRTaQkbRTT%?F^5l<6E9L$Osh=!dkc4OoM~hes-`x zm+wg56qdWP0qHbb-`^w^^lu#u%8lBF{-ul!tE6e7VTo3^8-Zb|y>IvRwE+TN>Tw~y zFU}XkGt*>}g8N;&`tqLDdBb$xeayo}9v8a0FUa5ZCa`GvEU3L_i^aF9K=Iv-(9e0o zY#ZMqo%fQ?%ST?%qdgB>Wz%&2eL6oC`DP8O>ABt_xkrn(K$FCQ#LD|l=gri4Ut+Qo z&PSWWXW`4dobu|{AF4@5YtogR7wNomop-p-d%%`=T%Ky(&+n_&-M}h%Jt_P#zy5_U z#5=4bGM1c8$~|V&8?uyu+OS?ln$I_0Q_o@O(&pWVHOPc&E;5_VxQ6#Zk-pP3MZ0Tz zHh>@Jcg2D4WpbATvdqob4sZRaA&r|M=z$qC!z;nD-})(Yjti#9sFpY*a-S2l*jYLA zb<|EghSAJ$6G}$rXRL>{7{Rm5NhAv%wY@}UoRL$GcFG>g{ujE+d>Cm@0PB(v!ti%n>FF7sPrGEr1h%px4hzMeI#VLGq| z1!_?M^0PppBDj* znuCZD<(jrR0ovfnI^ClAiY-Q+UUUe~E!}|)R?eHqDGcm2aKS|1`FRl-Pv9FKYQL;b z{vQLT=V{*6+0429mwlZl^QDA~(A@s}7Y6okr>p(j;nK{6*4a===qkNmI(SC-2jDq> zmHrI;lb2R8n}+rL$eTWV|Ijw~Zb{`o=Z5(_XWdwWaeopE`(R#D(~I2>d<+EUBET@Q zd}QaP;F0fa+5Z4cYL+)|0mZ`n)VzY$og4B~D_FpJi`-TV0c4RR1F0&eQka^{6gHt- zQky_aU^)TLkqE#iwWKZq;0e7XGc_Vd>iwWV50-__<%6+eKlzLi4M!S*4{1Cn?J(a) zPc>Z?H==>Q3A5c#`4CDDB#Y_^P<22Aj0YHZ-~$wf(0+IYT#2H8ba{dl7-4lJCx|Jg zTU|^vO7l}2FPfE~x!oE1K$a*&u(*iNX}TKBJfTxpJ6n4pb0Hi*xKtH#BFBH-)g3HI zl*!m-0B|iCf}zgdPSjG-+#UyhXkCT#=HnoN5xaHY!%F^p>u7IHSq^6_|J^EgkRigu zo+~ol0naeB%9CnXP@-q>Y!g!CM^BPBVrK4KjtY@~~)8@fluj|=xJ5+Nfm z)FmN{Nb9vL}4q|(-)xS9;iBShv0Cy3bXN?WO6?W+YHort@EZc$WJz;*2H7( z=3HLjfKhcywd|B6j5 zsQqT`J)TuXArza*CbHemmXD=fQsYaimcY%57y5}3@WUFhR;2*#1m)y)CG>6tH9v^h z#GKYm{*!SD66X=2A^SCDp3-noM;ye{xfs?)W+!5~E{5@Lt`#=AWo!(NeZ7K9NSfh7 zokI_h$4g=XQMNXXSn~PDE?R^dOc0o# zs=U9dJebzd`;f#K!4P79TEjgD!|9HA;fZsqiO0@)mpX5(tiZ?gC3^*mfNb7>PG zSb}OA4lr6ja<-l*_RwQA?}0sBl8QY8wGp@n%j`P`^%=u56rqs{8ydt?AUfhxG!jZB z(}-0-rz8-GXdTDVLM`hcFK5eNiTW*8!GHZzRULRq)dH&!d0-W+g&)?Wtzrqm%ni=g ze`7Ha<`?I})7_yfip=Avm3M%0$~`gQ@%pAewK08@Iaky^CAHDg5C@G^!>}ZqobhYO z3tFO$vzjpPD&Tz!u&6ZxZ-L9PnhU@dZNhgb6OM3M1sTW%J*!4#_%5uG)PCtEY<4vN z4BbD3Mmv|chyV+5uB7#8jXG=Sutz%LeX~q*(#NmRrx}`F*dGt}h)3{c+ff@WvgJWI0kzTb z0#JPmN#{>gQWJ368ffh-P|a?>kQGPv+vp^7Vf#j)Gu(_EUCtAHg4JLY;N28dlgkCG z|2jIni#Ruuc-qn24Vt)xGGhwAi3Ni>q@6vcb*fbs-Uz9rhBlyOe$O$q*V(WupH%*T8cq+dwXyLPsiO+}zM`YtQZB;j9J4 z+u0*{hdx7!(;cUoxthnTxJqZES#GB}O4B?78V_c}X}V$L=Wgyu@?`FEwjPg4!u;y_ z7}@eZ2PQ66f7$NYAkE)13@76z|P zlAZ>BqeX);L+?r2PK|RTfS8^@8N`U|@&wK@4@Hx3GWR*_#CfrrymF3plA6Zx*#;)= zPyuGH`Khy!$f3rCb1-2J3(~P6ufj|fOHJ4--0C@5rJqTn5NI(P(vGAEL|4gAq^uFj zgKAX#Ft4VaVplqyHz{C{H{*pPBUF%tw=NQ;6l#Md6hK`{ypDn-IAb``hcE^3Xp{8> zAC~!=Cv&&+rVm(I2$R;r)-Q|zmM}=e?to}Gpb&fIt+O$cEIjEecXO}P6Sx}N`ZxU@ zyZ|D&4M$47NGe8$(r7ZNMo#RA9}g+GB=U1Ai8pFwjYpMm9_!07?=CEO*DOi^U zu!O{pZjQ6%Dlqk&v_&1}VoMWCQ0Gl|!xq%!*@u+XtRhy5K&5`CP20@; z=RtZPe4I5;Kx~8|+Fk(;>1tmqf2}5_W`|)ICfa!( zubOpTjp~m>^?DzrujuubSYCVn0t8Ff55vfR(?+SQ7|ZP zGmk`jA`*r#5{>4=^$*r~sXlx)@MaTG*hQpUdG1+g0jsFJP$wJ|M&f$2`)j+ISSFHJ zDCQK2ru2Sh2a2;7SFt!WL#`Yk{=b+X!Ge13Ljbz+R1gjods^SK9EeIH{_|W^{NOg4 z9I(p?cA71Azz4{YxVce;MfZ16%te@7(nIPa_-*D-pAeys$>ogvLe?5eZh9XN?rUZo zjz*YoTg~5Q-Uun!{9i;#s;IP~%fH3^JzMaBw#(cLA^Bg@;4Q2d)e+;M0-KyIshfdF zZ#TD#RyLRm!2sdy5Dzp56Lz_?_EQ_QT_gx5q%b3teIdBwk}wdHhzO1}pArVLerfp! zIwNa8G*}gVSQSkmW!6Y%OADd_cnk**vr+QsgaaRPjErYPv-XQ@7*{Yql)dvb(0Zhb zuzH=J#64C;YgI*Aig*r*j}hW7Y^9p1;Hj5Zm7?A3Wu(}P(&!gb915M1;m>TnT7;cN z%3wnKxe%H+Bht*7;3%TuRVfm`$?~C$pGU!4V8dOLAvrrwjk3E4sJ;$qY>&_Zw?Id3YVrk1K=GoG4 z?74KA)Q3BBKWED`bbvrSPO?z;tq>U+HuorS5kPi!in*1D;Lz94VWT$~!F2P5r=`)6 zzTB*el%t6Gl1fj}#S_gZfPiU|SCV-iblSXoKWFPyj9(h;{b)e*ZoKV~uMsx`0gc`I zAH;Cc52gps))RDYwUEa-ZyArh1G(QwCW~1ol)nXX@wu7LOh(AoZ{W~Kl3^8e+6#Jh z!N014lrhmf&?5af0V;rJO_Ykv(VF;5T}z#=CB@mwH^w5;T!}YFR-wPCEpt8qypF%q zfby~iIa_Cvl9?2$h?$w|ovo9IMSvIV1PhxEEW2}vItF+M7gG2>*fw)FQiYjkV*YH| z&W*$b&-B1&xvcK515V-|_(tL$uE{!EdFs#5hYVZj(|#Ynjs(9(xeU3{fQgKI%$XF`L*s&6EDSL5@N-I+0|!h-Se%-FeNxr8Oem)CP)2PR zMfRm^^-b0F32{)>*DkDdeW;&n1Zu>prLW zPv0Zpt0Dc$-=x%Zl%id2Fn@teY8=u(B%$%)pMO@ejORcwNY+6EN-{^wXya54U=at! zZ1x`Rfk?o8M>I!s1G0!ZObQ>2DHewmtS}9O(()9Bh=E~j!;Jy1FzZZxd74*qwH*p9!7pX?JpxfDh-{EuyT~~{4Ti`MRm!X^g zseYk^IObtXuU$xAlC8z?LuiH+(3wkg=2ASv$;@nq9pRAc6K2?ldw2*Yq3+?|G@0ye>mll-14Bs6?}jcGkQe$Ch2r))lD5 z{ET^lPWnlh0GvX2FPkp)_B;=V{qqSB?7qj@aif67ZI+1q(}0zx!5YG;;D-FLhw%n* zf`TD_c(Q6elMP5Rmpj;CUUe*(_ks;2PDF+|h&p)?xqXWz?5|n(1$elu&XyOVQ*5(M zI6gU8rW1yl=IqAla-tiawkhjmv215cBLff&XZB;L4o+TfGs9e+H(kpdaq)M+ZQ_R8 zLOlvF9^lyNWN|4u4FQR_nUfgzxLF`WL;}Tz-sB^6xGG|$YUTMAl@^T0PC!l)F6b=D zb-|bhp$2#GlGMb)P(zX@K&Qzq-~^y4+}fpa9^Db3+7o?J(&SSmIJ5S|hF>j8}Ln|b2H1gb8r=~hvxRO@m z;9%l1iCZVWc$@hKFfx^S&`ZgguQyUa(#vJiUp zOsdQLBQT1wXu_`0Ii1<7kOu=PN2r1{ft1tO%ns}Kc7A+Rp?EafvIHEJzCI;7 zb2#*MGjKGU`=sotsrk^!FxcuH16-Pl&2SVb$Aj!7}EyxI=n5-%wZA_0J>_-SGsBhd2ezZ zAQorsiTM2v(6Hhx@CdI!b=LegBJX#irYL#wCWOd)ho~Cz(wi5B3oLon;K%h0T5BE) zCNdiIEUn?qo=zD?7%pN2GG^xeMN@I4;LV9`{x)-eawgc*0{M^C{N5&IU~^A9a7Ss} z%L>@URgTO#uC zK8BVjDrEjVNJ3_F=5Gs`=fX<*h0K$w1qqo=Oj#kbn8JjzQ;&wl&Q%$akU37}9ZmOl zuFiuLXv8-2NR^k!ywh|ZWS6|&C9Ip_@jVP&Wf>n4Cn0k+aO7|gIM_nwy;PpWcpWkq zum%a4T}&x~yUZIvXAhZmv`4t#zK}T^7+c89)r=z{)1^e8K+*jcYlgl_kP*5#i5W9g zMkHjOqVf)9-dQ@2x?njc^MVOYkQorI;`bD!Y+uMMIa8VEcS@hV$4?2$?^Mko$zpnP}NS zAu|U}`J|`W9xo@@<0XIK1`3>XIF$4uBQ~~Wff&nqumMH&CQk%pi?UOIN~Us9Pzk#HHXF=F9VlEH-^ zg=!=gZi4U-YN{|@YbLWOeZ;TELhpT=$`%XvNJ%ReCTo7{NvUJuAsY9z!VMG)_miv` zu`v9n7_qPu8NV?W!Wdg%pqK~g_hYT|dEO4SttMOK3MtR$@6+dwvUhYz87# zxR)a+W6kOMv|j$bJ-DrtfZz%HZdfZnIOL8KFArQI6ZkMcNvExIPMLTEpJLJA@X3;*Jtz;-+S3XvZp=OG|1KGeQYptta?w>S9fYUAdd}eM_|$xk5OgKpQFWm0O|rpUk+G+ zNsR4s045!)!AxBaI1MT56g@}A2vx(Q^zr(306C)JX>tY%->ojx^LoG@l;xQrwyaR2 zde3yg^Yp_pJ5apAH;pEl!`X^ea~c)Id<|vcL}J z{_tro8AWzmHk}yha#4Vm?+K{U=Lx(Nzl^#LCX8btJ$5*p%vE>Ma7`YY@G5%b*u>e& z$p~fwbT_ueGtRT6-QXb?2Hu8ORu={${l%bi^DkuE@4mpWyFs}Zcx(NQb+k*}7q|#g zQ{@1yjorTE{{NsH_YYm9C<>1j4*6^iw&C>T!u*})AS-d07;YByaagH2U?q-G5k49s zbDhsV9^k&fY<1{Az3>PozbL#I4BE^dc9{qd1Ge24ct_)!jFtMk-^@A*fsiQGvspxf>ip$&WrqM!cyo97W6%!HN?b z=;~w9S(iCQ`r=&@;TqtLgqqiJce~*WGZE5R<447N`2}tM90p%(YovKpxPzJejU5aU znUUg60}sTOw_p>^yQuJB7GAV6gOYqnA3?1x2T{JcA8NoP`m9r2aEi}IehwQeBhiJd zOb~a>gF#mqoDb!$=MGmFcyT;GgBZef9d{nFe!HEO;nJKD;BRL)~oQ(rX(2qwv}Uc)o>Q+><)-&o>MVkn%9eepr9cABSxyE0_XGJIG3Bbd{f`C!Fv3h!|!t5(VqUQ^QraC+v7kL|1amg zt8<+Nukq)SoKLMb7Q5rxJ^0zAO>Vpb4gnf&D=%>`}H*27axfGO1eu+^SmW8#iU=!XsE`OZR}>erprMBMMZ&S2fzAEBnQ2kJ#68wCEl1 zNgDD4+VVD=u?f~z?wso$k0!{`A18tAuzMx%7KJ8qpDGvL-@ERESS;Sy+-{Q;eA3|~$a!&e=38$B5RSq|9q zB&>A{WqXRh*-Hoa_ggV5IWB@eL8I>``u6~YABSvs27)uf+jDg9F9=ioS`a7@5Dq&C z#*N=U!Y-c~LP?}=Vc-ZP+)C<2w}9E_KT{sLmtKLNxNv8kk&467(@2dtHVzC?YPytK zYmP^Wo_BS;oF*qi@{g?Wr)=N_ zX)D8y6G;xo-(g_O)-JqBm-*_=Ai!y4(#rB0v?Q>_2n@%YcGo18Bhp zeXvi#%$sjd9PcnkK~m1dlDG^?S(ld;7-TJ@`^uJ4&YxW z^B=yC;V?fADUNCVIG(LI>J7>+^DZzIy-4_CQ_W)Nk{4Of%#ld4nQib9ea-wVh$y=} z;ro+SxXObLAYX6tRv;icvNr21-jS88w}3#xjnm#SY$<{hRqGDrCu1vVz5hU-xNS2H zNQ2qJcL(B*7Pva(DA2tDH?Ke&I2%8Ofm52u8V_OiuWSy&EH#GmjQ}a{GSG8hop;@K znWG`EM>??qUI7>*FdrFP%>@*NtCzrsq7btJ2yPkTpDlvIKZG=H{F6DJ%pNk@CjHCT z4e;wik>fGb*?K)&*}QJH<$&;uco_fWV4@EoPwNuOpZN)X^ar>{Z0f@?fR+nGJ{J#Q-~VNAbB>g9{T=L6Kk^PmTbd03Z!3`mskMdu^yc9vhD%3rCZRh#4`6uM`JPwh<^VnV1;O2}_E(X1*>vXnM z=)T37pHLbtg=kADLv&Bm#?wZSL?6DEwqsh)$^)ruBdj^ z`AX_*s!wwG8Y}&i&zn?ST2t$DOfIdhb=2X@TkrFuj*`-OUVmLlsSmuY*S7d_R$fUp zt3l0|`|A9z^17P&cqNPPVpaQH{spx@m%qkUS?%}D_0_q`YwBDL^*&cgJ>;sZESXc~ zb6rtc=C5!eqqfdhT3HXt`tj$>U;4E#RR<|t^|hcauk_W`gE!gEFY&vSBBJE8mAAyI zr$}gYS(Fg5>iv~+eKJ$UcD?v3e6A5Cu9E69*N9xLo`~r6P?1l%r@o=Kwx-VSlP=Gn@V&OY;V4HsJ3aI1BPLYURL{M%VSa5f z27G0;ugXzhIls2bS5Q}1Q|Iv4E%<}nc6DIzmQ~WuqQ=RTni|*qlIjJ~)*fqg*Dt6a zRa#P2Rqv=QV~zd=4j+8nWvnW%#9vzBptGr~YpC`2)miAPuP>SFi*lGHwY3-vjtMTu zD2`r-<7^oL^bX~f7$hgczxpOju<99k(us~TA7q_ddig|Dir#(PCgU6tj} z>Dx+c=GWqCvGSyFq;U6#nTr>E{9#Nm=^23k${7BB!%~g%Y_y5z~w};13 zT=~{KdNlS(WBDZu1I{b~ud`VcjV#MjV2NaDEDIZ1NWyUpH*GyFfuzw&BO{ZLC^i@; z2=L`;AOUaS11B*T!}4$z5)>d{5}ttt69*86fWsqb12Ia%diHm!PK~---PmUL+dpo9 z9-Z!AojP@@>eQ)6ca>|4sV2)ZT|*)hqHa`pypRgol+|hH;P7H}1n7daK6W0fNX1_& zE2^idtEjL>S6N$YLibOW=sNv+c=)!a?&g8^nib9QL_)K;h8H@+OY4&=zCCg#|JHCz z&ITcWND4nJkNaGi!({$KL z_~LXrHzGvGCF%4q=pKwmf>;&|V?dKK@C>>I1Fmt<-BW;L?PnN+>D^c3xDL?O^x_zD zIzZ0>9R%$K-2!?M=r+)8pgTZ!g6;*~0~*E$0B?epf(lGJR)f;VMO#1%L3e;wf$jxe z3tEa#CssqBZJ^)8^El|KcwUXqXSU&a5cD3f)ky(WD(w+ocMS9=n48)yvj$3YL^xdR`R(Yc3vK$G<>kAeUib!c zeaB{Jg=xNZB2avUCB`pIr^k>d5zSp1%|Crw&9Q?v~@r$CB7h%2l1%aVB|{cg zMRPa#*G5SX3i};?l)hWoIKqx|2_t=GA%Q}Kkqo+sE}%96>?8tYI;&+m#9xK5ItTx1 zVCT5-lPr3t)r+uI8Cm`s5=LQE_OtL!fLqv3gzZMy@f5Fu@_M$^f$Crn!X88z=_9J5 zxyunk_$csu95|&zFN=4&$f7!@&v|zt?0dKu-=cI*kH#9JxmAxu!wUcvBb>g;lHt%} zvmd6hDNq$H^q)O7nh&HR8eSF=#Jd)}7aajF=`;-91K?eRdsW6AGJh1d17X={-)>=# zBCOCOY!qRQ2x}#ovdsKn*Jb?)!ou(kZeg@Iat^}2PVr=Y`_I(ziV@ZenQmcq2)hUN zNVd=Et7QMkS#dt#(&R?iDPv;jv?$s{KlX!R)s`xB^=Dxs7}wEit1cS zZeRs|%kVoD`CW}N37!HTq;Q3QNi_fRz>;X}vfz^Ff?*>OE!-Tch?Z{3s*cwA-<}#R zt%w#vJ8*zg>2!KT#2{-Mekp+9~kSFb=M%f7%@J~F5-ynF-01eWyxsqp@;z2}; zJBr_};CUbIBlsS?hMQi{g7vD3ZuJFj@Vjb;(o+8Fqo5GwcdB$%9jX2<8{!9nuLYmd zrHXV}8jYO^vqqzWS+9)LT1>oFT6+~Ln96V%ZE73hZvqX{|HJ30YZa9tYBUh`K~*#y zMgaYO5cEZN79oOTwc`nULEim=kKNO*rnD(P^nXT4OyO* zs%?z;0z>XH;T4t4-H`da7t-l0k~wr_GW*?SCP*gLLpJQX^rCcn3TW_lCG&D>!Q zza5bI2Dscqb^WDfh`$+-Ti;ojCvmA(yL0KN^SUygXAh?jz2q*FD5+iiJ*Ef1sK zk*8&pC;#*Cj#ZZ2NZ&^7<3vOo#aQ+F9LTs157N(6=7g495lCDfOk8FphC|CYWd*(z zD9CLn5P^-psTc-Tk?tqr*NP}(8`9}Dq{CIpPi^uqkH!KA(MHHCBqCI9!{7_=Pp79+ zxm|=mj`p7*A8NI|UB198pHipPe~v=t3*X?hYuwWYoAlq=2;TN~I{h@tELft_zuaGK zum2iV|8SEvaFgVCtQH{4jPkr1(I3RvPv*JBJ6(*%vLh!6ahkDOoC{s?{Dqr*^BxY0 zr+uMIgJPTC7z>Dv0pszYcsO9(6%_9V2sso`%;6J^e-fqeb|mQYz`ke{ie)M+c}M6& z_`6%C8IMjA_fG@(&@^NJH1Wzbe4hWeX~y?+#TB$S>bhLxrd)9+-QJ&T{32KUhHihG zYy2%&{E==C<{CrO#l_R{apS*CH(s4C?wt}$M;>>!EngO^SiL~eY1_x*J>P zpGT^MxG6giw@@-CZY0GLW1~-;OGSI*4xjHP-xuc#(G|OHmbk@t(#{#;8Q=VUGsGvp zbyv(1@B59LXNji*215TFXt`&GcrEz2|D&1Wkt_obKMt25_=E7ApUo2ApOW+3Eb;hM z}A@c(3*0SQk`Gr+%*O8NE~2GV_T z1~EQCo4G%jN%5bYMI4xkUPy)v?R&TTLN`UlQ-0$=bHuBDBb@{H-~OGb_#kl7^-*zs zu<_S9;-g^dK#sUMWPFe<9twqkKbM7X;J1g3yQ1R3uyI>d90-Sg8WlrRj5l(`4O7m! zHCt?q7~6Bi!x7_I+-Dp2;GQ~tq`4kxzLO(H{KgYG;yJ&86!!VspUM&Y11J4HM+^rW zugei{2UCc0Wyp9STihKA0sjT1xiw6R?hG4$$`O0Rp^KwpENtAFBfdAq7>kI9ryB3c zl)s-%DG%j{cO!&fPAQ{qhJ3#Nb5`j6tjO1mN3$Zm#)uIaocl~xWRw3i{C&x9{KybD z1dQ=4sA60cjO-2Oyk>}xgEd0DV?e~_knw6rTpr4KASCVz0eG4I{3(?4P?oqc3vJmD z!3zTy`Y-fds7hGu%#Mle_zmN?6+dME^S|;LKg|*E`HknZ#lCC(NZ`|$oUE{yR7kXxP$*HGM~?$2MO2Y7~r`##{kd2=NJh6bB+PdYoZ2t?u{DY*&j9DLd%7hni1JQ z&A8>5$PaVj79#s|jkj<=-T1d-BKxNsSImw4aE7t**vL&Ya~{tVKb!dlA%2sGtdHdx zPs|bHc{#Vj3(NuVDE)b5PR@_>#pe9-$rbqC)R-8>?@9cg#Sg_kn~Z6@&)5e?9-1Ov2>invVR2PfQUHXBqc}#rMKc=2iORsS){vqquS#6mx93qC;mWI=WKPXX_RH|Ml-! zhsxPdH&cFg<~P=cKanaat6|wOv1}$IRdDS{~C(fXB`D^E?7P;9bl=*5%Gn z$LaPSVftH-GW}lWU(4w)b5~;$ilsYwX*!9)0j%B^lgJQL9qg6fu5%T~F0LP)j;1<( zF6&>;d|vW2zouG_rW*I6`M-DX%@8K1Sy98VM?jX|KF~#NiQsru786 zFzr;=XiQQJW#v`J1xRa#bm{VBF=F>=c6?7L{zyL`Ee*dCbB$8McD$77DyHpB z*E8MB^jfC3F};uJV@zLU`WDlJOheTw`8=lcnU*kJ%5)XecBbo@Zf1Hd)7zNd$Mi9# zFEV|L=|QHUWt=|K`AkcgE@irkX*<*POgA&Vmg#Lw?_>HH(-)b(#q=Q4Pz|TgbUxD( zrc0TwV%pAhJ=4uhuVs20)BBh{#`Hy|Z!tZ{G=yP2U3pCBGc94dl<6v_?M&A*-OTh_ zrnfP@kLhDfUu60g(}PSy%T@YY_Q%xC5&o^Nto(9p!CB1%J*k1%;)3FW!jt0z@+tnU zcxgdlaY6hG3O~HS0K7fDH*sD8;m`P6Ygho@Q~Oyv9u(Ri+i^^S*ss~~kkJ0nj%Nw& zSM7K>9gm8&+`rlJX(InIE28cL z*u@v#rB8IraIE14<5y z7h0b&2d?$m{$A&&Ykk(fV^3G>BmQ8=!>9+h z@}d9MgvzJfG2fld_>co%z_@dMd>U}77j5UdJn8?g^SQLg6wA@`aKYD*w7ysav-#aA#F(S~N&}9f{8p)5VTX$V2DnMh8Z{EAu#djCz=>x`FwQ$3%dh?tT^ za-w3_qYBVzeGa&kbFRYAVg5?S*WRiCdJjq0CQxjxn{?gB{Ca&!&(EF$uG6wcngM;? z@D|{?XkSiyI2SmT!*=e+^?dbW4>>>d!0+_HpYg!wPI1@sY!7@FaJTY#z=MCE#IuE7 zkDzz6bd3Y2`s)3qde$R`ys7T`{Jq4pgUD9Pr zIqKf}_ZyPmim3Z5Jn&~c@Yg(W`kqm?-}98*<2k7>0(Z0LQsx)Dt|q@vLY%Wa_*=F7 zbCqEEofmG0Jox{w2mYXz)2HNU|LlWZxas*-5Bxk2JPDl2VWd;ZKS>dBIdC^Sxxquu zEi6Z`_vW&kUwZIA4V>&`m!seL1LIp6kFlI!)IGl^0C&^#YaZo!j)$CYdf=CM;8%L! zHv@Om=Q$7lzj@%tO{4ryo|l~foa#%j>z~Z!bBg5W2Mb9szg`#C{X!>jD*ut+#w7O; z$GOBq&Rrh(11zV|(cWM3;C}-PrVxdZ6$*T~!TYgVx z{%6@=Y5A7|$G=I}jm%%ytAglzKf3M(PWqHP>iv1(Zt3o0{xPoK1xk$gQ1WA(#`D<| z7&joqP0lQUlwbW`N9+G(#$%$}Q@F|Dy-K5ID)->X2XS!G9bYIN6C_ztZKl9sE>2`hBae-(LXF#XM#A zyD?c_`p;v`Kg8vEDhC9RLDDD2`usiP1;E|%yMp=0cpjnk>G0t1^}sIxF6-Cf53lp! z|A`0wZ@@|aXIX#!{_JIzuit03a=!l4Lryy~BFmHQR_E(t+_;s))gJgba9RG2`D_qV zbdv9!_x&SqH~D8WzjHp-0o=`Q*K2;B_t7`TbPdBW$bOc+9Fx!JR^;_4_zCZ29BRvJ zzr<0VT)(h9dHt2;JMCvNCMItBH+bMF;8YI!eY{Th61MZ*4n3~|Kk1pz{f_QO?(~p< zKg-weyY;;0c@O^MU@&g_7kl8;&#sFW*H5-Q-`$a`bl?i@2Pxmg#1T zwg0LB`c{;#QQ$FoXaDnwIMJv1zsp4#wm0BI|ztnwExCJl@jL z)JNO$J5&8NEs2IE{7F<6CF+*9)Gu#XXclz#_cawPDi!9M74di_cB3YGa9TiHKY;30 zYl=#(`%aykOcSNFVb$#4&<*aR;dO5BZ0fJ6Z7@q#VKaKNuVzJOyf}g9_@cVLHtaoL z-_`&a8(*lQa8R@c=hDk+Z2TD%mQ#5PuJjqPdc!;uH6rgaIGw~{iImAEY{ zY74DB1Bc_NUa@jX#R}8xA80mfOxh5e7^JQF*lvuI5X?0VsTHNB*_rI@=}h5KWpYK_ zOiVm1s&7dSr0f8;%L1?qr+vezrU50`9^|Y9$`qM1TI)yQnuKg$;3wIPP3`6uW5scoymP0 zM0;OTcbnNd(A~WOQ7lg;W!%PUnoH{|&PbSv+A7mj)!?Z5c#-u?iNyNE!uEk4h^?+Y z%SUhy>SsdUDhZ7xOOg zjOBANQCK`~(GktMQ-`f7)qBzFmwyH)k&tqhERvqZn$d>TM0rAex8kIo*Ib<33^H4} zJu1BFYqb?;)KsG88_%fB*o!V2#{nagiov^p4r6Xq{GJ6!{!dxLr>C#0H)C&^%zEux zN2*T6+)Y`#g_TQ}Ho)^WRxDYO;D+=+S-aX)S{8M)^*`Hmu1N=_>`)m#ClQR zZg$~xgCeT-M)y{PpxSuh#5pA;t!*hBnotBQrJeZ=DI7Ag6i2`;l_#Kl^v&GntI z+d3F&p;@NzioSIT^iG{5O&&ckp;;mN1e_gnIu4XbWDJ26YEYaL#z}6raeblLQ&eoa zDq9stp3F51z6m0&6LI}gsYynGl^bya<`i( zFRxxvv!rtMYI9*hF`_tCWlhjbt3P{&+E;}FVtL#h;oWFnJo z7=}2QW%!Ymz@ko^^3ayRASPAPmBcX`Z8ptSQb#CipuH=3zN;=rJF+Un5liaO5E@~1 z+f<^$tT$^ryH?@<`fGEP1vM`l_otukF=6iULfe@a@>q&`QxBY0HlK?zl}(l!`B1&k z5hhEix7BUT-XT^MQtQCTvdL`5h@!`Ap~F+k%VEV8OKQw`!9rxGva6}TzhOgnbFz!% zxU~F<7G0J~;uwz_w5b{tP9lL;PtDIBNuAVGb2@sE&(di@a8)lgObCuIq#hg!;>-~Y zY$+`_8p5DG> zN*V-LdE9X#i>!0%X&M^qYig^_#EM3$A?1QC<^kc-Wx)WFbg0ij+gUy5!<@j}kSZ!P z%~VHU@_Z8ov@X@54upXc_RMOz#WaaPX?&LICgm*d?g-~UvX%DtcdqN9N+`Bh0yWd_ zwr+I!_PVXCph2iAipfPS3$EU*mP71BPo}G!aQ2hq6!&lFX`x;<(F^aKYKQp7Rjy?% zdoq^}SDCJCqU#2l`dZDVX3R$%)9W@> zl%T>OGjHkLFk$BLVx(5Ds#^CPnQ@gP-^%nd?ZHyGvfIc6ID046tLR{-`%wcj*G7hd zE=}SDLD^SRy}Gp(T^kq`ENY=~2bMcZG4;g}f^9vC7VS_|H0hkkH@S0E;iq<9@j~vI zGTqmDR4|+e{>hH4-dutfLc`{MoYN)ycgLYw81KM}xmCz8w=%7ta4mgRY_|%bX=mv;%#M*3bU7A?Q+&(qV9^%?`TTNF)_UA5gOTGKYCg8F-{{h*R*!k zIma|y@0QWL8tt~AlUD0KiHY$wCIG6MEgt^p2DI=kmT|Wn%~F4HWH~ECn`3ixWCQJ| zv!?^=OsRfqs;ZYGRe6N$R;>taw02eh-HcSqH4)uXCaO@*(_qZ1awjH3nX2U_jh2q{ zOqxGTv;uThymEySDlyHZnl5Oc<}z!5i%Z29sv2bO8cb3P%|t_qSx3u@TD+S06wp#! zt+j2&7s6@E#wpWwjDZW{n0{!Jv5Y>@%JqcFqZw%rKbWZ-x7KZK)8wwvIw#FSJE?|~ zJ?NjQqfjNI1~<)3tr%DKryRziM$9ze#CWT15XX>Wfz8p!!F5@CMH(>mcOI@|)DzrE z^Di<%%*th1H6$0E54&uWV?LYp(_np>9K*N?((e4QGudLYkwl$)Bf?szJKtaOr}i{Bj1dDS|Ji8E}gHuN-gcebd>k#rPtAg_n^8MR%P43R6(l?xh^a3||t z4U?Ar6N~%X+Rp1u_9q%r#qFRq@wjPrH)HgMX+;Vs8e$EmBTz1N!f|qi%8orV(T$hEIl???y$EaKYKZ!or!XG<1`lR3 zR#`aV*d9ozU~!YS*4nX_pmq@^@gR8Lfy#i`q@ z`Y?0oMJWfUoBs8tenC=?cQH7Smfm%FWl7J`sB>>fIuhVm z7F$uxvP7-1D49p{RtpOO)_;~|&6>tD%!xNTSR1GJIgN6h*Z8%%#LA@t(*v4$W86;5 zqR_Ci5{rz5a$tg9VPPh76KcZjc_3=Fp`(~Qw_qt3eN27#c)Hfjh>vuE!Lv8`!%v+*!Uo8m;?x#dg%#A&bB$|^awLWLjU zNYNf$nu_Hei6hip*PZN9u8ZD@Q3Y#zv{kQc;72%8?e9pQPe&wUHNl=*y}7h&puYo) zH2q1-)F9JSF2odgV!%zas>3gatNNQLK+A|-uZY_RTY5D`(Q!{9)Ptgv2 zhjW(O3fA=u6f_TXcD0_2^Zn&hM^k@?C}`c#gP4k9R4YR;v6ctxTOUk>Vd1BVBJf>r zS4tGfSwaDB3)Us^h{-hm%Egv~zNCCvP|()FFNQm?zLI%DQY~1{k>@w#^)^jbr7-o* zSjo`$B!^OmNU13er4T=vqJ)9aJBb1e0z?7cV$**mpudRtZ&7^iBhOyo&vG;Sll4C3 z6CKinj^_eAAH*37I=((PLenQXKz~;pcEm5ma|Ayfe}ogzG>q?2=^~guLDTW|_aNoK zNVblz_eE=}&!qs1y!bOV{V30|qcb9Oe7%2K(`z|`)?f2!dKM(nnG^b1@3Yo)gi|8@ zNru*+VbO{k3e)lRer!$kITt#89bf0ahvP3{f%@DGP4zh&6d#`|X42R32Z2!*`E#Gl z{%*a`o8nk6($#wPQcAzTVHS>3bAt zxt!_$2w~Rvdf&UI -#include -#include -#define NTYPE unsigned short //fine if N<=65535 -using namespace std; -using namespace Rcpp; -// [[Rcpp::export]] - -NumericVector getTau(const NumericVector ORIG_ID, const NumericVector x, const NumericVector y, const NumericVector onset, const NumericVector r, const NumericVector r_low, SEXP index){ //see how multiple vector arguments are parsed in here rather than the single posmat matrix - IntegerVector ORIG_IDint = as(ORIG_ID); - IntegerVector onsetint = as(onset); - NTYPE i,j,k; - float dist2 = 0; - float r2 = 0; - float r2_low = 0; - unsigned long num_cnt, denom_cnt; //counters for those filling conditions// - unsigned short serialintvl = 7; //mean serial interval of the disease - unsigned short r_size = r.size(); - NTYPE N = ORIG_ID.size(); - int *inds = INTEGER(index); - NumericVector tau(r_size, NULL); - float piInf = 0; - bool check = 1; - bool bstrapconflict = 0; - bool sameperson = 0; - bool iscasepair = 0; - bool temprelated = 0; - bool withindist = 0; - - if(*inds==-1){ //if index was set to -1 then it means we can turn off bootstrapping checks - check = 0; - } - - num_cnt = 0; - denom_cnt = 0; - - //calc piInf - for (i=0;i= r2_low) && (dist2 < r2)); - sameperson = (ORIG_IDint[i]==ORIG_IDint[j]); - denom_cnt = denom_cnt + (!(bstrapconflict)*!(sameperson)*withindist); - iscasepair = (onsetint[i]!=-999) && (onsetint[j]!=-999); - temprelated = (abs(onsetint[i]-onsetint[j]) <= serialintvl); - num_cnt = num_cnt + (!(bstrapconflict)*!(sameperson)*iscasepair*temprelated*withindist); - } - } - tau[k] = (float)num_cnt/denom_cnt; // pi(r.min,r.max) - tau[k] = (float)tau[k]/piInf; - } - return(tau); -} \ No newline at end of file diff --git a/anRpackage/src/get_tau.o b/anRpackage/src/get_tau.o deleted file mode 100644 index 771e71eee62427e1fde4abcbeae681439cd5fd21..0000000000000000000000000000000000000000 GIT binary patch literal 0 HcmV?d00001 literal 28144 zcmdsA4|rWgmA}bNX>3bw_@_bYl9m)(EYG|)EeRmL*W@+#ZPPX+X{nU*@{+tXku-@f z@0EZ^t7(xCBDeyw_^+&}UBPX2_k)VAq%F{bsNJ%=`UxVd?tVpO!SY8C*z7qo=e(TU zdnbRsf}gwRQ?6G<=GZdo^^?p$2TJFZ5X?WsBFaQwzeqV zd2o+1=fG7)`2~=kOZoRhx{&gVAUzMR55TpUu=61e!SzA7$_WcYx`gseA+4bNGDw$G z{sKrVDPINYg_OSt(u*noAxJNwe2mhSkk(S(r8ExdUsHY+rFD?5ru-U8*FsuP`AZ?i zuwdOA;cBA%I!M>Ubs1b6;Q9z$m&1j88{xVFu1#=V30DHHX1K5{NofnEDM;Jlk`%ox zx<0xtKKNF&Xw#D`+j%uy0VV=P?bBFDJyjco5~itgLP%)b%Ql!C3TQ1Jcb(F zp^M^!hba5~vV8~*9&87u8U@x3VVhY|=Yc=Q55)>GA0K@5p(0F&3d-VxF-!5S9@+we zkJJsmpl)}LI6J;m06Z5hF3S|Ugt|jDWgYI|Lg)S&sAvZ;Jp)WoJbDHKbwgcP{z!c2 z6JRB#J&X^nE_8=$9>7xFZ0BQ{W~@$0d}v!?_0ZhU)D1Nk9$(=U7yRUS zwNpIziFPOSv~yEyzGXQ(H^FcCsB^Xz9*I`Xy~DX}B>b4O=aGRI^MI@1r_NpbqfYE` zci=#tvuD41XFPADuIlmZ(fGjrJa8cU+4h+8!qIl;g%R|j`1oq)yeE&xo%2SH*Em-m z0bfR)9lru^-9t5H*z=)LYsv~=s&TGAS>Ht4hdK)JrW}KO1@wlY4eGuW@5>I>4Y}&J1fzfY{<@+5!sE5j4u{%noX|7QO_zgd zXJ-@q4mFh(hF=9!*pAM^gK#|o{i15|{l}dvYF>egT`}h+sLF_HCdfcj9XRf6*bi;> zNcvsj#~^*p`TYLyXnWKNo&JZ2bJrngzBiyX-GN7;1rJdp9*9@HnEg^5f~acm>^ygL zbQD?<6udBk?e~O=E9OA|m4`8y?(wK|<>SYlQ|8zZ{X!SQR*y|k>FnJ8Kck~4(BTd| z+`n;C^or=F=#|mLL2MVUkLq9dLyevUJ;ifghPwOJv_(YfmYxD4n4 zxAjOo??@PWSLnCS=MJdasa^w}60j$9$HZ#wM{ekw2jf+L$lif{^SQ)ZqoXR05R38g zJPPsr0mkz!70)+xJVCv)2f~j*4T|W9is7cy6bM1~Q513q4uR_joIQ^^dk(v;`{Q{J$E#j*b~Zv& zV9Pw0gs5VBK+Fz5w-rO28^N-rbyd)xYu^|hT@AySvtt?Lcz}CuehxlY!e_#b(ua@i za~t0Tg~K@8=&=cVbrEp5!|D-`MU;03RTb1o_S~=kg#H}Vtv*y4zv$w!%nViJ4*pW7 zS~)PO+u!{wnzF2e3qt4h>8$!~%xQlzJmL->cL$##H9y;U5IA7I(1R_M4m?~4Ve{M% zDX$=!mC!AV-ND+!&k6PLFsgExdfsq{x8-=))s5qr(cU|wkL*X>;rbHrGXj3XprBnJ zybQ-PwDp;@(>x|=mc%7mOTf$ZOO(>Fc@wJcFP?GU1|a! zJ>m|$mgm&`?sq>_l2~e4e`)@_GyO}4YVeoTn{OTK&u@2s8*d5Lz3<5LGTuu2pRGP^ z$xO00leQXKdwZ?X(Y!O3EwHRl+&t^4+13wc&v<>d_3~`!u~;w%(p`BoZY;2N+B3dg zVBJ@MZ!j1S8m9+-U#x2}q$UM-Sgm_$Iv8IuzCQA(E~msTer)<7gszI?j!l?^@G{-nmWJ;^@}V}UjJ z{|&lRtlmDHc7ZsZqL`2E?Y8%J>+?D60p1aNmA8y&Gjk(4%A?ynW?@6Ay86P>h3i_f z-I;9Z1)+*i`T60jx(i<)t_+n|gu;t7u0&5K1$oytSb6=$d1uX;S$HcffQkX1%Vm_7 zK#KcfNNRkJ78MWVZ!9W=d(?ybR7fam)f5%4f)d2r@E%Cwn7>Y;mW6teFuhgfENcO# zF6!T)^dtX$R3`aH2=nqI&(lpc^cjZhix;ZYqgWx1y( zD0_{{a1RjcF_53Xp5;!5{S+i`*<1}eWo66%7|yyNtRWB`4vUQ zx8$!VD*24PqNsGZAXc>S<{8nV@*7Ty7F7<+tS#D`XaB;URa6OBv}hsl058a_C@PH2 zCVl;+^K$S@=~Vu%Dk`;K&6}kHu2_j{e1iR^!1}}JcvRonx`lm3wsp1-L%_C54fs3* zKHq@j4z-_b-Rpe?wuL*_ezL7aeh?5_?P~jrY)fq``U-5dcjDI$hXFD>}yEa=ku&m z>)UjXd`Nh1fcRbyUP*W}G1gNWCma$*jN^@MRcp}mYg@sqC)>i413%eT@T_B7Xx&e? zg(nPtvU4N(_%hqV(+5A<7M@1<$+qxh!cVq3r3k{;`9VNzb%Nn9vek)35RRuCezMhR zM-ZOygMirTB*b52TP=Q2@wXaqHLCl{Y%Ar%5U{Nd1K#Pw@jS<}{C0JX{B||cYQZ!( zjA)L3R`XLbR=i7PXN&=qdq|})Sp<-6Z4N*&K7RekXu}g6 zL^hNJk*#MFUy*J3{X{0)GWj{vr$h65qpaC%>sUrd>8_@b*ceB&4Bajb9}LFskhtlAl%v3yL}h}w)GwZjw>ra z+18nU5D?pXuK_>Hfb)uItoRwo&bQ{}0Pd+LK_&J>k8&Sr0_*KJsiRRGUz!C9yylH2|=a|@cCkW!g{I|jdTt!uVVmN`mUU z!ho*>e7;yOvE7dWK2OYFtp8@j^T*~v#_v^n#@1KdPmTf}C)|Jb@qy{};1U;!P( ze+uC9#d?hSFHrc{dQ8tB@cn?o$JS$vZzH^y;;7zBApZr0Uu1dl2Nb^4Y9@a5&H?$S z6dv~Bqkx|(*6sM7fwT&S?0NYqeFMKd@@@t^sJ>|feysuDX~1tY;9oZ2djP*XKk4)H zhYD};;s0sS^LqpSrU5U8I0ePwYy-a3fUhv%4FSZR(LFj^URS{_{tp4wLDm!=UGY0e{_q*Ta%1 z$e;TR_-_pOS+Il(($i?bKW)Gt0(_o$?J)~iFjgj zDwD`0vmrddNnGBTsjEyRI(s_1J2P;rH=pamiP+f9vWC{4Y-T)=N~dF1sIIyKf0p6T za{O6>KTGkaJSSh-*_}zHyOUkv8J$<~zH)76$r*DQ-|&GC75=wZ`)qGH&< zA6Gn)l_kQ;M0{bv=l5oLe_wAl(;cf{)!LPd>DrYm8)Hq0rs#^gSi)+E^`@b=?Gx0g zqVJmSR!|MO?(6AFnFOGgEWmB!cdSO&uW5Z_CcLaA+1J^c7(41{1u)#& zkxVBt>11c7FWw3VAE2u-v%IUPHQAM_nYh-I6_=}+wRP1S8xmVewt%MsSUiCXSgdUX zo{NxZ)B3;|YAy?w#gdGjmO1Yq{0Xp6KsNb#KmeSn)MYvD#PzJmYIRyVfVWvMDT0cEJOut>bNH-mp4XLgabW&E9XzfoXHfNLR zwnVb^W7*C$xRR(%WvZ*eon#A+Jc1>yg{{Q8?yEby+dx@gcFTB{WJ^z4O%#op+!Ij| zhsmRw9%3)98zIrw*OTZ-cDLm$CNrsiJ>>|AY8ogeeI_e#zKYcck9ahD<%T$Njh{+S zzsILf4FrMD9&{sB+89{MvhfK^*RKheW0?DsEuCA#OCaCZnP}}w_Vpz)*Yu_mo!#v{ z%1Au8K7k{6Q?N(uZ?H#PO0u5OOZ%X!5?G{1%m885OID>Y6~TlvwaIB&CmhkIVhyWu z&4&ZlRHc{0QvfIYv3i`Nr!ICuUn+G~PkVc;3k)qu#9^p!1f;Dq4fjjc{o3AC_m-Zv z)U?(1_VmRXGZpQS#>3%6VoM8*$C*?j1Lz6puIlUPxq3^o`x;zLOrxqWlfZZg%b4j! zs4CkH3pO>SOv~J|?DS=-dRt(^!nq2kugja7%v~~|NOd?^1uv(SFuayqpM(oOKUEYp z;Jp8V{GSs1yrMo2umOQz4|ycy-vJlvA4U)e@ttrn{#68l5WfK~#=nCg5WZS|8ZO2U zA_#=|jc_skGX#N9kE~Z=o`MxV2wCQ@CLHzPE)C;r1&*;|oX>}lA9raOKNA5UOe8;_ zGkNuYSMWx4>@#LG}~(Ly(YP;$IQ^Wt{I3{6iqadj3)H zR|@<=f#WG5^FJi?+$?bV*38$wKN5Om`yLf~urIOx=LIh7^@_khBlvN@1qtIU>stW( z3`iKCTi{~-)GUhsLV@#6IF{qt4fB(Fr6)NB+y)*be2&`SVm(rSO5i4be!oY5B>z0P zM?$;yfp>!d{u7B0p8ZWHKGN=q#7EkF zyW)fKLc;i)p2sYdaooz083@}w0~g~l1cAJIc;^l6NIf_=At8P>T&zdx$2pPlspyyE zhDkr4BzpTTKi7!23A;F#A)%jg9JpAay1w#xSxGqRk^S~E!cmXxx7emgUVm_mMe_Rd z6%G6R;r`%#j`kA2*Ps1DkMyTm;BtI>+yw@L8JcEM*KiS_i8KOfa_4sXlTz{fl| zf2+97gAb4%9#2^R4#K%#F@7`Q=ra=IUnE@4ci6sM-@6sJZu2`J&-L9+IL2T0=RJhy z^i#ba^x$2D?-zR7AkX>_d-$&-{-1d8n+gAw2WS7E_uw4QmpwSg|8>GK2A_wE{ozkr zu)cQ+d?wlF`m+8a!jT`_l=Yk;_+>vV6*$^ve-;pqnC$ z>53aZ2NpQmD4~1#euOw)Nb5cPTrc*Q<8Th~;dc>8jF-`Up9kmqT`n7{i*Qisev-2O9=PcjT`va2>sITdVx#31p=4Powlz39>aV1 zMx+=P97~G{E2I?hyC{zdmm~fqfpc~uJzr)fNKZu1E`y%G6F6rl^5?6}1gZ2#>XC8y zf`R{|0_XQ!w9D^wSX_!fK!|@$<@C=&$|KHePK3$@{te1wQYrA=0*?s1Mc}T$|6btr z0{^DKHwgS*fj0~MK7n@#{96L=6?jJA{Q}=A@Iry}nj9h7F71L}#(A^AODX^7h;y^B z`j1R6` z7?=7HkHCe*<*X;7(DCuPhcHYeKhBR#rIbNE*D(lbrNE`0h`@Po93dH>YlufCsQ*?0 zlu|w)(isA;5co`i*9)9ug;0mU;gg&kZKoXSpOr)P+7fRO=XYnk--|y$h|f_u{qr#8 z5kFbr&kDRq;BN}t5jd|mQ4hB%7FS?kfDk`L<@8Upz)ux8zXu@yX#(FZ_<7BX#SeS% zUXpoM;HZc3HwBLVHxqve8U(@w>)Q+$(uW1UT_O6XN#K-K1^9=7kT7|-%ITlo0+;Iz z{&z3Te`mWiMNw!!a_zFg+9h`#$?qmoKo-Hw-*cJ_r%&=*4fG)#3WF3|o}%E?B>3Z?$<@P2_0*Z666!0L$<-n!eyuS%Cw}vx z<2P9(kkg>kbFi$wEeJ3z*-#sPnh?O|X*%GPO*b_O`N_=WCLmms2uRl_2(DxsiwgSpa{xJz zfnMDd{K*0FPofp7IxVfl@73xF%fu=Jm?zRVsf$vdnVG&+^Y=Z2v%&mk{Q;l5VH2qD z%L1D%;0WE^#sS~inVADzLBFTvin*czu4rZ$)1f<>zoNoH9~exZQt=be+Lx)!Jq~aa zUheJo>I2HqL?4BPT6?x^Np-_?)7_Ivg=*KVJ70Z_N19R{iS{&n6KjPwr+eUOd|T+B zP8z>nV5Jet@EQrH#^w_>1`hW2>uMCLaozsx736RJAYxbK{lmGnh^5=pK5{_TE|G{#cn&i94zLdv( z2A1R0B+q{bM8^Nq(~!TF>`VC@ry+l5K>W4`$d>|!tS}9ZqrD{0|1OUUi#bE-S&(8Z zB6QEmzOc@f+d9SgzYi|1Pl<+ef76-R2(_#S<)%dr)xL^b^*5hUr&oT2K z?~y-Np{2@67~_8qT;BLU0a%dz7sI`uzxNX53*lOjo}99ta27cgcP^$zlB$!R4(#j#EMQ%ZU+V zCdq4mmBD^L*>_18?SB9+uYEi}4YJ=%jFZ(rYp_2OQ2+A-?B8RszkQ1J$8*@A_`euX z{|^S($MaP9&-eW8CHvopGD+U}KW4CBxm-(?qH-XNKaMrt_`hPX-%FdYli7dKV1L5} zXw|Z${iSg4wLcfe4Rp~v{*=>Z?_}f8dq6%Y{{59`ee669>yJ(1wU6hlLH6s({!OSv zkT?E#ej8+eH`#9{S+tKK@!H>Ju>b57?O$uKUs$CLH|s>TKSEyn_&s<~{d;F?!5@=6 zC*JzsYp~xz_Hlec!usR;tJnTd4fZ3n@SDv3QG@+`7ot`8e;Ww<`2Od$|2owwrG+r@BL-^5(OD$u)ywkCd^X8@*FW_n&woFP`MvRfFUYI@#|)GfUaT{P zG>>BZqi}iazXas*8~~DWw3iM&n4b&zp40OU_KT~v1drcnA0d`w{1T8?@$&?5gzUeU z4O8l^e<#V?9{FcUp6#$2w$Jx%ARqMnE<|U6T%D-?yAtv&&w6h&_|Jz!lllLU!T(aK zDBB@7a`yTDL4*Ij0shAW{Qs@N{|MQaXEg?n%Rxy4|A7unyAN(8YyLpK@hyJWwE z1{{3;E`^KbxnDhGuwQzLwy<3%xnGuN;6a1^y)n%p?Kc9%@@)SFgZ*7(zn%?I%JK~V z)?mM+R&$gSBi4T%TrAJ_OW_46sQ|IZBek4&-tPa5nWCi^n}_^$4?j~gF|dB@LE`Y3?; zfklz`e%=YnM?fCipQ=k%EQ{?I=MbG?`<$*e$d{AC^SFfd?}Urxx&G*n$Ar=lp^rMYQ$uwD8ua|W3vgIM_-N!@4ZT9| YG0fntKMu1&@_ScliSI%EBzfik7cZ$g3;+NC diff --git a/anRpackage/src/rcpp_hello_world.cpp b/anRpackage/src/rcpp_hello_world.cpp deleted file mode 100644 index 98a959c..0000000 --- a/anRpackage/src/rcpp_hello_world.cpp +++ /dev/null @@ -1,13 +0,0 @@ - -#include -using namespace Rcpp; - -// [[Rcpp::export]] -List rcpp_hello_world() { - - CharacterVector x = CharacterVector::create( "foo", "bar" ) ; - NumericVector y = NumericVector::create( 0.0, 1.0 ) ; - List z = List::create( x, y ) ; - - return z ; -} diff --git a/anRpackage/src/rcpp_hello_world.o b/anRpackage/src/rcpp_hello_world.o deleted file mode 100644 index a446de8f072afeb5d71f08c8242f843b9a0a4044..0000000000000000000000000000000000000000 GIT binary patch literal 0 HcmV?d00001 literal 28504 zcmds94|r5pm4C@gAf;^*TWGTuw*f+>EyNiJkQCY`Bm*xcK$?Ux*?z`_XcZ3>N=jPN zHv`8MVIDTZu?ddNaBP9&-EeG$1M9ZI@g6v~!|`4?!f>?0fprn$ZNxi(N8ylEv+4O< zIIr1Q)3VV?{r-hwIB9l<|GTwe1XFHxI3s!UG1*S4u53$i%sEk4R)7`jMi!Peq!yN) z#*>EBYfkE@DjpmwI&dYJtY0y^ELq^BzV8Z*#>$GH>x`qFQ%-7jS^enk#hs~^veFmV zT)!BVw!7tbwT|| zUYT8gz;DXtmr`KJj*YE zQd+(gPB_AYV*_Q%cK07xb~BXHvUPA$bagf41Kr|Tx@tL`bQWEG0K8$ou>M9FRsc@4 z`e&s6XQaL}iZh1WCD;dfQ0a`~tYgV4Dmy9F%PK=>z5{*^Ng>q#EP#{L-w&d0XQqzV zr@ptYYN6_R?y!|*^($7E(dFc%`|)-dHOEDb>iEococih6V`u+19WEMHKF%J-Ou-}| z#sppa&Pc}Q=Uw5?IjQ4kc42q%#SK8_;zplE2f9Rep4mpi>biA>2VOq=tut<)A3yuG zGbnCZN9UJc*N@pHYEFmxGE82npF62nDJS3Fb`)e_sw+?c)D9mjIO;hcQ9a@epQv<( z$BLcQnlGPW+sGZazS&9r&KbQUEki!Ozb#8RBmZXGBmdRTXk!V)Sq*U(P-sqSGc?T& z-3@EhZBoi+iP85vu0*Zl2HaWBo7X{e(-oX*>LmeOzg#* zzdyhDCmLTH{|ScWdAy^RWZg;&f@4)bO9mmksb5szlPrTISGPJVRL>0D;EcAEO;?jI z3`=(gq3g^UX(%g5PqeAi1M^_mmTcP&v+Iz`=P1~E4)0^s0QvA5FJ@mpkp50_$!&|AmXJ zZ`jxUSAq56yemFbVErO*0a$o6e*qvb=idzSPfg1Q$WJn{$OQ!prSl@t&1y_nWUJrugS_8 zEY7*=lIaEboR0e%v$>9V3Ge`w)jV8NSUjA&t*`(vw1fLKNT_So6&A0B8pLrQ1xeyq zzD=Q)g?5pUA5taDnvdK;`&*QKlz)=yr2H6Rx;)CfLUp(o-YetO{f6prFAUo;oSVCe^=9G@7m!{z+l8&1eT@xjv`qgvB3JR4 z!LeONGB4#>7W}iURfWYL%w1DhGHS0XEWIzUrf|Vs`Spc?;c0b+hjQ%t+!=)dz-kIh zR~43k?5e_o)fbbUKC-nH;!(C#oNEh9?O)~0P>Bwqs4wAo$+P|ko^ZWnyFKfvu&t{- z7y`CcYQX0h@VN$jo&i^**HdR(xO3|z+rnL3FWG8W*AK@XUN6~d$JYGyx5~EKyr7bAH{cxxe1`$=^5CvzuYFh2JN|o=g`2Gc z%~C?0&(%=N+M^H}tfKZ)s+=znxA?qbYAfASXIq^f3;`Ri7(cS%%JU-|H2lbR?*g7G z+w#VVLbhd!bGpZdE3aqG6|=324DvG!xc_a*wk|QqUuwYdZHkv{tH=uiV!JCWPn~VO z&4VFeTbCJdd|~1xTdn5&aP`#bsj#gpJQxDDb)^Bvm5`Tg>o2?@Ahz{(1AdhO=M_w* z`g>5EYt2amJW}(a;D;|#_`I~dS|h?)3*fls&}pXu?^ky4F5?pGK46f?HQ-#a9$@(w z6`olSa0Y*DuyYw)a;(d~kaDG-(ny+bAfB6<{~5murV_vQ))?>>z~_o}4*Tr`e2%#9u>BFl zb2p~rXZ&-@PG()i!ZGvu2L&vGT7S3&rx3pxm{37J#kdw+tdBRxD{U0uwhkALpa<&&=#%_Ctz1Q+!c+*;^95< zM69D6-}r>LHYXb@!{P3DcW-wRPTk$)hG009Ia%7&9v?_%1JOhxWR*A6R^Xo?{#lHF z0%@(v?%rfH(Hn^cm#vTVL^~Q|k@lz)k99;7bxT9L({kEIJM3MLR}K5{S=HltSuDJa zCl?fQ-e{H$_V*1Wdqa(D+hggJZdkp#In)wvsae$!3R_K~z67+kbDTC+^zZF$2R}{W z)v@g7DL^Y(fqTRoSe@IxmQBsc;L^58e|LL0^Twhbz+iheptfFw-SiLFA8hlOT zk1{Y<8Q&dEbjIR)LfsZz_Q^>5u5f$Tu5f3h8^$hO+nN`+;Ih=Vw(Jgsdly%Pee*}1 zVwS;dm^yB}I=t>(J+p0;amS2V27S0P!Daf4x7(`QT%L$_mWRVT1_#4^(L{e7CaCV@ z-tg`KT$<~8}|N4Pg3}hM;(f(*+cXUJBPN>AH zo?W$Yomn^^G`BR>uU``mHMBss*KFFloZhZCE^iKPX>1DBhL&^=^g^Af2Fwk{qP;tk zT~<%Dr$3rhW;cauTQ)Q)vrXZqXe4@TO$cLx_kxxYKy9*gw%hm(8zqT%k|&bVvFzdgQlb&G#Q9It;w zTuZvzJSCgrTjkD;KFVyDyh>%Nf?LbvZbM7E;hlOk)U-C;eK<8uR(ly-6!1zu)QETE z$*V2zk4AUJJ3B)$aI`pF50iQ`ARXNaIA5a9H}pk&d*U6@DVyz!_lKI36`jEAgTZjP zrwwM_WHg)v^a5;G^>@Yh^hA31;<{l9O@({Ljwmizrk0_qfnHdeshh`?+$|lLx=vMJ z8{C@ku7o$Ltt~CydktT-K2>VKW#08(2P@zGaJ=EZu0OwuR>Jt$euu!<34BQ4?E=43 z;D0Uf4+#7jfky@YC4mnM{G7n0ooTS&g@kcRJ68x?+PPNXvRxGd9|BvPhg$@Gr@(I| z93N>92)s?;2L&DzxRk$L;8LF73uC+P5_b5VFyeO${0U(P-&Jw^-y|I6WqW@l@F5}p zTY=vx@XN8GKxpS4IM{xf0beZel#s{I{gKd)#FGLaf;zVUpuna4V*;1*Cj>6#zaelb z{|kXj`9d3l0Mg?uH{k0Hcst?f<^UYrt}Y>eufX>RdE7;0`TYW4Ch%clXGGxl3VG@G zeu2yJ_EBNyK4IrkAusJ55%Mzre-QFA{^tcQ<3B3w$oRi16N7{8 zqYM(_(l5^KNV?w#Tv+8z=JOMT>wX_6T=$FLKOkZJJ#cV5&mst<EzZa4p>G#db59SLA^KZHy@m)XXEr7y6*e~ycgYgi8K)Rg<7gqfu?cmx23H>MG zU^~)2{vS)mCu3jE8z%b~gH0rie*JIEuaq7S&eceIJnnZCZr)@(JRbCTo_5W6;`tKcdOWQHm-E|K4Y9 z=53JT*{R_jIWB;67ZUsBZ)@>c28r=uvcvNU;|B=mamDz3g!BE5@lOyg?{^p{xAzg{ z)$McC<@O#S9P=;7^WPDkwomPPPQ&^6@Jqr@1j=mx%bNUcWanQs{64~esNo#{Pc@wL z`AZGw{J%~(=HR1na6J5J7q<67flnv@++Ma{NI1&(3i(+=UXH_3funDZXFlPG$??BX zldpq%q@a+Owj``mS2gmsl%?{5CpVIJ4KnLk@4WI48)>9h3knk^RILCQH z!#SRpG@SGNGU0eG$KiJg$96pg2iy69up`G6{-yvD%6G%T@~;Vb^uajqhj71eJQpJX z#Q0;RKZ|g*FV_{N2Ko7fqg}a9UL^2QsN;AlG&^-rhjbI+y5D+({CZ(u`rRaO={HZ{ z^15@qyT7AxC7%e6Y56!D#e@|QNBm<{M#lAs^V}<{)pTRB^GP9iA$EAphUTQ5M}&PY zj@PbFv(PwpWF8(i$nOw1Kj)&~E`gU)8S%#izChrABk+L0c`b@irNBQ!Wn|SQazD3}l6L_n@pAa~=5&b?X@IE2WYXgJ^1@tRK~nXoJA1F^%Cb@_80awtjJ=>{PXm&h@2}363aUH@`N(Hnt$RO}afgculwZPvma7W;K1TN#3d64!c&Z-#C zPB=JjnIBxYFfQ#QUJVBl*R!2!g=XjHQNoZ(dAvU|l~Mui>}3#mrNE_~YJu;gGBTN; z-2&$#+W!!r057F-F7SLRBdee?;?t>&tdYuybFL8T68J?_Mm9ucl%GLmWVlvDLi}P? za{oL{WyE<+gwQd8UrJ?Ur>Knbg;YlN29*&nqB1gGZz9fpia5WYLi}wA5Vc~$fDpe- zmE1r49Dw-c0zV?;XA1n7hWAn3DS@LM#@`S)wyTxoOVA+@GHh2X97s0{{B{>oe_I4j zMTKt{_!UC_u)yVd#jI>$h@L8A3ySN13aYV}D`--Ghj-QxVfi;y25VGO}?_|7lwd6J5g)8Mb&CN>Z6g(oizZ;I7t zTN9f_N&Nh1VoChK$ju-3fhiLxX|Sn%KcPOM@o_E1&>DV$<0GbZ9Z1Tqo1BLH25e#% zkgf?7Wa#50SBhow<&I1qay#1F!_mR^XkRi6mjQfVB;tcFEfUFrK7#O@p>QOTi0lP` zhqvK78t?`mu&Uv#t^+864)kS<|Cco+KY>+f>XfXKylz!X%lIaJgvYZtVTjV6S=j$; zDsxu-<<0hezTZUS>Vrw&Zu9w4WWs1JU*r^9eGe?~|Grf&)P z3TW?7R_X&?eFI!R-q(=j?eU(TXfNDgdgIAx`I_|`7pd<_$x^f{+?jySx~=jZi8x$; zZvs6&Eqi^yCabByYa(1z;@T~baONg`Ubpd@2ler}fc2Q0^jpY})W_!>)?;qcZ}ri~ z=TX*UZqkpDKCgK&{;DbH5BccdFa>>HPs;Y=`kLcrZi@df`Iq|m{K9(7P5Mvz=;QM< z>oGU!KkK89&)KZU+@ycZM<1UbS&z9%|D=yTuJKrpxk>*OAN|28=$|HidHsBB3i{`K z{69PeeSTh${&C-n+s)k6{ybW*%khKzA*{#Tq+d+>(m(E(u^w}iK0nt;|A(ib&->)k z|EH&*&->}pKR&N=+{{h!SNh}+-+QngbCZ6ZkN!hb(C6Qem<>8eIDg^33HxPk@*gAr zQXg%w9&?jEe~&2he_#swhkVA*!71oJFvjRD#|P-d(7@6|8fro>)0N1oYSy8P4j<<^f70WbpQCyO66bE zKfKJzKVk5HlKgW&WBQsljNV{2USIS`{#hNU;ga--{^u}gIS`1@FX?GMaw1Aob7>A!mWpD_3@UE!+W z{DOq-M@aYo9|r%ti9eD5e>eE=TZvxb{etj+Bb@90Unl?i{qrSyi6Gw%pgl5|<^^!$ zM-0ql(SMRYe#;<<^%=XG^l7?f*y^K&FaORNSLpeFJLs$N#{$I1u(3b|G>>5ZYw#RF zAf8tjgFe0w0G3535d>m+_bhYH_#%V<9hV3VV#-1o|p}S^`gz|F^Sa;(GhL zN#EA=`6UheVKW>zpLcNln*J$@-yuPae>EJe&vrjzh`*Fx5=|8U3x@csNuT|2yzHOP zpEJb2AmmEO`0IV*|A`^~R`P!(J0h;<|1ubS*nfKd2PYZ-rJxW0tZFU)eLnFwfE>4< z^Kt|E&(-t~O)~yAL;Oc38UJ=e{LlKtzriQ|4;lQAk^dpkk;Lt1>|TR^{t3b^ns9LZ z-U;xDwKTGeE4*{$DfrA0z*b?1(t)GyD^S|C5d@ z5g_#YzwGXMCjuIC>&J`mIA&uaQOfaR$fX?;GQbAkHRppR1o zxnr}|ME*IhCeqjAcP1JCI?(sKeu~$+q40Hr*G~r=?3Z!er@{DnoVi-C?j!#T86Zy6 zq0&JJf47P;DSDfp5KW(DKSBE1?Z+A-{V}bA<7fQ`x%~_fUqJ2W`y#f#3l7%j_G9dt z3uU2tmSR>{vCDS>+X9EKkB=F!1n|(2`7ZS8jO;2sS9qED>36Mj^`3(ENz(QI4=~ox A$p8QV diff --git a/get_tau.cpp b/get_tau.cpp index 1c0e28d..2d12c6d 100644 --- a/get_tau.cpp +++ b/get_tau.cpp @@ -1,6 +1,6 @@ //Credit: Inspired and amended from https://github.com/HopkinsIDD/IDSpatialStats/commit/2782d6dcc9ee4be9855b5e468ce789425b81d49a by @gilesjohnr @jlessler //Author: @t-pollington -//Date: 29 Jan 2019 +//Date: 6 Aug 2019 #include #include #include @@ -10,66 +10,68 @@ using namespace Rcpp; // [[Rcpp::export]] NumericVector getTau(const NumericVector ORIG_ID, const NumericVector x, const NumericVector y, const NumericVector onset, const NumericVector r, const NumericVector r_low, SEXP index){ //see how multiple vector arguments are parsed in here rather than the single posmat matrix + IntegerVector ORIG_IDint = as(ORIG_ID); + IntegerVector onsetint = as(onset); NTYPE i,j,k; - double dist2 = 0; - double r2 = 0; - double r2_low = 0; - long long num_cnt, denom_cnt; //counters for those filling conditions// + float dist2 = 0; + float r2 = 0; + float r2_low = 0; + unsigned long num_cnt, denom_cnt; //counters for those filling conditions// unsigned short serialintvl = 7; //mean serial interval of the disease unsigned short r_size = r.size(); NTYPE N = ORIG_ID.size(); int *inds = INTEGER(index); NumericVector tau(r_size, NULL); - double piInf = 0; + float piInf = 0; bool check = 1; bool bstrapconflict = 0; bool sameperson = 0; bool iscasepair = 0; bool temprelated = 0; bool withindist = 0; - -if(*inds==-1){ //if index was set to -1 then it means we can turn off bootstrapping checks - check = 0; -} - -num_cnt = 0; -denom_cnt = 0; - -//calc piInf -for (i=0;i= r2_low) && (dist2 < r2)); - sameperson = (ORIG_ID[i]==ORIG_ID[j]); - denom_cnt = denom_cnt + (!(bstrapconflict)*!(sameperson)*withindist); - iscasepair = (onset[i]!=-999) && (onset[j]!=-999); - temprelated = (abs(onset[i]-onset[j]) <= serialintvl); - num_cnt = num_cnt + (!(bstrapconflict)*!(sameperson)*iscasepair*temprelated*withindist); + for (j=0; j= r2_low) && (dist2 < r2)); + sameperson = (ORIG_IDint[i]==ORIG_IDint[j]); + denom_cnt = denom_cnt + (!(bstrapconflict)*!(sameperson)*withindist); + iscasepair = (onsetint[i]!=-999) && (onsetint[j]!=-999); + temprelated = (abs(onsetint[i]-onsetint[j]) <= serialintvl); + num_cnt = num_cnt + (!(bstrapconflict)*!(sameperson)*iscasepair*temprelated*withindist); + } } + tau[k] = (float)num_cnt/denom_cnt; // pi(r.min,r.max) + tau[k] = (float)tau[k]/piInf; } - tau[k] = (double)num_cnt/denom_cnt; // pi(r.min,r.max) - tau[k] = (double)tau[k]/piInf; -} -return(tau); + return(tau); } \ No newline at end of file diff --git a/run_get_tau.R b/run_get_tau.R index f5cc411..8bbb278 100644 --- a/run_get_tau.R +++ b/run_get_tau.R @@ -1,10 +1,16 @@ -rm(list = ls()) #clean up from before -library(Rcpp) #sourceCpp() -setwd(FolderWhereYouStoreYourData) -load(file = "YourData.RData") #say you have called your matrix object "X" -setwd(FolderWhereYouStoreYourTauCFile) -sourceCpp("get_tau.cpp") - -r.max = seq(50,500,50) #upper bound of each bin, in metres -r.min = rep.int(0,length(r.max)) #lower bound of each bin -getTau(X[,"ORIG_ID"],X[,"x"],X[,"y"],X[,"onset"],r.max,r.min,as.integer(-1)) #note that as we want the point estimate of tau we set the last argument to -1. You could adapt this function for bootstrapping by parsing in a vector of indices for the last argument. +rm(list = ls()) #clean up from before +library(Rcpp) #sourceCpp() +setwd(FolderWhereYouStoreYourData) +load(file = "YourData.RData") #say you have called your matrix object "X" +setwd(FolderWhereYouStoreYourTauCFile) +sourceCpp("get_tau.cpp") + +r.max = seq(50,500,50) #upper bound of each bin, in metres +r.min = rep.int(0,length(r.max)) #lower bound of each bin +ptm <- proc.time() +for (i in 1:10) { + getTau(X.NW[,"ORIG_IDenum"],X.NW[,"x"],X.NW[,"y"],X.NW[,"KA"],r.max,r.min,as.integer(-1)) #note that as we want the point estimate of tau we set the last argument to -1. You could adapt this function for bootstrapping by parsing in a vector of indices for the last argument. +} +timing = (proc.time() - ptm) +timing/10 #2.0578s for commit 5dceac9 get_tau.cpp +#2.03s for NTYPE change to unsigned short \ No newline at end of file diff --git a/run_get_taupersonal.R b/run_get_taupersonal.R deleted file mode 100644 index b856ec7..0000000 --- a/run_get_taupersonal.R +++ /dev/null @@ -1,42 +0,0 @@ -rm(list = ls()) #clean up from before -library(Rcpp) #sourceCpp() -setwd(FolderWhereYouStoreYourData) -load(file = "YourData.RData") #say you have called your matrix object "X" -setwd(FolderWhereYouStoreYourTauCFile) -sourceCpp("get_tau.cpp") - -r.max = seq(50,500,50) #upper bound of each bin, in metres -r.min = rep.int(0,length(r.max)) #lower bound of each bin -ptm <- proc.time() -for (i in 1:10) { - getTau(X.NW[,"ORIG_IDenum"],X.NW[,"x"],X.NW[,"y"],X.NW[,"KA"],r.max,r.min,as.integer(-1)) #note that as we want the point estimate of tau we set the last argument to -1. You could adapt this function for bootstrapping by parsing in a vector of indices for the last argument. -} -timing = (proc.time() - ptm) -timing/10 #2.0578s for commit 5dceac9 get_tau.cpp -#2.03s for NTYPE change to unsigned short - -ORIG_IDenum <- file("ORIG_IDenum.bin", "wb") -writeBin(X.NW[,"ORIG_IDenum"], ORIG_IDenum,size = integer()) -close(ORIG_IDenum) -x <- file("x.bin", "wb") -writeBin(X.NW[,"x"], x, size = double()) -close(x) -y <- file("y.bin", "wb") -writeBin(X.NW[,"y"], y, size = double()) -close(y) -KA <- file("KA.bin", "wb") -writeBin(X.NW[,"KA"], KA, size = integer()) -close(KA) -r.maxbin <- file("r.max.bin", "wb") -writeBin(r.max, r.maxbin, size = integer()) -close(r.maxbin) -r.minbin <- file("r.min.bin", "wb") -writeBin(r.min, r.minbin, size = integer()) -close(r.minbin) -indx <- file("indx.bin", "wb") -indxvector = as.numeric(sample.int(16221,16221,replace = T)) -writeBin(indxvector, indx, size = integer()) -close(indx) -indx <- file("indx.bin", "rb") -readBin(con = indx,what = 4,n=16221) -close(indx) diff --git a/tau-statistic-speedup.Rproj b/tau-statistic-speedup.Rproj index 8e3c2eb..5bf29f4 100644 --- a/tau-statistic-speedup.Rproj +++ b/tau-statistic-speedup.Rproj @@ -11,3 +11,8 @@ Encoding: UTF-8 RnwWeave: Sweave LaTeX: pdfLaTeX + +BuildType: Package +PackageUseDevtools: Yes +PackagePath: anRpackage +PackageInstallArgs: --no-multiarch --with-keep.source diff --git a/testloadbinary.cpp b/testloadbinary.cpp deleted file mode 100644 index 2e9f890..0000000 --- a/testloadbinary.cpp +++ /dev/null @@ -1,328 +0,0 @@ -#include //printf -#include //cout -#include //streampos -#include -#include -#include -#include -#include -#define NTYPE unsigned short //fine if N<=65535 -using namespace std; -using namespace std::chrono; -int main(){ - ifstream fileorig_idenum("ORIG_IDenum.bin", ios::in|ios::binary|ios::ate); - double* orig_idenum; - int idx=16221; - orig_idenum = (double*)malloc((idx+1) * sizeof(orig_idenum[0])); - if(orig_idenum == NULL) - { - printf("not enough memory"); - return 1; - } - if (fileorig_idenum.is_open()) - { - streampos size = fileorig_idenum.tellg(); - cout << "size=" << size << " bytes\n"; - char * memblock = new char [size]; - fileorig_idenum.seekg (0, ios::beg); - fileorig_idenum.read (memblock, size); - fileorig_idenum.close(); - cout << "orig_idenum is in memory \n"; - orig_idenum = (double*)memblock; //reinterpret - } - else - { - cout << "Unable to open file"; - } - printf("orig_idenum[0] = %f\n", orig_idenum[0]); - printf("orig_idenum[1] = %f\n", orig_idenum[1]); - printf("orig_idenum[2] = %f\n", orig_idenum[2]); - printf("orig_idenum[16220] = %f\n", orig_idenum[16220]); - vector orig_idenumv; - for(int i = 0; i xv; - for(int i = 0; i yv; - for(int i = 0; i KAv; - for(int i = 0; i rmaxv; - for(int i = 0; i rminv; - for(int i = 0; i indxv; - for(int i = 0; i tau(r_size, 0); - double piInf = 0; - bool check = 1; - bool bstrapconflict = 0; - bool sameperson = 0; - bool iscasepair = 0; - bool temprelated = 0; - bool withindist = 0; - - printf("N = %d\n", N); - printf("r_size = %d\n", r_size); - - if(indxv[0]==-1){ //if index was set to -1 then it means we can turn off bootstrapping checks - check = 0; - } - - num_cnt = 0; - denom_cnt = 0; - - //calc piInf - //#pragma omp simd collapse(2) reduction(+:num_cnt, +:denom_cnt) - for (i=0;i= r2_low) && (dist2 < r2)); - sameperson = (orig_idenumv[i]==orig_idenumv[j]); - denom_cnt = denom_cnt + (!(bstrapconflict)*!(sameperson)*withindist); - iscasepair = (KAv[i]!=-999) && (KAv[j]!=-999); - temprelated = (abs(KAv[i]-KAv[j]) <= serialintvl); - num_cnt = num_cnt + (!(bstrapconflict)*!(sameperson)*iscasepair*temprelated*withindist); - } - } - tau[k] = (double)num_cnt/denom_cnt; // pi(r.min,r.max) - tau[k] = (double)tau[k]/piInf; - } - auto stop = high_resolution_clock::now(); - auto duration = duration_cast(stop - start); - cout << duration.count() << endl; - for (k=0;k