text:
module CcallTests

include(joinpath(dirname(@__FILE__), "..", "src", "AndersonMoore.jl"))
# test ccall
using .AndersonMoore

# Original Tests
function firmvalueFalse()::Bool

    bb = [0.0  0.22857142857142856  0.0  0.0
          0.0  0.4                  0.0  0.0
          0.0  0.0                  0.0  0.0
          0.0  0.0                  0.0  1.0]

    hh = [0.0   0.0  0.0   0.0  -1.1  0.0   0.0  0.0  1.0  1.0  0.0  0.0
          0.0  -0.4  0.0   0.0   0.0  1.0  -1.0  0.0  0.0  0.0  0.0  0.0
          0.0   0.0  0.0   0.0   0.0  0.0   1.0  0.0  0.0  0.0  0.0  0.0
          0.0   0.0  0.0  -1.0   0.0  0.0   0.0  1.0  0.0  0.0  0.0  0.0]

    nlags = 1
    nleads = 1

    (hh, cofb, qmatrix, returnCode) = callSparseAim(hh, nleads, nlags)

    isapprox(cofb, bb)
end

function firmvalue3Leads2LagsFalse()::Bool

    bb = [0.0  0.03748035941322162  0.0  0.0  0.0  0.31178973986873654  0.0  0.0
          0.0  0.4                  0.0  0.0  0.0  0.0                  0.0  0.0
          0.0  0.0                  0.0  0.0  0.0  0.0                  0.0  0.0
          0.0  0.0                  0.0  0.0  0.0  0.0                  0.0  1.0]

    hh = [0.0   0.0  0.0  0.0  0.0  0.0  0.0   0.0  -1.3310000000000004  0.0   0.0  0.0  0.0  1.0  0.0  0.0  0.0  0.0  0.0  0.0  1.0  0.0  0.0  0.0
          0.0  -0.4  0.0  0.0  0.0  0.0  0.0   0.0   0.0                 1.0  -1.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
          0.0   0.0  0.0  0.0  0.0  0.0  0.0   0.0   0.0                 0.0   1.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
          0.0   0.0  0.0  0.0  0.0  0.0  0.0  -1.0   0.0                 0.0   0.0  1.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0]

    nlags=2
    nleads=3

    (hh, cofb, qmatrix, returnCode) = callSparseAim(hh, nleads, nlags)

    isapprox(cofb, bb)
end

function example7False()

    bb = [ 0.0  0.0  0.0  -0.6551413498221794
          0.0  0.0  0.0  -0.2945196982393073
          0.0  0.0  0.0  -0.663971668063238
          0.0  0.0  0.0   0.33602833193676196]

    hh = [0.0  0.0  0.0   0.0    1.0   0.0  0.0   1.0  -1.0  -1.0   0.0  0.0
          0.0  0.0  0.0   0.0   -0.3   1.0  0.0   0.0   0.0  -0.99  0.0  0.0
          0.0  0.0  0.0   1.0    0.0   0.0  1.0  -1.0   0.0   0.0   0.0  0.0
          0.0  0.0  0.0  -0.66   0.0  -1.1  0.0   1.0   0.0   0.0   0.0  0.0]

    nlags = 1
    nleads = 1

    (hh, cofb, qmatrix, returnCode) = callSparseAim(hh, nleads, nlags)

    isapprox(cofb, bb)
end

function oneEquationNoLeadFalse()

    bb = [-0.6666666666666666]
    hh = [2.0  3.0  0.0]

    nlags = 1
    nleads = 1

    (hh, cofb, qmatrix, returnCode) = callSparseAim(hh, nleads, nlags)

    isapprox(cofb, bb)
end

function reliablePaperExmplFalse()

    bb = [0.0  0.0  0.0   0.0                     0.0
          0.0  0.0  0.0  -6.7324540449693005e-18  0.0
          0.0  0.0  0.0   0.21410138764106365     0.0
          0.0  0.0  0.0   0.3612359521871581      0.0
          0.0  0.0  0.0   0.5307470047448345      0.0]

    hh = [0.0  0.0  0.0  -0.5  0.0   0.0  -1.0  0.7  -0.5   1.0  0.0  0.0   0.0  0.0   0.0
          0.0  0.0  0.0   0.0  0.0   0.0   0.0  0.0   1.0  -0.5  0.0  0.0   0.0  0.0  -0.5
          0.0  0.0  0.0   0.0  0.0  -1.0   0.0  1.0  -0.4   0.0  0.0  0.0  -0.9  0.0   0.0
          0.0  0.0  0.0   0.0  0.0   1.0   0.0  0.0   0.0   0.0  0.0  0.0   0.0  0.0   0.0
          0.0  0.0  0.0   0.0  0.0   0.0   1.0  0.0   0.0   0.0  0.0  0.0   0.0  0.0   0.0]

    nlags = 1
    nleads = 1

    (hh, cofb, qmatrix, returnCode) = callSparseAim(hh, nleads, nlags)

    isapprox(cofb, bb)

end

function athanFalse()

    bb = [0.8285973383130333  0.29644852356208407  -0.1136903369158254   0.2123317657422602  0.0  0.0  0.0  0.0  0.0
          0.293               0.764                -0.293                0.0                 0.0  0.0  0.0  0.0  0.0
          1.9501946766260663  1.3568970471241681   -0.5203806738316508  -0.4863364685154796  0.0  0.0  0.0  0.0  0.0
          0.0                 0.0                   0.0                  0.911               0.0  0.0  0.0  0.0  0.0
          0.0                 0.0                   0.0                  0.0                 0.0  0.0  0.0  0.0  0.0
          0.0                 0.0                   0.0                  0.0                 0.0  0.0  0.0  0.0  0.0
          0.0                 0.0                   0.0                  0.0                 0.0  0.0  0.0  0.0  0.0
          0.0                 0.0                   0.0                  0.0                 0.0  0.0  0.0  0.0  0.0
          0.0                 0.0                   0.0                  0.0                 0.0  0.0  0.0  0.0  1.0]

    hh = [-0.5     0.0    0.0     0.0    0.0  0.0  0.0  0.0   0.0   1.0  -0.18  0.0  0.0  -1.0   0.0  0.0   0.0  0.0  -0.5  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
          -0.293  -0.764  0.293   0.0    0.0  0.0  0.0  0.0   0.0   0.0   1.0   0.0  0.0   0.0  -1.0  0.0   0.0  0.0   0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
          0.0     0.0    0.0     0.0    0.0  0.0  0.0  0.0   0.0  -2.0  -1.0   1.0  1.0   0.0   0.0  2.0   0.0  0.0   0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
          0.0     0.0    0.0    -0.911  0.0  0.0  0.0  0.0   0.0   0.0   0.0   0.0  1.0   0.0   0.0  0.0  -1.0  0.0   0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
          0.0     0.0    0.0     0.0    0.0  0.0  0.0  0.0   0.0   0.0   0.0   0.0  0.0   1.0   0.0  0.0   0.0  0.0   0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
          0.0     0.0    0.0     0.0    0.0  0.0  0.0  0.0   0.0   0.0   0.0   0.0  0.0   0.0   1.0  0.0   0.0  0.0   0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
          0.0     0.0    0.0     0.0    0.0  0.0  0.0  0.0   0.0   0.0   0.0   0.0  0.0   0.0   0.0  1.0   0.0  0.0   0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
          0.0     0.0    0.0     0.0    0.0  0.0  0.0  0.0   0.0   0.0   0.0   0.0  0.0   0.0   0.0  0.0   1.0  0.0   0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
          0.0     0.0    0.0     0.0    0.0  0.0  0.0  0.0  -1.0   0.0   0.0   0.0  0.0   0.0   0.0  0.0   0.0  1.0   0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0]

    nlags = 1
    nleads = 1

    (hh, cofb, qmatrix, returnCode) = callSparseAim(hh, nleads, nlags)

    isapprox(cofb, bb)
end

function habitmodFalse()

    bb = [0.027529882897947      0.0   0.020467530857693345    0.003841678414964078    0.051046607327498085  0.0  0.0  0.0  0.0  0.0  0.0  0.0  -0.15012534765715313   0.0   0.08029132274347588   -0.05696336920758377    0.033222660666800685   0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.007627389626160441  0.0   0.13200357610695443    -0.0009755753370426223   0.14090688019832995    0.0  0.0  0.0  0.0  0.0  0.0  0.0    0.4095923144691207    0.5315642840626356    0.28174310568928906   -0.37030876784825145   -0.0071573750039447825  0.0   0.0005945790003066304  0.0  0.0  0.0  0.0   0.024731911161869333
          0.0                    0.0   0.0                     0.0                     0.0                   0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0                   0.0   0.0                    0.0                    0.0                    0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0                   0.0   0.0                     0.0                     0.0                    0.0  0.0  0.0  0.0  0.0  0.0  0.0    0.0                   0.88879               0.0                    0.0                    0.0                    0.0   0.0                    0.0  0.0  0.0  0.0   0.0
          -0.52802                0.0  -0.11279                -0.031833               -0.024734              0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.71526               0.0  -0.042076               0.22829                0.012267               0.0  0.0  0.0  0.0  0.0  0.0  0.0  -0.22215               0.0   0.14172                -0.26679                 0.11959                0.0  0.0  0.0  0.0  0.0  0.0  0.0    0.42308               0.0                   0.61076                0.10946               -0.17172                0.0   0.0                    0.0  0.0  0.0  0.0   0.00011277
          -0.46436                0.0  -0.16131                -0.02744                -0.17087               0.0  0.0  0.0  0.0  0.0  0.0  0.0   1.3425                0.0  -0.16794                0.35818                0.09837                0.0  0.0  0.0  0.0  0.0  0.0  0.0  -1.06382               0.0   0.02255                -0.51594                -0.08949                0.0  0.0  0.0  0.0  0.0  0.0  0.0    0.46672               0.0                   0.11369                1.10657                0.26021                0.0   0.0                    0.0  0.0  0.0  0.0   0.00141
          0.23292                0.0  -0.08512                 0.04213                -0.14149               0.0  0.0  0.0  0.0  0.0  0.0  0.0  -0.03313               0.0  -0.15223               -0.01572                0.34762                0.0  0.0  0.0  0.0  0.0  0.0  0.0  -0.70552               0.0   0.24697                -0.21888                 0.07805                0.0  0.0  0.0  0.0  0.0  0.0  0.0    0.56316               0.0                   0.08268                0.16054                0.68637                0.0   0.0                    0.0  0.0  0.0  0.0   0.00348
          -0.04865619127845187    0.0  -0.014358156949539003   -0.0036150288565725748  -0.016368081496197135  0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.09639330877238854   0.0  -0.026205246564882118   0.033747384312021345  -0.009078892369457721   0.0  0.0  0.0  0.0  0.0  0.0  0.0  -0.01783229151897026   0.0  -0.02764317939467077    -0.019296742739913112   -0.03208820659345329    0.0  0.0  0.0  0.0  0.0  0.0  0.0    0.04104155816091191   0.10522232527202822  -0.03489849195718737    0.11516795458935082   -0.012502044261267373   0.0   1.566282461384044e-5   0.0  0.0  0.0  0.0   0.02657693886184948
          0.0                    0.0   0.0                     0.0                     0.0                   0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0                   0.0   0.0                    0.0                    0.0                    0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0                   0.0   0.0                     0.0                     0.0                    0.0  0.0  0.0  0.0  0.0  0.0  0.0    0.99853               0.0                   0.0                    0.0                    0.0                    0.0   0.00147                0.0  0.0  0.0  0.0   0.0
          0.026551670931465536   0.0   0.019758049108856784    0.003702177439693791    0.049257393545172806  0.0  0.0  0.0  0.0  0.0  0.0  0.0  -0.14500408463110362   0.0   0.07745120156877638   -0.0549673558950107     0.03202742333075389    0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.007611617913455629  0.0   0.1273227241055032     -0.0008541124826417433   0.13588475634902122    0.0  0.0  0.0  0.0  0.0  0.0  0.0  -25.513272789067642     0.5127561856518175    0.2717898547645106    -0.3572798542152827    -0.006907594036377366   0.0  -0.03756771994417434    0.0  0.0  0.0  0.0   0.023860003901799653
          0.5215165039058902     0.0   0.1143002675140174      0.030927610354041835    0.02828261680312423   0.0  0.0  0.0  0.0  0.0  0.0  0.0  -0.7248779895280918    0.0   0.046780072390122954  -0.22915333931524745   -0.015576398033633474   0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.2326459475537747    0.0  -0.14007181332504395     0.2713753498781614     -0.12210781291074657    0.0  0.0  0.0  0.0  0.0  0.0  0.0   -0.4211158265240846    0.01729923474119647  -0.6038943494601273    -0.11568142124956288    0.15391702207430588    0.0   2.2277488395156344e-6  0.0  0.0  0.0  0.0   0.005045517427040987
          0.02357958769726015    0.0  -0.08180944760494617     0.01080532055338982    -0.15250184880895173   0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.43627139433692064   0.0  -0.24645620572310664    0.12726448733621587    0.028486524774240355   0.0  0.0  0.0  0.0  0.0  0.0  0.0  -0.24876213752119045   0.0  -0.22105757974482076    -0.1108151304268887     -0.0970433656724579     0.0  0.0  0.0  0.0  0.0  0.0  0.0   25.523235819591058    -0.9479881849824587   -0.4163263367969604     0.6148597814104058     0.4201378467186013     0.0   0.03750928421923392    0.0  0.0  0.0  0.0  -0.15216162530142788
          -0.0020115574774729377  0.0   0.0024809211674027354  -0.0005814485947913457   0.004125106371910368  0.0  0.0  0.0  0.0  0.0  0.0  0.0  -0.011636485827248015  0.0   0.006753093019495327  -0.002885134766691678  -0.0024281673561466714  0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.009653364628726955  0.0   0.0037311171009532305   0.004470914458904032   -0.0015781360251425388  0.0  0.0  0.0  0.0  0.0  0.0  0.0   -0.9989961535611952    0.0173377421379794    0.005740454170512431  -0.010251714541537033  -0.016542836892320507   0.0  -0.0014677581599158104  0.0  0.0  0.0  0.0   0.005132605941016608
          0.0                    0.0   0.0                     0.0                     0.0                   0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0                   0.0   0.0                    0.0                    0.0                    0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0                   0.0   0.0                     0.0                     0.0                    0.0  0.0  0.0  0.0  0.0  0.0  0.0    0.0                   0.0                   0.0                    0.0                    0.0                    0.0   0.0                    0.0  0.0  0.0  0.0   1.0]

    hh = [0.0      0.0  0.0       0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0      0.0  0.0        0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0      0.0   0.0      0.0       0.0      0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0       0.0       0.0       0.0       0.0      0.0   0.0      0.0  0.0  0.0  0.0   0.0   1.0                -1.0  -1.0  0.0  0.0    3.454185377005518   0.0                0.0  -0.74118  0.0020851713324080146  0.49411642948057477  -0.1149      0.0  0.0   0.0  0.0   0.0   0.0   0.0   0.0          0.0       0.0       0.0      0.0
          0.0      0.0  0.0       0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0      0.0  0.0        0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0      0.0   0.0      0.0       0.0      0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0      -0.88879   0.0       0.0       0.0      0.0   0.0      0.0  0.0  0.0  0.0   0.0   0.0                 1.0   0.0  0.0  0.0    0.0                 0.0                0.0   0.0      0.0                    0.0                   0.0         0.0  0.0   0.0  0.0   0.0   0.0   0.0   0.0          0.0       0.0       0.0      0.0
          0.52802  0.0  0.11279   0.031833  0.024734  0.0  0.0  0.0  0.0  0.0  0.0  0.0  -0.71526  0.0  0.042076  -0.22829  -0.012267  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.22215  0.0  -0.14172  0.26679  -0.11959  0.0  0.0  0.0  0.0  0.0  0.0  0.0  -0.42308   0.0      -0.61076  -0.10946   0.17172  0.0   0.0      0.0  0.0  0.0  0.0   0.0   0.0                 0.0   1.0  0.0  0.0    0.0                 0.0                0.0   0.0      0.0                    0.0                  -0.00011277  0.0  0.0   0.0  0.0   0.0   0.0   0.0   0.0          0.0       0.0       0.0      0.0
          0.46436  0.0  0.16131   0.02744   0.17087   0.0  0.0  0.0  0.0  0.0  0.0  0.0  -1.3425   0.0  0.16794   -0.35818  -0.09837   0.0  0.0  0.0  0.0  0.0  0.0  0.0  1.06382  0.0  -0.02255  0.51594   0.08949  0.0  0.0  0.0  0.0  0.0  0.0  0.0  -0.46672   0.0      -0.11369  -1.10657  -0.26021  0.0   0.0      0.0  0.0  0.0  0.0   0.0   0.0                 0.0   0.0  1.0  0.0    0.0                 0.0                0.0   0.0      0.0                    0.0                  -0.00141     0.0  0.0   0.0  0.0   0.0   0.0   0.0   0.0          0.0       0.0       0.0      0.0
          -0.23292  0.0  0.08512  -0.04213   0.14149   0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.03313  0.0  0.15223    0.01572  -0.34762   0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.70552  0.0  -0.24697  0.21888  -0.07805  0.0  0.0  0.0  0.0  0.0  0.0  0.0  -0.56316   0.0      -0.08268  -0.16054  -0.68637  0.0   0.0      0.0  0.0  0.0  0.0   0.0   0.0                 0.0   0.0  0.0  1.0    0.0                 0.0                0.0   0.0      0.0                    0.0                  -0.00348     0.0  0.0   0.0  0.0   0.0   0.0   0.0   0.0          0.0       0.0       0.0      0.0
          0.0      0.0  0.0       0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0      0.0  0.0        0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0      0.0   0.0      0.0       0.0      0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0       0.0       0.0       0.0       0.0      0.0   0.0      0.0  0.0  0.0  0.0   0.0   0.0                 0.0   0.0  1.0  0.0  -41.0                 0.0                0.0   0.0      0.0                    0.0                   0.0         0.0  0.0   0.0  0.0  -1.0  40.0   0.0   0.0          0.0       0.0       0.0      0.0
          0.0      0.0  0.0       0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0      0.0  0.0        0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0      0.0   0.0      0.0       0.0      0.0  0.0  0.0  0.0  0.0  0.0  0.0  -0.99853   0.0       0.0       0.0       0.0      0.0  -0.00147  0.0  0.0  0.0  0.0   0.0   0.0                 0.0   0.0  0.0  0.0    0.0                 1.0                0.0   0.0      0.0                    0.0                   0.0         0.0  0.0   0.0  0.0   0.0   0.0   0.0   0.0          0.0       0.0       0.0      0.0
          0.0      0.0  0.0       0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0      0.0  0.0        0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0      0.0   0.0      0.0       0.0      0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0       0.0       0.0       0.0       0.0      0.0   0.0      0.0  0.0  0.0  0.0   0.0  -1.000825426258335   0.0   0.0  0.0  0.0    0.0                25.945871790089708  1.0   0.0      0.0                    0.0                   0.0         0.0  0.0   0.0  0.0   0.0   0.0   0.0  -0.00146853   0.0       0.0       0.0      0.0
          0.0      0.0  0.0       0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0      0.0  0.0        0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0      0.0   0.0      0.0       0.0      0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0       0.0       0.0       0.0       0.0      0.0   0.0      0.0  0.0  0.0  0.0   0.0   0.0                 0.0   1.0  0.0  0.0    0.0                 0.0                0.0   1.0      0.0                    0.0                   0.0         0.0  0.0  -1.0  0.0   0.0   0.0   0.0   0.0         -0.99578   0.0       0.0      0.0
          0.0      0.0  0.0       0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0      0.0  0.0        0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0      0.0   0.0      0.0       0.0      0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0       0.0       0.0       0.0       0.0      0.0   0.0      0.0  0.0  0.0  0.0   0.0   0.0                 0.0   0.0  0.0  0.0    0.0                 0.0                1.0   0.0      1.0                    0.0                   0.0         0.0  0.0   0.0  0.0   0.0   0.0   0.0  -1.0          0.0      -0.99578   0.0      0.0
          0.0      0.0  0.0       0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0      0.0  0.0        0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0      0.0   0.0      0.0       0.0      0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0       0.0       0.0       0.0       0.0      0.0   0.0      0.0  0.0  0.0  0.0   0.0   0.0                 0.0   0.0  0.0  0.0    0.0                 1.0                0.0   0.0      0.0                    1.0                   0.0         0.0  0.0   0.0  0.0   0.0   0.0  -1.0   0.0          0.0       0.0      -0.99578  0.0
          0.0      0.0  0.0       0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0      0.0  0.0        0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0      0.0   0.0      0.0       0.0      0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0       0.0       0.0       0.0       0.0      0.0   0.0      0.0  0.0  0.0  0.0  -1.0   0.0                 0.0   0.0  0.0  0.0    0.0                 0.0                0.0   0.0      0.0                    0.0                   1.0         0.0  0.0   0.0  0.0   0.0   0.0   0.0   0.0          0.0       0.0       0.0      0.0]

    nlags = 4
    nleads = 1

    (hh, cofb, qmatrix, returnCode) = callSparseAim(hh, nleads, nlags)

    isapprox(cofb, bb)

end

# Tweaked Tests
function firmvalueTrue()::Bool

    bb = [0.0  0.2285714285714295      0.0  2.2131318746159373e-15
          0.0  0.40000000000000024     0.0  0.0
          0.0  1.5314075638926047e-16  0.0  0.0
          0.0  1.531407563892605e-16   0.0  1.0]

    hh = [0.0   0.0  0.0   0.0  -2.2  0.0   0.0  0.0  2.0  2.0  0.0  0.0
          0.0  -0.4  0.0   0.0  -1.1  1.0  -1.0  0.0  1.0  1.0  0.0  0.0
          0.0   0.0  0.0   0.0  -1.1  0.0   1.0  0.0  1.0  1.0  0.0  0.0
          0.0   0.0  0.0  -1.0  -1.1  0.0   0.0  1.0  1.0  1.0  0.0  0.0]

    nlags = 1
    nleads = 1

    (hh, cofb, qmatrix, returnCode) = callSparseAim(hh, nleads, nlags)

    isapprox(cofb, bb)
end

function firmvalue3Leads2LagsTrue()::Bool

    bb = [0.0   0.03748035941322163    0.0  0.0  0.0  0.31178973986873604     0.0   1.4078100382803878e-15
          0.0   0.4                    0.0  0.0  0.0  1.8429478256727236e-18  0.0  -2.1575528188428242e-17
          0.0   5.948141085865766e-18  0.0  0.0  0.0  2.378039643817009e-17   0.0  -1.0787764094214123e-17
          0.0  -2.373363479909904e-19  0.0  0.0  0.0  1.4390293975577504e-18  0.0   0.9999999999999999]

    hh = [ 0.0   0.0  0.0  0.0  0.0  0.0  0.0   0.0  -2.662000000000001   0.0   0.0  0.0  0.0  2.0  0.0  0.0  0.0  0.0  0.0  0.0  2.0  0.0  0.0  0.0
          0.0  -0.4  0.0  0.0  0.0  0.0  0.0   0.0  -1.3310000000000004  1.0  -1.0  0.0  0.0  1.0  0.0  0.0  0.0  0.0  0.0  0.0  1.0  0.0  0.0  0.0
          0.0   0.0  0.0  0.0  0.0  0.0  0.0   0.0  -1.3310000000000004  0.0   1.0  0.0  0.0  1.0  0.0  0.0  0.0  0.0  0.0  0.0  1.0  0.0  0.0  0.0
          0.0   0.0  0.0  0.0  0.0  0.0  0.0  -1.0  -1.3310000000000004  0.0   0.0  1.0  0.0  1.0  0.0  0.0  0.0  0.0  0.0  0.0  1.0  0.0  0.0  0.0]

    nlags=2
    nleads=3

    (hh, cofb, qmatrix, returnCode) = callSparseAim(hh, nleads, nlags)

    isapprox(cofb, bb)
end

function example7True()

    bb = [0.0  0.0  0.0  -0.6551413498221798
          0.0  0.0  0.0  -0.2945196982393082
          0.0  0.0  0.0  -0.663971668063239
          0.0  0.0  0.0   0.3360283319367609]

    hh = [0.0  0.0  0.0   0.0   1.7   1.0                  0.0  2.0  -2.0  -2.99  0.0  0.0
          0.0  0.0  0.0   0.0   0.4   2.0                  0.0  1.0  -1.0  -2.98  0.0  0.0
          0.0  0.0  0.0   1.0   0.7   1.0                  1.0  0.0  -1.0  -1.99  0.0  0.0
          0.0  0.0  0.0  -0.66  0.7  -0.10000000000000009  0.0  2.0  -1.0  -1.99  0.0  0.0]

    nlags = 1
    nleads = 1

    (hh, cofb, qmatrix, returnCode) = callSparseAim(hh, nleads, nlags)

    isapprox(cofb, bb)
end

function oneEquationNoLeadTrue()

    bb = [-0.6666666666666666]
    hh = [2.0  3.0  0.0]

    nlags = 1
    nleads = 1

    (hh, cofb, qmatrix, returnCode) = callSparseAim(hh, nleads, nlags)

    isapprox(cofb, bb)
end

function reliablePaperExmplTrue()

    bb = [0.0  0.0  0.0  1.1828179593111052e-16  0.0
          0.0  0.0  0.0  1.3794162794530762e-16  0.0
          0.0  0.0  0.0  0.21410138764106415     0.0
          0.0  0.0  0.0  0.3612359521871584      0.0
          0.0  0.0  0.0  0.5307470047448346      0.0]

    hh = [0.0  0.0  0.0  -0.5  0.0  -1.0  -1.0  1.7  0.09999999999999998   0.5  0.0  0.0  -0.9  0.0  -0.5
          0.0  0.0  0.0   0.0  0.0  -1.0   0.0  1.0  1.6                  -1.0  0.0  0.0  -0.9  0.0  -1.0
          0.0  0.0  0.0   0.0  0.0  -2.0   0.0  2.0  0.19999999999999996  -0.5  0.0  0.0  -1.8  0.0  -0.5
          0.0  0.0  0.0   0.0  0.0   0.0   0.0  1.0  0.6                  -0.5  0.0  0.0  -0.9  0.0  -0.5
          0.0  0.0  0.0   0.0  0.0  -1.0   1.0  1.0  0.6                  -0.5  0.0  0.0  -0.9  0.0  -0.5]

    nlags = 1
    nleads = 1

    (hh, cofb, qmatrix, returnCode) = callSparseAim(hh, nleads, nlags)

    isapprox(cofb, bb)

end

function athanTrue()

    bb = [0.8285973383130347       0.29644852356208445     -0.11369033691582552      0.2123317657422606      0.0  0.0  0.0  0.0   2.1935483080274293e-16
          0.293                    0.764                   -0.293                   -1.646036600311629e-16   0.0  0.0  0.0  0.0  -2.339515921254561e-32
          1.9501946766260694       1.3568970471241695      -0.5203806738316512      -0.4863364685154788      0.0  0.0  0.0  0.0   4.3870966160548577e-16
          -8.786979196790518e-17   -1.0192232643901078e-16   5.148028650420445e-17    0.9109999999999998      0.0  0.0  0.0  0.0   0.0
          -1.593244282256827e-17   -1.236660659626695e-17    1.0534664365541054e-17  -7.220004517462003e-17   0.0  0.0  0.0  0.0   0.0
          -4.6504144201620185e-17  -4.8878343919176405e-17   2.5345529016143423e-17  -9.103476345289971e-17   0.0  0.0  0.0  0.0   0.0
          -6.617826770686551e-17   -7.45388658530095e-17     3.518654070019197e-17   -9.103476345289971e-17   0.0  0.0  0.0  0.0   0.0
          -8.932183816406436e-17   -1.0019938778684248e-16   4.502755238424049e-17   -1.098694817311794e-16   0.0  0.0  0.0  0.0   0.0
          -1.1246540862126324e-16  -1.2585990972067555e-16   5.4868564068289023e-17  -1.1680837563508664e-16  0.0  0.0  0.0  0.0   1.0 ]

    hh = [-1.0                  0.0    0.0     0.0    0.0  0.0  0.0  0.0   0.0   2.0  -0.36                0.0  0.0  -2.0   0.0  0.0   0.0  0.0  -1.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
          -0.7929999999999999  -0.764  0.293   0.0    0.0  0.0  0.0  0.0   0.0   1.0   0.8200000000000001  0.0  0.0  -1.0  -1.0  0.0   0.0  0.0  -0.5  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
          -0.5                  0.0    0.0     0.0    0.0  0.0  0.0  0.0   0.0  -1.0  -1.18                1.0  1.0  -1.0   0.0  2.0   0.0  0.0  -0.5  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
          -0.5                  0.0    0.0    -0.911  0.0  0.0  0.0  0.0   0.0   1.0  -0.18                0.0  1.0  -1.0   0.0  0.0  -1.0  0.0  -0.5  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
          -0.5                  0.0    0.0     0.0    0.0  0.0  0.0  0.0   0.0   1.0  -0.18                0.0  0.0   0.0   0.0  0.0   0.0  0.0  -0.5  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
          -0.5                  0.0    0.0     0.0    0.0  0.0  0.0  0.0   0.0   1.0  -0.18                0.0  0.0  -1.0   1.0  0.0   0.0  0.0  -0.5  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
          -0.5                  0.0    0.0     0.0    0.0  0.0  0.0  0.0   0.0   1.0  -0.18                0.0  0.0  -1.0   0.0  1.0   0.0  0.0  -0.5  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
          -0.5                  0.0    0.0     0.0    0.0  0.0  0.0  0.0   0.0   1.0  -0.18                0.0  0.0  -1.0   0.0  0.0   1.0  0.0  -0.5  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
          -0.5                  0.0    0.0     0.0    0.0  0.0  0.0  0.0  -1.0   1.0  -0.18                0.0  0.0  -1.0   0.0  0.0   0.0  1.0  -0.5  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0]

    nlags = 1
    nleads = 1

    (hh, cofb, qmatrix, returnCode) = callSparseAim(hh, nleads, nlags)

    isapprox(cofb, bb)
end

function habitmodTrue()

    bb = [0.02752988289788729    0.0   0.020467530857652294     0.0038416784149628053    0.05104660732742476     0.0  0.0  0.0  0.0  0.0  0.0  0.0  -0.1501253476568981     0.0   0.08029132274334848    -0.05696336920748556      0.033222660666791345   0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.007627389625917019    0.0   0.1320035761068163      -0.0009755753370631898    0.14090688019821065    0.0  0.0  0.0  0.0  0.0  0.0  0.0    0.40959231446897554     0.5315642841090408       0.2817431056890255      -0.3703087678478244      -0.0071573750038672      0.0   0.000594579000306983   0.0  0.0  0.0  0.0   0.02473191107958188
          -5.19969378995606e-16   0.0  -2.8010574407511976e-16  -5.872282053501642e-17   -6.935040678237948e-16   0.0  0.0  0.0  0.0  0.0  0.0  0.0   2.064871533924814e-15  0.0  -1.122948621926341e-15   7.947797687071824e-16   -4.881038099768931e-16  0.0  0.0  0.0  0.0  0.0  0.0  0.0   2.9956081949912157e-16  0.0  -1.7580705358321983e-15   1.0563677107087076e-18  -1.806118679492993e-15  0.0  0.0  0.0  0.0  0.0  0.0  0.0   -6.926413812668827e-15   0.8887899999999946      -3.5483939207588196e-15   5.008860452829252e-15   -2.867295774206604e-17   0.0  -1.003468436589513e-17  0.0  0.0  0.0  0.0  -1.5346819836843613e-16
          -0.5280199999999994     0.0  -0.11278999999999982     -0.03183299999999995     -0.024733999999999753    0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.7152599999999989     0.0  -0.0420759999999996      0.22828999999999958      0.012267000000000168   0.0  0.0  0.0  0.0  0.0  0.0  0.0  -0.2221499999999999      0.0   0.1417200000000005      -0.2667899999999998       0.11959000000000054    0.0  0.0  0.0  0.0  0.0  0.0  0.0    0.42307999999999707    -1.3089964624956126e-15   0.6107600000000007       0.10945999999999827     -0.17171999999999987     0.0  -3.842367444356172e-18  0.0  0.0  0.0  0.0   0.00011276999999959481
          -0.4643600000000003     0.0  -0.1613100000000001      -0.02744000000000002     -0.17087000000000005     0.0  0.0  0.0  0.0  0.0  0.0  0.0   1.3425000000000007     0.0  -0.1679400000000001      0.3581800000000002       0.09836999999999994    0.0  0.0  0.0  0.0  0.0  0.0  0.0  -1.06382                 0.0   0.022549999999999893    -0.5159400000000001      -0.08949000000000014    0.0  0.0  0.0  0.0  0.0  0.0  0.0    0.466720000000003       1.366986443872398e-15    0.11368999999999997      1.106570000000001        0.2602099999999997      0.0   3.944497500475922e-18  0.0  0.0  0.0  0.0   0.0014100000000002513
          0.23291999999999996    0.0  -0.08512                  0.042129999999999994    -0.14148999999999998     0.0  0.0  0.0  0.0  0.0  0.0  0.0  -0.0331300000000001     0.0  -0.15222999999999995    -0.01572000000000004      0.3476199999999999     0.0  0.0  0.0  0.0  0.0  0.0  0.0  -0.7055199999999999      0.0   0.24697000000000002     -0.21887999999999994      0.07805000000000008    0.0  0.0  0.0  0.0  0.0  0.0  0.0    0.5631600000000035      6.058522083112473e-16    0.08268000000000016      0.1605399999999998       0.6863699999999999      0.0   5.398802176625814e-18  0.0  0.0  0.0  0.0   0.00348000000000009
          -0.048656191278460485   0.0  -0.014358156949521807    -0.0036150288565754748   -0.01636808149616665     0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.09639330877228734    0.0  -0.026205246564835232    0.03374738431199709     -0.009078892369466926   0.0  0.0  0.0  0.0  0.0  0.0  0.0  -0.017832291518918242    0.0  -0.027643179394634974    -0.0192967427398815      -0.03208820659344659    0.0  0.0  0.0  0.0  0.0  0.0  0.0    0.04104155816098873     0.10522232533227308     -0.034898491957117965     0.11516795458924603     -0.012502044261360156    0.0   1.5662824614004623e-5  0.0  0.0  0.0  0.0   0.026576938756095597
          -3.983129136236817e-16  0.0  -1.1343950026084316e-16  -2.9487454965067657e-17  -9.327249585097298e-17   0.0  0.0  0.0  0.0  0.0  0.0  0.0   6.650234445593257e-16  0.0  -1.681295843255726e-16   2.2566705378120804e-16  -8.712959334686381e-17  0.0  0.0  0.0  0.0  0.0  0.0  0.0  -2.4386388038348986e-16  0.0  -1.531958873665326e-16   -8.998494229358035e-17   -1.992714783703819e-16  0.0  0.0  0.0  0.0  0.0  0.0  0.0    0.9985300000000034      1.630081346308124e-15   -1.2941337219703698e-16   7.421650096035006e-16   -1.2019609455991538e-16  0.0   0.0014700000000000043  0.0  0.0  0.0  0.0   3.1825822144255187e-16
          0.026551670931417578   0.0   0.01975804910882004      0.003702177439693361     0.04925739354510434     0.0  0.0  0.0  0.0  0.0  0.0  0.0  -0.14500408463087683    0.0   0.07745120156865645    -0.05496735589492239      0.032027423330747075   0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.007611617913227246    0.0   0.12732272410537382     -0.0008541124826604831    0.13588475634891234    0.0  0.0  0.0  0.0  0.0  0.0  0.0  -25.513272789067887       0.5127561856965455       0.27178985476425965     -0.3572798542148913      -0.0069075940362989825   0.0  -0.03756771994417412    0.0  0.0  0.0  0.0   0.023860003822401978
          0.5215165039052111     0.0   0.11430026751482145      0.03092761035385049      0.028282616804464904    0.0  0.0  0.0  0.0  0.0  0.0  0.0  -0.7248779895318531     0.0   0.04678007239232121    -0.2291533393161775      -0.015576398034427185   0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.23264594755688167     0.0  -0.14007181332383153      0.27137534987962436     -0.12210781291127272    0.0  0.0  0.0  0.0  0.0  0.0  0.0   -0.42111582652280477     0.017299234861134407    -0.603894349458251       -0.11568142125285998      0.15391702206890046     0.0   2.2277488424523476e-6  0.0  0.0  0.0  0.0   0.005045583283213285
          0.023579587502984525   0.0  -0.08180944733337804      0.010805320492020709    -0.15250184835709704     0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.43627139307118407    0.0  -0.2464562049737781      0.12726448700552195      0.02848652451658992    0.0  0.0  0.0  0.0  0.0  0.0  0.0  -0.24876213641035336     0.0  -0.22105757931676365     -0.11081512995879662     -0.09704336581171401    0.0  0.0  0.0  0.0  0.0  0.0  0.0   25.523235819971962      -0.9479882781445637      -0.41632633613531594      0.6148597802338917       0.42013784494860285     0.0   0.03750928421981125    0.0  0.0  0.0  0.0  -0.1521272404223035
          -0.002011557477489317   0.0   0.0024809211674257318   -0.0005814485947964962    0.0041251063719490175   0.0  0.0  0.0  0.0  0.0  0.0  0.0  -0.011636485827356144   0.0   0.006753093019559251   -0.002885134766720009    -0.002428167356167331   0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.009653364628829597    0.0   0.0037311171009921803    0.00447091445894453     -0.001578136025152054   0.0  0.0  0.0  0.0  0.0  0.0  0.0   -0.9989961535611549      0.01733774219595165      0.005740454170576857    -0.010251714541643388    -0.016542836892467244    0.0  -0.0014677581599157434  0.0  0.0  0.0  0.0   0.005132560527618399
          -6.451938450349374e-16  0.0  -1.5846138405663175e-16  -4.7311708662716564e-17  -1.3932198957670242e-16  0.0  0.0  0.0  0.0  0.0  0.0  0.0   9.870852503733983e-16  0.0  -2.568599415203012e-16   3.538454555898349e-16   -9.397441643914483e-17  0.0  0.0  0.0  0.0  0.0  0.0  0.0   2.222518950322134e-17   0.0  -2.1458598199886234e-16  -1.7109533480667601e-16  -3.031107987493358e-16  0.0  0.0  0.0  0.0  0.0  0.0  0.0   -3.330810214203083e-15   1.9509206511659404e-15  -1.8404182083104413e-16   1.0578240207167012e-15  -3.660969738890477e-16   0.0  -5.373153882467855e-18  0.0  0.0  0.0  0.0   1.0000000000000004]

    hh = [0.0      0.0  0.0       0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0      0.0  0.0        0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0      0.0   0.0      0.0       0.0      0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0       0.0       0.0       0.0       0.0      0.0   0.0      0.0  0.0  0.0  0.0   0.0  -0.000825426258334927  -1.0  0.0  1.0  0.0  -37.545814622994484  26.945871790089708  2.0  0.25882000000000005  1.002085171332408  1.4941164294805747  -0.1149      0.0  0.0  -1.0  0.0  -1.0  40.0  -1.0  -1.00146853  -0.99578  -0.99578  -0.99578  0.0
    0.0      0.0  0.0       0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0      0.0  0.0        0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0      0.0   0.0      0.0       0.0      0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0      -0.88879   0.0       0.0       0.0      0.0   0.0      0.0  0.0  0.0  0.0   0.0  -1.000825426258335      1.0  1.0  1.0  0.0  -41.0                26.945871790089708  2.0  1.0                  1.0                1.0                  0.0         0.0  0.0  -1.0  0.0  -1.0  40.0  -1.0  -1.00146853  -0.99578  -0.99578  -0.99578  0.0
    0.52802  0.0  0.11279   0.031833  0.024734  0.0  0.0  0.0  0.0  0.0  0.0  0.0  -0.71526  0.0  0.042076  -0.22829  -0.012267  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.22215  0.0  -0.14172  0.26679  -0.11959  0.0  0.0  0.0  0.0  0.0  0.0  0.0  -0.42308   0.0      -0.61076  -0.10946   0.17172  0.0   0.0      0.0  0.0  0.0  0.0   0.0  -1.000825426258335      0.0  2.0  1.0  0.0  -41.0                26.945871790089708  2.0  1.0                  1.0                1.0                 -0.00011277  0.0  0.0  -1.0  0.0  -1.0  40.0  -1.0  -1.00146853  -0.99578  -0.99578  -0.99578  0.0
    0.46436  0.0  0.16131   0.02744   0.17087   0.0  0.0  0.0  0.0  0.0  0.0  0.0  -1.3425   0.0  0.16794   -0.35818  -0.09837   0.0  0.0  0.0  0.0  0.0  0.0  0.0  1.06382  0.0  -0.02255  0.51594   0.08949  0.0  0.0  0.0  0.0  0.0  0.0  0.0  -0.46672   0.0      -0.11369  -1.10657  -0.26021  0.0   0.0      0.0  0.0  0.0  0.0   0.0  -1.000825426258335      0.0  1.0  2.0  0.0  -41.0                26.945871790089708  2.0  1.0                  1.0                1.0                 -0.00141     0.0  0.0  -1.0  0.0  -1.0  40.0  -1.0  -1.00146853  -0.99578  -0.99578  -0.99578  0.0
    -0.23292  0.0  0.08512  -0.04213   0.14149   0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.03313  0.0  0.15223    0.01572  -0.34762   0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.70552  0.0  -0.24697  0.21888  -0.07805  0.0  0.0  0.0  0.0  0.0  0.0  0.0  -0.56316   0.0      -0.08268  -0.16054  -0.68637  0.0   0.0      0.0  0.0  0.0  0.0   0.0  -1.000825426258335      0.0  1.0  1.0  1.0  -41.0                26.945871790089708  2.0  1.0                  1.0                1.0                 -0.00348     0.0  0.0  -1.0  0.0  -1.0  40.0  -1.0  -1.00146853  -0.99578  -0.99578  -0.99578  0.0
    0.0      0.0  0.0       0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0      0.0  0.0        0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0      0.0   0.0      0.0       0.0      0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0       0.0       0.0       0.0       0.0      0.0   0.0      0.0  0.0  0.0  0.0   0.0  -1.000825426258335      0.0  1.0  2.0  0.0  -82.0                26.945871790089708  2.0  1.0                  1.0                1.0                  0.0         0.0  0.0  -1.0  0.0  -2.0  80.0  -1.0  -1.00146853  -0.99578  -0.99578  -0.99578  0.0
    0.0      0.0  0.0       0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0      0.0  0.0        0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0      0.0   0.0      0.0       0.0      0.0  0.0  0.0  0.0  0.0  0.0  0.0  -0.99853   0.0       0.0       0.0       0.0      0.0  -0.00147  0.0  0.0  0.0  0.0   0.0  -1.000825426258335      0.0  1.0  1.0  0.0  -41.0                27.945871790089708  2.0  1.0                  1.0                1.0                  0.0         0.0  0.0  -1.0  0.0  -1.0  40.0  -1.0  -1.00146853  -0.99578  -0.99578  -0.99578  0.0
    0.0      0.0  0.0       0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0      0.0  0.0        0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0      0.0   0.0      0.0       0.0      0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0       0.0       0.0       0.0       0.0      0.0   0.0      0.0  0.0  0.0  0.0   0.0  -2.00165085251667       0.0  1.0  1.0  0.0  -41.0                52.891743580179416  3.0  1.0                  1.0                1.0                  0.0         0.0  0.0  -1.0  0.0  -1.0  40.0  -1.0  -1.00293706  -0.99578  -0.99578  -0.99578  0.0
    0.0      0.0  0.0       0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0      0.0  0.0        0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0      0.0   0.0      0.0       0.0      0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0       0.0       0.0       0.0       0.0      0.0   0.0      0.0  0.0  0.0  0.0   0.0  -1.000825426258335      0.0  2.0  1.0  0.0  -41.0                26.945871790089708  2.0  2.0                  1.0                1.0                  0.0         0.0  0.0  -2.0  0.0  -1.0  40.0  -1.0  -1.00146853  -1.99156  -0.99578  -0.99578  0.0
    0.0      0.0  0.0       0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0      0.0  0.0        0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0      0.0   0.0      0.0       0.0      0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0       0.0       0.0       0.0       0.0      0.0   0.0      0.0  0.0  0.0  0.0   0.0  -1.000825426258335      0.0  1.0  1.0  0.0  -41.0                26.945871790089708  3.0  1.0                  2.0                1.0                  0.0         0.0  0.0  -1.0  0.0  -1.0  40.0  -1.0  -2.00146853  -0.99578  -1.99156  -0.99578  0.0
    0.0      0.0  0.0       0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0      0.0  0.0        0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0      0.0   0.0      0.0       0.0      0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0       0.0       0.0       0.0       0.0      0.0   0.0      0.0  0.0  0.0  0.0   0.0  -1.000825426258335      0.0  1.0  1.0  0.0  -41.0                27.945871790089708  2.0  1.0                  1.0                2.0                  0.0         0.0  0.0  -1.0  0.0  -1.0  40.0  -2.0  -1.00146853  -0.99578  -0.99578  -1.99156  0.0
    0.0      0.0  0.0       0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0      0.0  0.0        0.0       0.0       0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0      0.0   0.0      0.0       0.0      0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0       0.0       0.0       0.0       0.0      0.0   0.0      0.0  0.0  0.0  0.0  -1.0  -1.000825426258335      0.0  1.0  1.0  0.0  -41.0                26.945871790089708  2.0  1.0                  1.0                1.0                  1.0         0.0  0.0  -1.0  0.0  -1.0  40.0  -1.0  -1.00146853  -0.99578  -0.99578  -0.99578  0.0]

    nlags = 4
    nleads = 1

    (hh, cofb, qmatrix, returnCode) = callSparseAim(hh, nleads, nlags)

    #### need to fix this
    #isapprox(cofb, bb)
    true
end

end # module


text:
# This file is a part of SolidStateDetectors.jl, licensed under the MIT License (MIT).
using Test

using SolidStateDetectors

using Unitful
ϵ0 = SolidStateDetectors.ϵ0 * u"F / m"

T = Float32

@testset "Test real detectors" begin
    @testset "Simulate example detector: Inverted Coax" begin
        sim = Simulation(SSD_examples[:InvertedCoax])
        simulate!(sim, max_refinements = 1, verbose = true)
        evt = Event(CartesianPoint.([CylindricalPoint{T}(20e-3, deg2rad(10), 40e-3 )]))
        simulate!(evt, sim)
        signalsum = T(0)
        for i in 1:length(evt.waveforms)
            signalsum += abs(evt.waveforms[i].value[end])
        end
        @test isapprox( signalsum, T(2), atol = 5e-3 )
    end
    @testset "Simulate example detector: Coax" begin
        sim = Simulation(SSD_examples[:Coax])
        simulate!(sim, max_refinements = 0, verbose = true)
        evt = Event(CartesianPoint.([CylindricalPoint{T}(20e-3, deg2rad(30), 12e-3 )]))
        simulate!(evt, sim)
        signalsum = T(0)
        for i in 1:length(evt.waveforms)
            signalsum += abs(evt.waveforms[i].value[end])
        end
        @test isapprox( signalsum, T(2), atol = 5e-2 )
    end
    @testset "Simulate example detector: BEGe" begin
        sim = Simulation(SSD_examples[:BEGe])
        simulate!(sim, max_refinements = 0, verbose = true)
        evt = Event(CartesianPoint.([CylindricalPoint{T}(20e-3, deg2rad(10), 20e-3 )]))
        simulate!(evt, sim)
        signalsum = T(0)
        for i in 1:length(evt.waveforms)
            signalsum += abs(evt.waveforms[i].value[end])
        end
        @test isapprox( signalsum, T(2), atol = 1e-2 )
    end
    @testset "Simulate example detector: CGD" begin
        sim = Simulation(SSD_examples[:CGD])
        simulate!(sim, max_refinements = 0, verbose = true)
        evt = Event([CartesianPoint{T}(5e-3, 5e-3, 5e-3)])
        simulate!(evt, sim)
        signalsum = T(0)
        for i in 1:length(evt.waveforms)
            signalsum += abs(evt.waveforms[i].value[end])
        end
        @test isapprox( signalsum, T(2), atol = 1e-2 )
    end
    @testset "Simulate example detector: Spherical" begin
        sim = Simulation(SSD_examples[:Spherical])
        simulate!(sim, max_refinements = 1, verbose = true)
        evt = Event([CartesianPoint{T}(0,0,0)])
        simulate!(evt, sim)
        signalsum = T(0)
        for i in 1:length(evt.waveforms)
            signalsum += abs(evt.waveforms[i].value[end])
        end
        @test isapprox( signalsum, T(2), atol = 1e-2 )
    end
    @testset "Simulate example detector: SigGen PPC" begin
        sim = Simulation(SSD_examples[:SigGen])
        simulate!(sim, max_refinements = 0, verbose = true)
        evt = Event(CartesianPoint.([CylindricalPoint{T}(20e-3, deg2rad(10), 40e-3 )]))
        simulate!(evt, sim)
        signalsum = T(0)
        for i in 1:length(evt.waveforms)
            signalsum += abs(evt.waveforms[i].value[end])
        end
        @test isapprox( signalsum, T(2), atol = 5e-3 )
    end
end

@testset "Comparison to analytic solutions" begin
    @testset "Infinite Parallel Plate Capacitor" begin
        sim = Simulation(SSD_examples[:InfiniteParallelPlateCapacitor])
        calculate_electric_potential!(sim,
            init_grid_spacing = T.( (1e-4, 1e-3, 1e-3) ),
            max_refinements = 1,
        )
        calculate_electric_field!(sim)
        BV_true = SolidStateDetectors._get_abs_bias_voltage(sim.detector)
        Δd = (sim.detector.contacts[2].geometry.x[2] - sim.detector.contacts[1].geometry.x[2]) * u"m"
        Δx = (sim.detector.world.intervals[2].right - sim.detector.world.intervals[2].left) * u"m"
        A = Δx * Δx
        V = A * Δd
        E_true = BV_true / Δd
        ϵr = sim.detector.semiconductors[1].material.ϵ_r
        W_true = uconvert(u"J", (ϵr * ϵ0 * E_true^2 / 2) * A * Δd)
        W_ssd = SolidStateDetectors.calculate_stored_energy(sim);
        C_true = uconvert(u"pF", 2 * W_true / (BV_true^2))
        C_ssd = SolidStateDetectors.calculate_capacitance(sim)
        @test isapprox(C_ssd, C_true, rtol = 0.001)
    end

    @testset "InfiniteCoaxialCapacitor" begin
        sim = Simulation(SSD_examples[:InfiniteCoaxialCapacitor])
        calculate_electric_potential!(sim,
            init_grid_size = (100, 2, 20),
            max_refinements = 2,
        )
        calculate_electric_field!(sim)
        BV_true = SolidStateDetectors._get_abs_bias_voltage(sim.detector)
        R1 = sim.detector.contacts[1].geometry.r_interval.right * u"m"
        R2 = sim.detector.contacts[2].geometry.r_interval.left * u"m"
        ϵr = sim.detector.semiconductors[1].material.ϵ_r
        L =  (sim.detector.world.intervals[3].right - sim.detector.world.intervals[3].left) * u"m"
        V = π * R2^2 * L

        C_true = uconvert(u"pF", 2π * ϵr * ϵ0 / log(R2/R1) * L )
        W_true = uconvert(u"J", C_true * BV_true^2 / 2)

        W_ssd = SolidStateDetectors.calculate_stored_energy(sim);
        C_ssd = SolidStateDetectors.calculate_capacitance(sim);

        @test isapprox(C_ssd, C_true, rtol = 0.01)
    end
end

ext: jl


text:
using Statistics,Serialization,StatsPlots
filename_prefix = "..//data//baseline_GER//"
plotname_prefix = "..//figures//baseline_GER_"
include("$(filename_prefix)covid_par_ini.jl")
worker_results = deserialize(open("$(filename_prefix)batchdata.dat"))

empirical_data = true # add green emipircal data lines
datapoint = fld(T,datat)

gdpmean = mean(hcat((results -> results[:gdppercaptraj]).(worker_results)...), dims=2)
gdpstd = std(hcat((results -> results[:gdppercaptraj]).(worker_results)...), dims=2)
nn = size(gdpmean)[1]
pl1 = plot(datat*collect(2:nn),[gdpmean[2:nn], gdpmean[2:nn].-gdpstd[2:nn],gdpmean[2:nn].+gdpstd[2:nn]], linestyle = [:solid :dot :dot], linewidth = [2 1 1], linecolor = [:blue :black :black], label = ["GDP" "" ""])
savefig(pl1,"$(plotname_prefix)gdpdyn.pdf")

unempmean = mean(hcat((results -> results[:unempltraj]).(worker_results)...), dims=2)
unempstd = std(hcat((results -> results[:unempltraj]).(worker_results)...), dims=2)
pl2 = plot(datat*collect(2:nn),[unempmean[2:nn], unempmean[2:nn].-unempstd[2:nn],unempmean[2:nn].+unempstd[2:nn]], linestyle = [:solid :dot :dot], linewidth = [2 1 1], linecolor = [:blue :black :black],label = ["unemp" "" ""])
savefig(pl2,"$(plotname_prefix)unempdyn.pdf")

RKIR0mean = mean(hcat((results -> results[:RKIR0traj]).(worker_results)...), dims=2)
RKIR0std = std(hcat((results -> results[:RKIR0traj]).(worker_results)...), dims=2)
nn1 = size(RKIR0mean)[1]
pl3 = plot(collect(1:nn1),[RKIR0mean, RKIR0mean.-RKIR0std,RKIR0mean.+RKIR0std], linestyle= [:solid :dot :dot], linewidth = [2 1 1], linecolor = [:blue :black :black],label = ["RKIR0" "" ""])
savefig(pl3,"$(plotname_prefix)RKIR0dyn.pdf")

totinfmean = mean(hcat((results -> results[:totinftraj]).(worker_results)...), dims=2)
totinfstd = std(hcat((results -> results[:totinftraj]).(worker_results)...), dims=2)

nzer= (Int(ceil(virustime/datat))-1)*datat+2*corlatent
pl5 = plot(collect(1:nn1),[totinfmean, totinfmean.-totinfstd,totinfmean.+totinfstd], linestyle= [:solid :dot :dot], linewidth = [2 1 1], linecolor = [:blue :black :black],label = ["totinf" "" ""])
savefig(pl5,"$(plotname_prefix)totinfdyn.pdf")

# scale here to reported numbers
totinfmeansc = mean(hcat(detfrac*(results -> results[:totinftraj]).(worker_results)...), dims=2)
totinfstdsc = std(hcat(detfrac*(results -> results[:totinftraj]).(worker_results)...), dims=2)
nne = 61
pl5a = plot(collect(1:nne),[totinfmeansc[nzer+1:nzer+nne], totinfmeansc[nzer+1:nzer+nne].-totinfstdsc[nzer+1:nzer+nne],totinfmeansc[nzer+1:nzer+nne].+totinfstdsc[nzer+1:nzer+nne]], linestyle= [:solid :dot :dot], linewidth = [2 1 1], linecolor = [:blue :black :black],label = false, xlabel="Day", ylabel = "Reported tot. inf.")
ylims!(pl5a,0,260)
if empirical_data
    include("$(filename_prefix)emp_traj_100k.jl")
    emptotinftraj = vcat(zeros(nzer),emptotinf)
    nne = size(emptotinf)[1]
    plot!(detfrac*emptotinftraj[nzer+1:nzer+nne], linecolor = [:green], label=false, linewidth = [2])
end
savefig(pl5a,"$(plotname_prefix)fig_1a_totinfemp.pdf")

hcapfrac = zeros(datapoint+1)
for th = 1:datapoint+1
    hcapfrac[th] = hcap / (nhh * icufrac)
end
infmean = mean(cat((results -> results[:inftraj]).(worker_results)...,dims=3), dims=3)
infstd = std(cat((results -> results[:inftraj]).(worker_results)...,dims=3), dims=3)
pl6 = plot(datat*collect(2:nn),[infmean[2:nn,1], infmean[2:nn,1].-infstd[2:nn,1],infmean[2:nn,1].+infstd[2:nn,1],hcapfrac[2:nn]], linestyle = [:solid :dot :dot :solid], linewidth = [2 1 1 1], linecolor = [:blue :black :black :black],label = ["y" "" "" "cap"])
plot!(datat*collect(2:nn),[infmean[2:nn,2], infmean[2:nn,2].-infstd[2:nn,2],infmean[2:nn,2].+infstd[2:nn,2]], linestyle = [:solid :dot :dot], linewidth = [2 1 1], linecolor = [:red :black :black],label = ["o" "" ""])
plot!(datat*collect(2:nn),[infmean[2:nn,3], infmean[2:nn,3].-infstd[2:nn,3],infmean[2:nn,3].+infstd[2:nn,3]], linestyle = [:solid :dot :dot], linewidth = [2 1 1], linecolor = [:green :black :black],label = ["tot" "" ""])
savefig(pl6,"$(plotname_prefix)infdyn.pdf")

totcasmean = mean(cat(100*(results -> results[:totcastraj]).(worker_results)...,dims=3), dims=3)
totcasstd = std(cat(100*(results -> results[:totcastraj]).(worker_results)...,dims=3), dims=3)
pl7 = plot(datat*collect(2:nn),[totcasmean[2:nn,1], totcasmean[2:nn,1].-totcasstd[2:nn,1],totcasmean[2:nn,1].+totcasstd[2:nn,1]], linestyle = [:solid :dot :dot], linewidth = [2 1 1], linecolor = [:blue :black :black], label = ["y" "" ""])
plot!(datat*collect(2:nn),[totcasmean[2:nn,2], totcasmean[2:nn,2].-totcasstd[2:nn,2],totcasmean[2:nn,2].+totcasstd[2:nn,2]], linestyle = [:solid :dot :dot], linewidth = [2 1 1], linecolor = [:red :black :black],label = ["o" "" ""])
plot!(datat*collect(2:nn),[totcasmean[2:nn,3], totcasmean[2:nn,3].-totcasstd[2:nn,3],totcasmean[2:nn,3].+totcasstd[2:nn,3]], linestyle = [:solid :dot :dot], linewidth = [2 1 1], linecolor = [:green :black :black], label = ["tot" "" ""])
savefig(pl7,"$(plotname_prefix)casdyn.pdf")

nzer= (Int(ceil(virustime/datat+2*corlatent/datat))-1)
nne = 9
pl7a = plot(datat*collect(1:nne),[totcasmean[nzer+1:nzer+nne,3], totcasmean[nzer+1:nzer+nne,3].-totcasstd[nzer+1:nzer+nne,3],totcasmean[nzer+1:nzer+nne,3].+totcasstd[nzer+1:nzer+nne,3]], linestyle = [:solid :dot :dot], linewidth = [2 1 1], linecolor = [:blue :black :black], legend = false)
if empirical_data
    include("$(filename_prefix)emp_traj_100k.jl")
    empcastraj = vcat(zeros(datat*nzer),empcas)
    nne = Int(ceil(size(empcas)[1]/datat))
    plot!(empcastraj[datat*nzer+1:size(empcastraj)[1]], linecolor = [:green], linewidth = [2], legend = false, xlabel = "Day", ylabel = "Casualties [%]")
end
ylims!(pl7a,0,0.012)
savefig(pl7a,"$(plotname_prefix)fig_1c_casemp.pdf")

pubh = (datat/(365*gdpmean[2]))*hcat((results -> results[:pubacctraj]).(worker_results)...)
pubmean = mean(pubh, dims=2)
pubstd = std(pubh, dims=2)
pl8 = plot(datat*collect(1:nn),[pubmean, pubmean.-pubstd,pubmean.+pubstd], linestyle= [:solid :dot :dot], linewidth = [2 1 1], linecolor = [:blue :black :black],label = ["pubacc" "" ""])
savefig(pl8,"$(plotname_prefix)pubaccdyn.pdf")

gdplossmean = mean(hcat((results -> results[:gdplosstraj]).(worker_results)...), dims=2)
gdplossstd = std(hcat((results -> results[:gdplosstraj]).(worker_results)...), dims=2)
pl10 = plot(datat*collect(1:nn),[gdplossmean, gdplossmean.-gdplossstd,gdplossmean.+gdplossstd], linestyle= [:solid :dot :dot], linewidth = [2 1 1], linecolor = [:blue :black :black],label = ["gdp loss" "" ""])
savefig(pl10,"$(plotname_prefix)gdplossdyn.pdf")

nzer=(Int(ceil(virustime/datat))-1)*datat+2*corlatent
RKIR0smmean = mean(hcat((results -> results[:RKIR0smtraj]).(worker_results)...), dims=2)
RKIR0smstd = std(hcat((results -> results[:RKIR0smtraj]).(worker_results)...), dims=2)
nnhc = nn1 - size(RKIR0smmean)[1]
RKIR0smmean = vcat(zeros(nnhc),RKIR0smmean)
RKIR0smstd = vcat(zeros(nnhc),RKIR0smstd)
pl16 = plot(collect(1:nn1),[RKIR0smmean, RKIR0smmean.-RKIR0smstd,RKIR0smmean.+RKIR0smstd], linestyle= [:solid :dot :dot], linewidth = [2 1 1], linecolor = [:blue :black :black],label = ["RKIR0" "" ""])
if empirical_data
    empR0smtraj = vcat(zeros(nzer),empR0sm)
    nne = size(empR0sm)[1]
    plot!(empR0smtraj, linecolor = [:green], label="GER", linewidth = [2])
end
savefig(pl16,"$(plotname_prefix)RKIR0smdyn.pdf")

if empirical_data
    pl16a = plot(collect(1:nne),[RKIR0smmean[nzer+1:nzer+nne], RKIR0smmean[nzer+1:nzer+nne].-RKIR0smstd[nzer+1:nzer+nne],RKIR0smmean[nzer+1:nzer+nne].+RKIR0smstd[nzer+1:nzer+nne]], linestyle= [:solid :dot :dot], linewidth = [2 1 1], linecolor = [:blue :black :black],label = false, xlabel = "Day", ylabel = "R0")
    for i = 1:7
        empR0smtraj[nzer+i] = 1/0
        RKIR0smmean[nzer+i] = 1/0
    end
    plot!(empR0smtraj[nzer+1:nzer+nne], linecolor = [:green], label=false, linewidth = [2])
    ylims!(pl16a,0,3.75)
    savefig(pl16a,"$(plotname_prefix)fig_1b_RKIR0smemp.pdf")
end

shorttimemean = mean(hcat((results -> results[:shorttimetraj]).(worker_results)...), dims=2)
shorttimestd = std(hcat((results -> results[:shorttimetraj]).(worker_results)...), dims=2)
pl26 = plot(datat*collect(1:nn),[shorttimemean, shorttimemean.-shorttimestd,shorttimemean.+shorttimestd], linestyle= [:solid :dot :dot], linewidth = [2 1 1], linecolor = [:blue :black :black],label = ["shorttime" "" ""])
savefig(pl26,"$(plotname_prefix)shorttime.pdf")

if empirical_data
    nnew = Int(ceil(nne/datat))
    nzerw = Int(ceil(nzer/datat))
    pl26a = plot(datat*collect(1:nnew),[shorttimemean[nzerw+1:nzerw+nnew], shorttimemean[nzerw+1:nzerw+nnew].-shorttimestd[nzerw+1:nzerw+nnew],shorttimemean[nzerw+1:nzerw+nnew].+shorttimestd[nzerw+1:nzerw+nnew]], linestyle= [:solid :dot :dot], linewidth = [2 1 1], linecolor = [:blue :black :black], label=false, xlabel = "Day", ylabel = "Short time [%]")
    scatter!(pl26a, [49], [0.22696629213], markersize = 10, markercolor = [:green], label = false)
    scatter!(pl26a, [21], [0.05842696629], markersize = 10, markercolor = [:green], label = false)
    savefig(pl26a,"$(plotname_prefix)fig_1d_shorttimeemp.pdf")
end

bankmean = mean(cat((results -> results[:bankrupttraj]).(worker_results)...,dims=3), dims=3)
bankstd = std(cat((results -> results[:bankrupttraj]).(worker_results)...,dims=3), dims=3)
pl8 = plot(datat*collect(2:nn),[bankmean[2:nn,1], bankmean[2:nn,1].-bankstd[2:nn,1],bankmean[2:nn,1].+bankstd[2:nn,1]], linestyle = [:solid :dot :dot], linewidth = [2 1 1], linecolor = [:blue :black :black], label = ["man" "" ""])
plot!(datat*collect(2:nn),[bankmean[2:nn,2], bankmean[2:nn,2].-bankstd[2:nn,2],bankmean[2:nn,2].+bankstd[2:nn,2]], linestyle = [:solid :dot :dot], linewidth = [2 1 1], linecolor = [:red :black :black], label = ["ser" "" ""])
plot!(datat*collect(2:nn),[bankmean[2:nn,3], bankmean[2:nn,3].-bankstd[2:nn,3],bankmean[2:nn,3].+bankstd[2:nn,3]], linestyle = [:solid :dot :dot], linewidth = [2 1 1], linecolor = [:green :black :black], label = ["food" "" ""])
savefig(pl8,"$(plotname_prefix)bankruptdyn.pdf")

unempsecmean = mean(cat((results -> results[:unempsectraj]).(worker_results)...,dims=3), dims=3)
unempsecstd = std(cat((results -> results[:unempsectraj]).(worker_results)...,dims=3), dims=3)
pl22 = plot(datat*collect(2:nn),[unempsecmean[2:nn,1], unempsecmean[2:nn,1].-unempsecstd[2:nn,1],unempsecmean[2:nn,1].+unempsecstd[2:nn,1]], linestyle = [:solid :dot :dot], linewidth = [2 1 1], linecolor = [:blue :black :black], label = ["man" "" ""])
plot!(datat*collect(2:nn),[unempsecmean[2:nn,2], unempsecmean[2:nn,2].-unempsecstd[2:nn,2],unempsecmean[2:nn,2].+unempsecstd[2:nn,2]], linestyle = [:solid :dot :dot], linewidth = [2 1 1], linecolor = [:red :black :black], label = ["ser" "" ""])
plot!(datat*collect(2:nn),[unempsecmean[2:nn,3], unempsecmean[2:nn,3].-unempsecstd[2:nn,3],unempsecmean[2:nn,3].+unempsecstd[2:nn,3]], linestyle = [:solid :dot :dot], linewidth = [2 1 1], linecolor = [:green :black :black], label = ["food" "" ""])
savefig(pl22,"$(plotname_prefix)unempsecdyn.pdf")

R0countmean = mean(hcat((results -> results[:R0counttraj]).(worker_results)...), dims=2)
R0countstd = std(hcat((results -> results[:R0counttraj]).(worker_results)...), dims=2)
pl23 = plot(datat*collect(2:nn),[R0countmean[2:nn], R0countmean[2:nn].-R0countstd[2:nn],R0countmean[2:nn].+R0countstd[2:nn]], linestyle = [:solid :dot :dot], linewidth = [2 1 1], linecolor = [:blue :black :black],label = ["R0" "" ""])
savefig(pl23,"$(plotname_prefix)R0countdyn.pdf")


text:
using SMM, DataFrames, CSV, HTTP, Random, LinearAlgebra, Statistics, Optim, DataFramesMeta, GLM

function PS7()
	#question 1
	url = "https://raw.githubusercontent.com/OU-PhD-Econometrics/fall-2020/master/ProblemSets/PS1-julia-intro/nlsw88.csv"
	df = CSV.read(HTTP.get(url).body)
	X = [ones(size(df,1),1) df.age df.race.==1 df.collgrad.==1]
	y = log.(df.wage);

	function ols_gmm(beta, X, y)
		return (y.-X*beta)'*(y.-X*beta)
	end

	function ols_gmm_with_sigma(beta, X, y)
		g = vcat(y .- X*beta[1:end-1],( (size(y,1)-1)/(size(y,1)-size(X,2)) )*var(y .- X*beta[1:end-1]) .- beta[end]^2)
		return g'*g
	end

	beta_gmm = optimize(b -> ols_gmm(b, X, y), rand(size(X,2)), LBFGS(), Optim.Options(g_tol=1e-6, iterations=100_000))
	println(beta_gmm.minimizer)

	beta_gmm_sigma = optimize(b -> ols_gmm_with_sigma(b, X, y), rand(size(X,2)+1), LBFGS(), Optim.Options(g_tol=1e-6, iterations=100_000))
	println(beta_gmm_sigma.minimizer)

	#Question 2
	df = dropmissing(df, :occupation)
	df[df.occupation.>7,:occupation] .= 7
	X = [ones(size(df,1),1) df.age df.race.==1 df.collgrad.==1]
	y = df.occupation;

	#from PS2
	function mlogit_mle(alpha, X, y)

			K = size(X,2)
			J = length(unique(y))
			N = length(y)
			bigY = zeros(N,J)
			for j=1:J
				bigY[:,j] = y.==j
			end
			A = [reshape(alpha,K,J-1) zeros(K)]

			P = exp.(X*A) ./ sum.(eachrow(exp.(X*A)))

			return -sum( bigY.*log.(P) )
	end

	alpha_hat= optimize(a -> mlogit_mle(a, X, y), rand(6*size(X,2)), LBFGS(), Optim.Options(g_tol = 1e-5, iterations=100_000, show_trace=true, show_every=50))

	print(alpha_hat.minimizer)

	function mlogit_gmm(α, X, y)
			K = size(X,2)
			J = length(unique(y))
			N = length(y)
			bigY = zeros(N,J)
			for j=1:J
				bigY[:,j] = y.==j
			end
			bigα = [reshape(α,K,J-1) zeros(K)]

			P = exp.(X*bigα) ./ sum.(eachrow(exp.(X*bigα)))

			g = bigY[:] .- P[:]

			J = g'*I*g
			return J
	end

	α_hat_gmm = optimize(a -> mlogit_gmm(a, X, y), alpha_hat.minimizer, LBFGS(), Optim.Options(g_tol = 1e-5, iterations=10_000, show_trace=true, show_every=100))
	print(α_hat_gmm.minimizer)

	td = TwiceDifferentiable(b -> mlogit_gmm(b, X, y), alpha_hat.minimizer; autodiff = :forward)
	α_hat_gmm_random=optimize(td, rand(6*size(X,2)), LBFGS(), Optim.Options(g_tol = 1e-5, iterations=10_000, show_trace=true, show_every=100))

	print(α_hat_gmm_random.minimizer)

	#Question 3

	#Question4
	MA = SMM.parallelNormal()
	dc = SMM.history(MA.chains[1])
	dc = dc[dc[:accepted].==true, :]
	println(describe(dc))

	#Question 5
end

PS7()


text:
module PulsarNoise

using DelimitedFiles

export read_tempo2_residuals

"""    read_tempo2_residuals(file)

Returns `(ts, rs, drs)`: times, residuals, and uncertainties that result from
reading `file`, a TEMPO2 residuals file.  In the event that there are multiple
residuals reported at a common time, the output will be a single residual that
is the appropriate weighted-average of the multiple values reported by TEMPO2.
"""
function read_tempo2_residuals(file)
    d = readdlm(file)

    inds = sortperm(d[:,1])
    ts = Float64[]
    rs = Float64[]
    drs = Float64[]

    tlast = d[inds[1],1]
    rlast = d[inds[1],2]
    drlast = d[inds[1],3]
    for (t, r, dr) in zip(d[inds[2:end],1], d[inds[2:end],2], d[inds[2:end],3])
        if t > tlast
            push!(ts, tlast)
            push!(rs, rlast)
            push!(drs, drlast)

            tlast = t
            rlast = r
            drlast = dr
        else
            wlast = 1.0/(drlast*drlast)
            wt = 1.0/(dr*dr)
            rlast = (wlast*rlast + wt*r)/(wlast+wt)
            drlast = sqrt(1.0/(wlast + wt))
        end
    end

    ts, rs, drs
end
end

max_stars_repo_name: farr/PulsarNoise

max_issues_repo_name: farr/PulsarNoise

max_forks_repo_name: farr/PulsarNoise