From 92a77a942133d35ae24541c781725111025e7905 Mon Sep 17 00:00:00 2001
From: SebastianBruijns <>
Date: Mon, 27 Feb 2023 11:47:09 +0100
Subject: [PATCH] new code

---
 .../mcmc_chain_analysis.cpython-34.pyc        | Bin 0 -> 6282 bytes
 __pycache__/pmf_analysis.cpython-37.pyc       | Bin 0 -> 7266 bytes
 __pycache__/simplex_plot.cpython-37.pyc       | Bin 2989 -> 3017 bytes
 all_first_pmfs_typeless                       | Bin 0 -> 45403 bytes
 behaviour_overview.py                         |  64 +-
 canonical_infos.json                          |   2 +-
 dyn_glm_chain_analysis.py                     | 615 +++++++++++-------
 index_mice.py                                 |   3 +-
 pmf_analysis.py                               | 343 ++++++++++
 simplex_plot.py                               |   2 +-
 10 files changed, 754 insertions(+), 275 deletions(-)
 create mode 100644 __pycache__/mcmc_chain_analysis.cpython-34.pyc
 create mode 100644 __pycache__/pmf_analysis.cpython-37.pyc
 create mode 100644 all_first_pmfs_typeless
 create mode 100644 pmf_analysis.py

diff --git a/__pycache__/mcmc_chain_analysis.cpython-34.pyc b/__pycache__/mcmc_chain_analysis.cpython-34.pyc
new file mode 100644
index 0000000000000000000000000000000000000000..1f4261481b77817f0a5299f8eaec91dab934679b
GIT binary patch
literal 6282
zcmcgwTW=&s6|U}{i#;Ci#cS`)o83+}kWAv(1Q9_lWJ5OE08%Vv7Zk}#(rS0tc-q^y
zRP}h*6Z>hyMLhFB^2Qr^LE;}E#1G&Jp?V_X0r>^qc?sV+)ziLg2tfp9rfaI}Ugv!0
zJE!Vz{rTQgfB$Rr#|uLIMI8CGkbeb_{xgzThz+EUh+UC7V#DRUCsI#jEwRxOSzBzh
zg;Tfah>Z@nXo<Kj;*R*h6&qa<;}z+Q_y9x9h<H}SJ@EmCMrKaLea`npJTKx?oJZkl
z5uf3FU&Locyuf)hI49yo&Yu$Tl8Dc99*r-E_#)@ei1?C-FLNFpmqmPq^9v$=M8sD`
ze2tfZaa?g{;7@jYK0dd`c~m7uuKlX;<^4*9QRV9@tP))%k@knG$o$aPvhs_eAC1B!
z*AI7$imV)0(!c%gD6CHOX$JFmWK<PO_YNRnu|9Ckj0*F+aae_>l@}_z6}eP{k9XwJ
z&*Ejz5qeqtM2I0&!dm&U5c{sk7p8?C?aLPK(FR+S-`tkfSYc<OKdg7@SF}D=msca{
zU#D$VDhX4x)4GU~u#zzegKG0?<foz9lG+MI`bkdG2(vPkECz@$3Uhx`Vh)8a%Gafg
zl3`*8De`pJPmAKNk9q5wsP3;<Fh`*!tswNXB9@Uuo5zR#+j!35;r&-aY&&8PDUj%i
zp9*mi+l3?^xDK$2-;eOx;)<awF5MZl*AHkE(+ck=`ruQ1R0r<^&$xMMX3J1j2`S#p
zW82AA9z?}Buk^q%&V5SqkHtS;T)VrwxfX7f@tWRL<77M6Yso0ff=xLJcLLI^)N5Ik
zMFArt2=g%A)rnp$cg=YSFQ|oQn+BwYQ?w=4vom>?M+kH>k-<o&rBtggr9}i?yuA9^
zy&Htk2#tiE!^=?#0WFulsFYL-ofQqm@~MTSx>PVj$t;qAt9tk|UXsTHN6k@rpOSe>
zXiI1D)Z5~kfyb@PJn1!<YMz_g`g#-^z?^+Z4m6vl>4n#Emb?X5(FD!zzZCoU-50hG
zt7!k;D2YZsG%tZxC8uTsn<Art_9B^A;T9@JO6pOO!jArr21AfwgLwoQG8LOTVj1Rg
zx5cIdgK<C(?mDX-N+Ru=ejX4Dl|=#zFrFS~P%^22452#Iv^ixJm5(#V@?|v9c*Nw%
z!-|tYKCB`5H!7snHY<J<FZGH6UCRCYRjUPbQncY3o$55@9kod1hv;0w>lHk_&?z>#
z%^7M))l)=rDI@aC;i~O^&eI4DAl(H=QOPmFBD%Rr%}8=~H`vCVfkRY32hWo0AckAX
zO?xZMvapu1U23yl)8k5(uTP`l5$Y0>FFcuoYZs;+rse2gPQsNu$#pcglPsT{gv^Y0
zL_(_d!@8NXgbSMG&GoE92+~xW**vgspi;)r$^CWHl;{0wGw;=-^OEq0jd<Qqkg@Z=
z?TUSeJYMAxkWpye-~-3m_ojs&p?kk2-~d5^a1a-lsktZi@N2<kc;XydEbh0t2MRGf
z94Eq*R_#E@d2qXCcnroP?RTcrba<LBO+z*KTi)x6#l1F_EMM545qL!$bk`N?_Y@7?
zp;cg6T(L}RoE3Zc%?M0)XV9Je>ZT1KbR9tyAQb$R3U|pn4rK^WF5xhTMVc1(l6=d5
zU8Wh@l&k*Nzx9pb)E8Syt4B#HEz{u&hKl{2P*H61ZkdiTjAgP!`?-`cJoT_p?C<4p
zh!`aS!C4*p<eK2%(`cL`Br0<2D9B|Tknyf70zK@Hoe1aKrmaUI$dEjUX=O6Z*+@-G
zzo)7}&va}MWt_}<h}lvnL0nY9wwaq|@oJn`&72Ol1O47OR5CWsTc(d)1lUDDyD)9e
zpM=*L+k#~WpPrQitV5Hr3h%VD<eYY{JI^7#?DXB1bIwsKsEG)1^da<7n84`!2$|Zg
zPNOb@8`KP4a=Zn6aYzERLHmFgeiVp!q*&mnCC)WCS|k9m^$>0to-`U@1o#n*S~ZM1
zRKhT#<^btY*i@K$){=O8^5zfov4;9~q|#6h`5LQPU1a2ifD`TClW<jafExSZNQN<E
zstE}TC>t8Dl1#21>UUjT!}bSl8wos$JgZSXhMe)T@cy8uu2bn(D0z~S$0_01fNWbm
zK?zA;Jw?gWlrYU5a^0fb3DY2iH?07)IAGnqfQFjd2zSYO9xA&`YFfp6!+<nNqoz;s
zQj>&Nlf)q<9Oa4D)PNyy$lrq)S(^|s>#Si5$Qasuh!tF`P436qwA_N(xoN{8O=jp3
z2<mPwO-HctKub`2O=Mt-W-U;(t;gRo)V{a}EFgOUd{I*fvj$U0Pwe40Ybm5XdAg>M
zDT56CY@AjJ;`D|-C<ZU`XsTt>FB|SUwCvLh2Yts?J(y!H8z(rTURL>`j#>Giujr+d
zRvFOhUc;oC+_->4UBdTy=jExszlPFA-=yotDf|L2wISIIy92Q5MPgP*RLtwh6Xn8!
z+DzB*94P1bhssfF%i%maV%w5=8&-%8IH`8T1r{<OV{4(8E^5&YwKH|CWdf#ssA<>&
zm(NZ<oH8?!pnn@Dyl3ce7pBRCDU?Y(hG|-2=6JuZdcpq=j5sS*u>-5VO;pF!j1!8)
zAqA`-6*-K3OG1le#Bbbu<K~q5=`2Z2?j<<1!JrjIqEu+vJ(2prgX9h<h+(xK5NU7W
z!RJpm>dYUenIk2GSxe?#xT#GilW_u0C;B&CigS6NcvCdI@yc=tu7^1fIt-nSc{nXN
z8GA&+0JJ9X`LwuxJx%{XEh&(iIbU}>QGSv_E^&R5@6>$*4z%FVvNy3D%Ol<u=c04b
zz3N<b)VEQ33|9@0+(r$g00iB?LBM(Koo{l$0sAJ~6>ql$cm?4FJODWXU;$tQ`0xX?
zTH*mQgK*xP7jHFKLxeQ><<vBBJV-#s^<@HI7CeC~1~}goap4h9Z(TneHMwqkT=zk_
zNjAo1M5GFhUUf{cI@NzY36Ak;T9ZKS7ib*q7iEdVHXrSobv1yNN=vze7R$wS!m^Zk
zY|fCsU=vMj1N97qQ!5l2BxJ$?vnN+CqDs9)$r>dmD89#&agTRU&A~Ap>gJ)6M_lzP
z${KYLfCyJNso@nQFmot^WPc?5vqnsU{q0n++he+*;0@WW{}05<R~d9EQqNQEmTZHp
z#uYy=pvaW$+J`;~lqh1eC&OKV`q(;?3cQR{<S^8gREV62=d5%+>rY>yg3prV*@x6N
zLQ_pvEIb@YNNN#@kUm7^RKudYh0HNgGK{H&`H+MQNQ3JdsH3FAB?xSQZJ-@NuSbWO
z$#<thAC3XE<5*ttw-7w$97_U9swi=@R7u<cHC71(+5m7Bn_t1g?B=_6^9`0vFAcLA
zoh%Krbb<&hFwK?ErimK1$HAy>V7>y6J;}l8(3eY1w^VX733Veul3wF;eumf<(y;~)
zhfQG_L<3&nGXPS66zl+D!Mju(3!(&(0IJCAH>LJ|6rlj>;8F<QuH(SVJixWre>4R;
zJTyLEzCOUc?0+{nU)nyjelGp8){ntBjy=8wF8$moVgVHP)+bMKcvIiIQ2?vga1mML
z2~H~k+!Dfz;Qf94R!ifmP4;-V6dn10VETM1g}baIl6vc-59#yqjji_=fBpI=lb_!h
z%$QD@M0ZnZT4@o+rd5CyO^2LHSef}cmdWD)UYReZOdCMYE4~u|mwpekn64ggZsRnC
zv7|~0jDpL;`Z9%KSzmqeT`FG?+AHUKgijxda36QhJC8YC=WztU*Z6<>PXFo-+6)#=
zH^A{J3Ifx{y-c}lef6A<l5%&It}!(GYWB|R`{=?e=B4trf}&%;qBFQPM+&qpE0L>Y
zf@W%AGvPUpQ2gOWW?q6&$Ekdo?1tChar#()-&^RP?)5$Xp6&Pg^ZkB*uHWs?{1dFA
Bh1viB

literal 0
HcmV?d00001

diff --git a/__pycache__/pmf_analysis.cpython-37.pyc b/__pycache__/pmf_analysis.cpython-37.pyc
new file mode 100644
index 0000000000000000000000000000000000000000..5d73cb225dcfed4bac0c68eecdffceec3453a4d1
GIT binary patch
literal 7266
zcmd5>O>A4o5q^K9C`uG%nUN^lap-@^Ka#B|j(_CHwrt0aW5<czSb2?0liC$ciNE_O
zn_}D|sC*~_I3RjTgP@>*AV4o^QM88uMUq1gZIJ>-4{h;K6bV`sWza)UJ)}k2Z+0Jl
zq{vNxB!Efq?au7%%<O!#yYr-zp-@o6=jZ$XJ#%JU)BZw4;}b&XIDXE@>okpNOwVYx
zYW1vPo4V>-wvQRiWES(iVcLG~t939x>tF%a$%2^G#kyID^{{ZNYgxAgte363quHH2
z$ol4u>N@c^uGPADH(SplO9muD9{rF<_oOsCT)x9Yr?hLg?Oy0#$NC}HSH8kSymwFo
z|Llu~z21&s>>(cE>jyQA_w#->FmH(ucpuYz+Mve5k87qjp%vf2NOWEkZ(t<6a^$MU
z!?gkEiSj<)R~uZ`HbU33u4~0R7&B`dc$g@@v9@VhLu=L^_S(&CP>z_jM~K5TJBHsD
z{I=q^t!k9}c+V-VYFrE05w-yyZ?|{IcBj2dwnO%C`CmN3!g5U?+q=;xzc&tQRb9No
zcT>;qhqgWCOMDNl$KKPB;d`Qxg<Q488E9}YQ;m}jzK8e&8tdN(ZI>{oSjF1*^6k8T
z310Ci-ywTDz}qQf>UlzOg~U!s>{^4r3-Uv-0PT4CrHyb@!;B%!7+!N^7}g80zMZhD
z8u=XGLD<;)_)d@WGT$XRBYcQ$k}}`t!|Y+mjPfmfA7ZeDZ{;IsxAJX#6zw*?uWC5&
z@)7wQMNH1}n2c7eJl>4&5aMgr_G29fRK%8b<+-A;Rzn9`#d{L>F}At(DBs6N_~^3f
z;XID)4wDo&=N@Al*&{4=S5XI%?@g+<nQg(`n93vAe4O$F>plz3hj?s>^63!Y;<a17
zb{n<G{ln#_WbU7|d&@_f89(A>d?$}A>8u|+$rP{h!^qtTMP98T3cO&+EWUb`o-fl%
z<~H!dY%ANQJgtqf?d@3$dq<F4wAaqU-tnfr@fG&8c@y?N<wsy|ud>C$$}-{YS^ZR6
z{ikR3$O@}F+O4vk(0<gOfOp6EQ9i+Tv7x&rcu&}q<)3nRv82PJx4=1p6W|0t%K9T3
z)@+H_FKK*|pJ2lgjZbXQ8ohpP))>_0PV!#1o9&^smRWNx)81Md^JEMM+UGYpo*tDY
z>%9z-Ck*#<DTV#H_b$%Ds#W|A)@v2NdsQp{0zOYlpMN5Ko+`fy&XnZ52F@w_G(SCe
zhV5e`Dfn4?l8wrEoUWa1<`k_Gy9bf>Ip5=_nCyuk^3%0*DTyn3PpMw+INNWZhsJZr
z=BL?#Wt|83c{au#y=x)2PR~6Bt#k&D;}rI__aU`Mdherbton)v39Ys>?XfS6(_fo7
z{c#meTGcPkRQ<)>xy=&0gTRh>$p(^9%|4McfMg3<zbHL0$oh+KA~x{<D@V4&NqWe3
zi1OEjY)5e`yzeOXu0r<lFNv%|F(6{qK+*D$JtUD0h=~TyE0CobMl~QVJ}}4zir+zO
z;Qv>SY^O8fA=@dQ`<jpq6#pZ=|9BO$hrcAU3dK(GNwu^1XAjvU64{{Wm&ih51+p~5
zs0PK*1A}a^IE~oA|F0a`E~no^woB~&nve|^e-H12#b2&M_Gmk@UY$JQ)yWq;`j|&w
z#BKVFyh*7B&{6pu=a={eesRg5d#I{ki_tlaU&MYEXZ~~SAbSG$XLvwc#(f6P7=M=R
zp$cf!q~laG&HbO6N#2UPmvr2UF?p|a)82xvR;|>ioTVi;iTkZv+=j&EWvzL$DXMY1
zG1$Dqj<ZR1lY{yOJJF(O=G-~lY8%<fl!<o_RAp$V)Xk<e-~R?|$|-+NA&<MnJygYA
zUgg)S?@zZ!pPjMKvnTE6y{Mn{qJEW~LqzXW%fK@&WvFM$u&3GioWY*9U*IpWr(riF
z1iQdP?82ay#(PlthOEP`qvDb8Ikt_;%$Ai|1T_z;4|-FeSg(PI9D@XEaa4n{R%u9k
z^-8Gt3UVw|ECZiOl549yCP{V$-M@;4{C%noIrqL}&~x80()qcS3VelcUr3H@eGz?R
z_p`N_va=SS^sdFH7=9Dzo|NwdDu(Se#ju?wyC>NR{=Aph5{0=J+3_al5TBNB18o%<
zf1y1Kz03)VxA?Us_0~C!{c((i*>JOp{+2yeuE|__4K;nWSqI-irqWxptYcXiJ25Q&
zQ4JS=hbl6zatF`mdmO6YP4jxSS8U+Z@@;CG_1o8TQJzq;tv7&qYe~1S^C^BED_o;?
zW*5C2ge~M?Se|oJt8!4$t$0a#jkqiR=i{v8)~=A%SUVO|t9F2*+r6P3;+;@jj;`q&
z<?8)x(pWp%Q;agC^#-8mcAGTbvblC1%31fZxkf&*XWH|L>6f+POLh7S&_&mm5}fB;
zE16l~t|>ShHieFnx_lVa`!7BzlGd~j-}5^68lCswOHCjC*{^QC`{QP3(vx|=(Ye>?
zd^jn5&@w~1=|e3srtouK#^W~r05wQYyf8&IX5tSxQD$7$N*<3{N`hWMnQz{5&a9f@
z$BeR?PvfH0nxW{D&N^zwZM<C3rX$_V!Y_cIhOl&A6dmwGh9LF&59sr2<#g)pCqJr0
z-u}a+$Np%tG8}Jxu890XGQ*vGmd7q#ICCWS(zjcKap8y3Vbe9UiL&b}h;**x`e*YZ
zn<x$Iu35;={zwzukojbMbbc{8nn)GcsFM^6>A9RUn!b_E#*_R;;%3|_B}&{G#e{ew
zm&hzSX=kLc=mv;US}bt-L!Z-&ho(n#%NWr&={>rIp6CJ34aDQwL^>Cby90?#CO(@M
zPAQJLj+_kl9oG=aN{|MLJyMXKH<z6H?FaPv(`04W$Y_@HQaYD{mehhEr$oMx7o~JQ
z7tfb&a1l?Y6V9+9V<~a36V{%tcN%eg=e^1n5~PUjjg@$g@!ZU!qB4ns!>Q|9NYh1x
z$~TkG<V7X+S?i=DK64|HOY!(i->whA6y&XBI_FeER}h-mSWM8UQ`vkroq~TcsUFM{
zkL71$6u#J&W3l~kxev4J5;LPt(X4kWbsbo%ldkZ(Nlg63RN5d8Z0R$uPe60YHB0%z
znL6QF@1n;(5;VlaXxxCv149JK+a-UWtaK)F>1-n7#17Ogh_vW2F|dMuv>sSNk5|HS
z#e{-_s}gBv#B+%(2Z(jb>3TmNcIQ++R}u-Qgs5B^lJ4G02<8meDRo4Y1UKVR2{jbO
zNOPUir@ny(x>UaKY+<B;MG>ZABOX`acP1ml=v2DSU*3OUToJBM60V-CL|A%umJ6OM
zr77w%<YU80*I!7_%xAc}4jZD7$(J(eB+{CC7*efS$dp`PZXt_wH*y6xNR!88I=H@s
z5Q#;X?8j%YKBt-O6aZwqmB}ZVYvl_(CkRhB(8zFm-7LR_TvZz@?xY!x!Fj&u`YK%H
z9Xq73k-5wrw}<iBbdIxlA&))m2!hwOvOGaYC^FFXCkimau+!5y*K`)Lu5}~rASIlf
z<90Ca6foWGE>WJuGl|9gLdorL5;vhX<yy{-{4G0FLWa;n7IGyn>~59W9@X2}A<glS
zOyef$!x0lt+CJDv(@j)k3Mep}R1b}*$pnktL!&L))eCkf<2Qj6vQ<##HGGk~jM;6z
zc$|hWigYcBpzjto8+U!BbSc9HrI~9jX3|-=BT>GYp>UvcLw0;KnY_aV?KV4zJd3Lx
zYx|{Jh;A+}^0#CQ<Q+HQEF{T><LYtOOq3nBlOshq?P*s}yFrJS;)|I?l4o24N3TBT
z>Iv5`L2%qI8oilz(m1%GN{YfQIjLMiMmnAX*|$5XL55T~P?nCl>tx2!Y8zY7?O)AO
zn-IsnJSJVUk}kLw<t1=e-n(XM240nwpMKBZpd^*O9#Yymt=*1Bas(-L;xe+XURH}J
zQ*6~Dlt!z#kVXblhM9?^WBV5z9!E4xq_Dho+8*@(AU80P&9j9Je}XP(XEUmRfZnYK
z%x*nq20}r-SMN0}J&11~u9pCw0bDnmKnL)}6Qk-&qx9{Oa+YqOza5e!71X2nZqQ?x
zPhTZTZI3<#ecf_aube5zq)oE8)!*Y;Q63C|M^?!T8aMQvXl<#L^h7b*joDFszqH(v
z^h3JGFdn2iCPf#KhWA1EH4e{J{CW+FVoV=GODn*NO^R|a))bLXS|wo+)Grc7_->Nl
wA;fb?Mmy^95pQUaIw0+bo(Nison^(q_Zv|Z52PWAnVTR@{T@&&`%1_E09gD_`v3p{

literal 0
HcmV?d00001

diff --git a/__pycache__/simplex_plot.cpython-37.pyc b/__pycache__/simplex_plot.cpython-37.pyc
index 8fa16cdde21c4a7d3d6daa3797043fc2e294e842..cc971e6effeb8c5cd96ac67921d02e35f9142f1b 100644
GIT binary patch
delta 444
zcmYjNF-yZh7`;oHv}p)dL24b0SWyC{YAxvG(0?F0bqQ&5wrSHOd^v3^f<Hh6A^w0m
zIygud|AIKW-7l!SgW%#_JJ=h(_wIW<-n)A#z88&#VdN1%gNtAH=heUx{+=vJpLoaH
z@XY%RMmL!U3{#56YflKHFvnR*Fb|(0Yjcvx4DqgC!3dpOGSgeFXlVPP+aU>)z1U|u
z1t&~m2S`ZMJLRj=%qPwjv8hi30u`Rj)ZV~>cEY6KpsB?Wb%{%f_l%&Hj@3bK1*<qO
zL49Ior+B|49LhC;C}B#{abm(0p3o++E)?12dx_`Di62@ysHR6+xlrR0%*LT$<+Q6=
zxST#}Z?~d>WwNYody#9~ObIXh@d#!y%YBE&K|}+;RgFi2F?EixD#8|Z>bpGxOr0Z`
zTyGW@{DB|xFHvRyMd3DtnttR@)jGfOQd_#UeJm6?AR#owa!-u^_VXvCk$&<Ar^be3

delta 289
zcmX>pzE+&qiI<m)0SGeAJW9U1k(ZH)@xx|*rpt_6y-bV@H7pYti=IsWz%0%vx0!>b
zf{|a7@fLSdQhr5zW?pheYVqXdtS|X)aTO<)rN)<JmgJ-sX-%HTCd+6zc@JAT2OlFF
zBNtPV@n&oG7mNntKwECHR;3o@7vEwn&PXgsy~SFTn3tYfBmt6PP6moI6=dFGF3!x;
z<S5da9L%XKU<i~d0ux3cmg3|#PEAI|$s0I#NN52CH05ux$H%ASC&$O%V$CbfEvPIK
zpWMJD&nP^30g#lPe1J=qQFQV>t|&HRh<-b6O%-z>^A>YXYF?2eNKh9<AS{B|Ho2SI
F4gfSmP}u+g

diff --git a/all_first_pmfs_typeless b/all_first_pmfs_typeless
new file mode 100644
index 0000000000000000000000000000000000000000..336970346beed422f9657b5d9db49519023836d9
GIT binary patch
literal 45403
zcmZTxcU(>X|5tWoQ^`n38L7;ymt>1Lwxnc5M)&Q=URl|j>?pG&Nhv}*N>h7}+a9vF
z-}9a0{yvYs{qcOh@B8(7yk5_zbMCq4p7TD}uKL?$^>pVou2-+#%(<pk21Wy9N|rKZ
z-6d-`tg+v$zhoU-roU#xYH{Ua#ujf@v0}@XtXsQYWE+-<Dzf^2vsBhD{Tr1vm{M<_
ztf8rrjjYiDSz|%gM0bkPlxF*6O{K)8;%570Dl(OM%IIffC~J1m>42=cAXC+yQh%Nj
z=C;XNOq(`s@*n)w!jyWVtfkFC(`HKbl>YwJQ<Ak3Wvy*&>)Y0LYh#O6Z3S67sa1R1
zztxuym383N)dZRPzv?4o8vj$hBdX6+>;8Sh6Wc8`R3rY@o}{fL(@r~}eCQypFP^I}
zdfp?}I$`UxW8u`}v7O<Byj%+VoU&=c*gVQv*Ut0I`8;}lMp>!f*F5?Yb>6zLXFh#?
zH;Qz7=aWoRly$P%FYEkIZ5OV#t02>o)ONE~tz%AGlyz5_>%r4|3Njrjy_c;z)XtPO
z@jd)0ls2VwSh{Lu80qykck13Xio$d<^G32{a{Qi~`tfcK+0<W|v)`tWibKnqb&Dw`
zpD#_Ct{7KF`lH`Z`rfsIY~~!8WPGuL7O#vRoN}XrWV)hE4{Ce=Q`?8P?<>gqN$vIj
zs#W4!svh4Uvi_pXK%v%<ryB{f0aE(FTD1%3IXykS@)>m*y}f-~gFqVI`S<av#ZX^r
z)7L^jnN(WcI@qdp9=YsKGHzlmr@<GSf1NU^h;Hm(IY)b5F=;<IZ8W@l38kK25#FqG
z32E(We#u|Ggk*z6*<h$0@=xth-riV{4U^gr|EE^kbChHwMA=A%+EF~61ledQ-K18n
z%_}!`f4wkLJ8t*Fzg`%vbBb*=W=#l1-|&((%zsI(?(})NFy%Qd9DGtc_}(*W-+HLV
zHjPATx~6@R@`pTf{O$W`-}y3nAp0ZpDyX2xgO<-oZC**TF`~>AYRCRlJC3&>FUTfH
z?I-?^+DW2pvO?_?o^B?{rb_A4YSnh!*lm%+nPgfMmRN2nOQw4ro;9;-m_##A`#27p
zo<Qo))qnS&7Dq3qy=>O#Of>1QvTu?g&!&JeBkwmGQ%rlU7j4Q&te_9KzeI#+RZ*Do
zgJh+lRU|VPWz(T{#y_<JZ$DFz&63*Bmekhn2z*zN%@JjD6>8`4^!b8pfs}4htJbxt
z#)>IdUs7h3Vfo#Zmo)Cxf?XHALnv>Df8H&qTWq;0bc9|a`MGH2j7d(VVGTZbh8%)B
zb+f^8*;_h($zXxkoCKQP@^{Ka_jnqSGW&4K)p(LE6lIoByXc=<E8gB(kS&(l+x(B(
zC8BJpLhUl1zFd&4kkV~y)mpUqb*8*JnKU20y;u4*nM%H9ubLvJ(5LC2o5@R3$zs0m
zsEl#xr29Rf>zl_J<Y=pQd}v_?T?;aQ*lkcIWp+%LC8lN42aD0Jr^aQG%Jr#%JEv!n
zY^5k$1+{j|<~77W-;|jKhJ$6REoEzT|H)s=^{f+Q_L3giU-|zz0y0syULjxP=^F%@
zgOt9pR=(Nb(7c=ux%9%p!sy_@Tx#ljZl=?^JUaB~+BUt20{Rd(BxHq!oJv-lt+Mtn
zA}ZERGyPCP>SLqN-?(2!TC-f%E!3}|xzCESZg^CX>I<*lVZJy>Y!YRgA%Dw1wOe`n
zZGvpO)PBeRsNE^bb}7{E=IM@tY>$+_w^psktaQ6AVWG6{r-Q-cUSZTIaJQ{(%`58M
zu4eF~mx*LGZja~5G3hjKQi4iobS|kZZ~4JUT~6};Yg0d)m5}*`Yu)6H%W0xpmkwWQ
z%E@)$lJt)apmv`q+Yhw|{;751?GFmFLsI+0|Fdfx5oJdeYLD^s<AUsjlzy^S?aV6~
zRi=ANN!M-p;xF4uY171T<DWXEw5G+AZktp~$m3F%PP@~JsCm05vj^MBX>W_}0V54c
zXy7%a`}t$a$aC~s%e3rrdbRA<-QP1SXx4%`@f8c9_LL|)4Yg<fsXfcvpA%%~rS{JM
zqxOO*yQolmiKkx{WLKp0tF>x}Pw+6D5|l|pc25p#<cb~B$9ed=H5ufjJmhNSvNRIS
zT@MGFrqZ@iy^Yj2q)<EA0gVd>QlRe1u2GpO^khxuIOo_@a(ff|<@13wQtJ5AZ|0>m
zl3f#JE>L^@pV}L|y{jN|liJ_>AGNnc*=>c|J3RfaAiF1}yVt7Ky6>i|y}yvY)@#4I
z*N#H!v}u1@YR5wI&u{ScQ`-W1d@OscXP<mh*Ap+^rabDi{n4zlFZq=6=Kh=8o`sYZ
zG|hg)l_K)g>$zYUK0av6NcoNx#U%3(W%r@h^Pk!WyuFtodnmQ{{vWlEM4695?PH$q
zE6Dt$^e44yyALUwroRSTvT2KxCs&pe9diqh=~qrK#R3z1qcXZQ?eO~_UrT6{%lzWZ
z<;C=Luj<5`FG}gN>sH5&teigkDBT+$R6(WhcV3R~R7vx;n|p8Qj*km}Q5FETPyeZX
z#@h!9vgcC!7yqO7r6>zhs14@nA%ZMaN)P+1HZihK&v&I1xxHpXd_gIF8Pci2({ZKr
z#6RovT3Ip0jL`p};v%P3^QUeOzgj>}N5{68Tgs_N=OsfOJ{6JbzAKY9o+zQbe&4<Y
ztC!NRHKDWnbuT4ZxG0N&+DOGA&B#dh%2M`PaY&2edZGo{8%a;hU-|zzq{WJ|IEDOp
zo}M7c5~cLFwenNjtv<9X<UU1c&tEw{%9E@EzCCZ;(w8jf`+Byx98QxbwkT@3E0OFv
z{d%B2ESH#}HFJq8r-fH`<}|KXO+S}U_Z*y1P48m1ZMa@uO}h*&4{BCclPpP;B}0D7
zzklph-abu`rAzHIq<`#xClX4sOi`AlP@B!ua|Btgl%7|s*6(V+pqPj5)bsqP^)WT>
zWPgq9PAfm4jWdiRdUW-n*3p(qYkv9DIqmb;;?sht{6ga{s}@Akj_vav9q%7Ug<C8i
z20p^)&$FiWr>;#Wt^H3v_K8X-S-vPMfZD=;YUR9rksvFU+L!!OTl?Y6&(caoS(!p@
zIZv+;WR+5SRjpd>PS4t!U4KiJQK}gOg5J`X)eft=8YR>H&`<j#uBXw2_)`0)Co*Y^
zVe$M<=W=Lp-lV1{WO;N#-^}!aW&wR#w9%nNTTav<F6y3A5uNn7o!Gd25y`4W**mDM
z`KR_hZ~sA%eU#dN`kx*7vnczbQ2Ujqe-mWirSu=QYJc4cJ)gG#2ZNq(igs8R(A5P!
z2fjT~K(Wh~=jJ@bYv1m!KQ2Pmw-Y~S*%naCE`tVaYAL4*ZJR0b(?!%vN$!8XO9{!t
znpzp(ETIFx9yIWKfOE>9qU;yc{{E--58AVOINuc*KkNM~f+=HIS2nAUNs%=G{u9p{
z@;nr11Qu9hpj4m<5AaiY=i$Le^F}4ni}4SX?u<;LvkjZnbe^6>&tg5RLfa=((>G4H
z7OAI@QHog$JM&b^_%NaMs-Y>QYSqACie)m9`<cAF_#_%WKWN*6>?At1;N-&t`ANi@
zqJzj(fHwOXzCG2G9I@tn3{9$Ffwcfinp*Ne?Y9vOq+N%#!dPUjfeKS?cpeJ01q-Yl
zP%6-#2e>J>oqmBHKe5N&KH>iTb2{}J8{gaEbvmtGSQ+~<G@W!BX|#9zkxuqL4xL8L
z${@WHraxG91{qge@0%Z(N!{LfI9%wMO~=e04&0xYO>cC29=rcBn^*^Q5Sbbfrqn^p
z<~4&-)0q>k7MP61m)Ai|1GAQ_Bk<p`hiT#k7u*Rfu+BVG`>m99;ek5ezgSmHicAZr
zaNLdOp@246VBLXIfgU`-9hVuNJ9}>TDGLAL8?<8LX?kkE*)%Bs3cdPrH+bB3Z%SQw
z%<u865Q;jsvG1OTZ%O^YvE?UK^XTR#wTMG|ig2=J^C6^TIh7j4em{D+f+E_CaK3z}
zf>=*<5Sb1Tj(hQTf4}xJT|R~;J+Q!f10_v;c;N45#0Y;l0I<FoimV?{VM?Fpp+JAI
zzzl#=0Ye_(rtbM0YJ9qtNHdxoPAom2NTk=lcl~CGB%V3S{;W<Qy)iTMy6#D!J#N2L
z&2$pU_odR*)xIgz-|79xN6)ip(WS|js^$f>@XfD%HIwBei@M)ppM#v35ju!$01&1I
z@^+G`L3|8NgTVqD0+cij<$=1Uj4>&)VL*kc;XDroMt}u25-1fI#RJ^bz|?-ZK@D=L
zSMkoMZ{WBE%c7$j=aO!lBN4mW<kGvE1J7)<a>?~ezvM4QxwONoV9AkLc~qJkw&8+P
zKDiCvnR}x(PBUU!j;IbOq;c1#H4<MH5+igF*=Qh4necX!sWE&EO{QRhjRi`Y#_>Si
zH{5tkifjT<VQM1JLxD+PflUTV1*Y%-Hx;85cc?h0iiR#P%I}d?MaBm{j9=ncMT1JZ
zOj*$fU*}xzHDuP6#5PTtKd4>>Sswg7{(yNUJ$zQObW5`;8hI|z@T*f5*=yPL8g;jd
z4AvC;>~XImW`+(Tn+k-fX}q0e%AAj(X*yV7Gk}sNfd}fEnu$q~%>pV+&E|P1Fb6EK
zxj?A^zYJ7s>Y_>flEZsTDM5RY#|4~Z*PpvN$fZpw85~a53<)i!fN`5AZn7?-ia<5C
zC$NyrKUW1@TwF|@KZlHVJXA{U7H%1BZHH5`tf<F}{*==S`^pGk6}+eBql3s60Ab34
zx06gQ<YQ>E1Pg2tP|{?@19eSVV^U;`feKSLJP!qyfCaV`C>2=71AI#zn)Iww*KrkO
zc#8HMA5}r8ZDzfVmzPtARu;oupOsOyM}5=kWu^4&P3o&3?!}aGH9<w^YzeK}ZK+>i
zP)aFYkB_y#R7!1*Oz(E(b}3b9EMCUIY&klJYy}XeY(dNBsp8GQ(q^^Q#+9PDY^{~y
z5VjKYmTVPJad5Q53odyzSYT^-=<mUCEf3Us{{NFSwhogbvj-|X%Xl6NFtEVZ1Em5Y
z4{*;*4C+}n)DEJxt9$ovyC#T~%Qo5iC5O_}oS@@Fqux-7#qkZ(0upIZu---wSvt-#
z-A>Ay=F_J#>**b`<m6}CHgUq{Qi@;s%;$Pc8KrI;(|=ic8L<uMATkFaJa6Rfqz$==
zkD+NZSYTU#lBTUZ@OMLEBpn>LVJNcgK!vFtJP!qSf(5n<C>7Yv1KiY^)P_SXqVwr%
zgGTEw<>2c`WUErcTLt7j?cNHZwVZxlR?hcoT11Z*%FX9h7E$y#y#oPzi|C?zkCu}$
zuh!wjt2|RV<%d4m=@eT?PJNQanh+dE9MM5!dw?*tm$#El?c-x;+7A}k0idMGi3jR#
zse_ml*&(39)M1{70!P3CI|`Hv9OD6Q>Sw@yjrXV1D17kYS8f$)q^Z8EIwl~CQXiTH
zEP9(yy#fY?8ivTpx8Bu_xf@Domh}ST2@gxCS1Ui=UEw7(FR{^->_;WkI%{r^5a$wV
z6O-KcD)wG>934b<0tizlc{|C}DL#g#(_n#}0ZN+A@<3fv=P)U<^FW0uXP$=w7r+9$
z2$Tw3;sI{Tf5fr8W`WUkuG6NkZYpnRORUw?@^NwW^m4<17g2A?$*omP^`SIcB3O1k
z9hyVIUF$n-Jyb}BKQ8FbNhzl0zv8YvxKc)U)|T9BeX5+!FEp(XZj}?ej1D5Z0)(ln
zyq#p~8XrTG3s_*+fs&>hJW$t^D<(ze22_~3$@5U)7Fb}nfl`4xJitxqWaKSNHgTtQ
zj=|Ypf(Jw=4moW8n4)sm_;_1J(*CuoGumBBrqHbskInQ8=$6WFyOB8FXnh@5?v!3d
z>(&|Bj(%88mBR+NnVMKl9usFQ-j!5M>@GTp>>d!N+<806lm{O}(|xeOJb{v?2R!ii
zJyrXb;%6Mp3zH&y2vnHz=6NXa2rMukpj6;74{%fKn`F1Oa7-raF40|{Tu3J8h2JfQ
zd8Lq^uJ)9t)6>a6I#}zFT_&A>#_HEwo=tPw?lM0yCzHfS7b36pN~hgAul@wG6nfk`
zvF*b-$<!n<Z{N33$;5ooL1cbFn0f+QHm?b^;=Iw!%D`|a^T(7W3jivPi%;=_t9=F*
zSRfDmJuW`yfx6r91tvxI5~#2n#Pd)f7%Z?5pj04~2e{poHbv2`MqHpRV+&f(OTR!j
zWz8nmJa?g8JGzvvO!Xiq&&`wL3Vo?iKBK<9eF%-c|7c$4C&_fZ%f6bcFZ1c{-t+qF
z^Giu{{%UWf78Ud-b83rTL-7F@h7KYN2f}UyZzpZTNIr(9S73p?21=Tuc;N3g#0Y=<
zq%j&pk-Y&bOvUg#6o>^2EDk6Yi01)rDt24b*qh68Y5T83^#|_CrDmrOe*bkbkBaPP
zm3O~hKs$b)V%O6OsnrvEy;<RMYNO*i@q@gScE0|e-tb5THBV2r(jQetJ4Sq2-v3?|
zjgGqLpy`cMkpy%QSt1ao-tu;msU$vzrev_dQh<`CR34~nDh-n&O9v`UW$-){$OH>4
z3n&%H<^gVM$UXZn^TJ|it?`%Ecb>*lzQ#zGC;my~*FMH|Sympdr4IjXai@qjM=i5*
z*Ds?j2aU3x=T?$Z`i+BH1FEU5^7hLceXA*bux8Zk*lIF={MuVJ9`CjsbP!oC5T^2Y
zJIPc&A45|CSYU-fNt2uh>b~2GFe$QPpu$uM&qIMyu)xZIQh{<F;HEZgQ+Zw_OQj~W
z2TW<UDwSp~Id&nxe=5CftT(;IkW{*-@9?(fKq{%mo(xp`l}eX~<a9c?KZEKOf6Fx4
zl1;81j*eKpIFGt)9uob<KA$4qG=AcKFrQcjI*6<i2vb$Oon)$-kD=)uSYS0kNz;2C
zsB7v2CPnrUs4(@3=b^x7u)w|mr2=1hfSa<2?Bdz)bqOUKehM!LE}>=vRwPbdR6@O<
zc&>5iQcPvv4o|Vll+*W{v98q{3Tc2rn@jmkis(n59(T5^!_~qMNoOkJOUU(nQMUEi
zQtFj!!Ya*he)SC<MD`sBQ$Ik<YVVmYm}+G>)OZN{i5W}w3#d4P{Kg9|_77NK^>D{i
z603yazx`54y06KUF)6b8K!w`|JP!pLf(6zHC>3bT1Ke%0wYIZr`y?8FVRw{8ha^&s
zy3ksCY!YdHR&8&un@qPe+A0NgPa(4bZ`))|Or^Zg(bM&C&GeX7)2f%gMKpb)>9>x_
z<<w`3$<XV=D#>Nsy88CB@I9dkI*6<(5N=g?J82U(<6~%Q4i=ayP}0<b2WmIre@+)$
zVp3$SfC^Ksc^(S10Sl}xP%6-l2e_#f9VebGd=)`$a;nPC`9)A}(T<zZQzK}C_^#ow
z)sa-=Q+y|{V-yW;662F=7E8CMJTVFPOCwk9>w(wS7tj&qMLxzJW%OaU^P?Ny735Ix
z#^g*24vX#4L1Z0(Fr~)ZNv72K7@9P|0_zBrG->icT~nPfDYDK$g{dw)4+Xk{1*Qd*
z3UuQEZmMxy`|R{aXKD9{)Q^>xXDNKV-O8TX7wF!PYwfnSyh+zOm6z>^@S<_K3u4XI
z1=6C)<;}Zg$I{A8!ya9v3~D#eLp$(mAuVc6FO{o{DcMiQ&t4yAR@&$wvhF~b>cQJd
zrh4)*H0gi^)(a?U(&d4=Z#F$limW$KVX6<$LxH|vf%OAQ1@w7<n-Z+v&x_1Sq3`v_
zAL(70Le`TP+%GtdBgWR%Lz_=Zr!fs`{<sEb(8{xJ4Mp4mJMJ43yh1sPz8TmK8}K%h
zG=FW;Z5x?M)!7&Qn`LIwH{oYxLUtyx{^%ew10YNp@^+FbBR+<v0bqd*1WK9)@jzWu
zgE1+xAwY$xp*#--jKKmM29yd6=K*diY{}@Qv%Tch%kI;6V;?!G1vGx587-%2lLpG(
zT$PiV`R7kBhRUfOog1X@QApP=8C~eQpoDCCv}o6}PdVK+w5n$}wUSad*!Q-~t0W^`
z566WScu$Q$2a$~g!qh0xe@{Y;3>aa`l8pu`4k0FZ!PSld3(S;<{+@)4<$=08<~U4>
zY&=k5cLL8tfr(&&O#(^<Ci4KdJ9fkj>(-_(sa?Co$kXRvlB}ub^XL0R$>I3?S?AWr
z(k$gJiFY!SaMLF_WtEUkBQ@Lkd>v3s7Q3RJx`vjMm*M;HHkMVS?xkn0@ui9u1_yra
zri_oaDd-?FGa&3v<?W<xIE{~?$s8=O=|D--3?BHq4HXyKnSe==%>*h;&Ek0|FdHne
zIY6nvTpr-2%zMTZ+pJHcq^B(+b_b-<+ab@quE%FmO{-pWKb_5^4fnNf-}WmY-B_ot
zt5?d&#a2^IX<7+YZ0*0^bVwQHj#+m+(6EA<Tx>AB>~jSzbPB$q+ZfxBuUiq>d>~9M
z;O!(+7JLj%3&8@j1WK9~@j%@zWrazRSpyZO7V|t5umKBf2~aAqln1z}0cZ7EKdOqS
zuVV_53au09-mwJ1d~XUFZZ<GkxHOOI`A$98ZJ3-q%RSGH99}}#=Q_#G1(lNL3!_d)
zYs$#NFWKjVMg>h?+_!J__6oYv<xiyRJ{%X9p@Ybl17T_fZzq|u<zr}C2^QEYprpx;
z2c#{fI9*(gNs+AqDom~Ac_^?BEHHbZR6xc9+?3U$7LOhJzNNl1+qPA``j(E4d{ej$
z2lM2k9TpA1NkvP|Z|e2O7E!AgPkS~uE~S#Q@|AbAODU?*>Gi@-C6qL`v3I}YC6s3C
zcyRdW64JA|H2K6392Xfnh-^I&rbOOOGPQw^p~(R(u#G@T(<UCMyQMZ`Qe<0z3R7Eo
z9tvy&3v4@3DzJkGYE9{P-JTkLjI>g``Yyh8oZ3bCEsS1xmfD3js_(1oM&=L4I9uC1
zqy?@GoA;dalDggh=@!{Ap2kgN1N#2Vpt;}NPrK&h8)!fC`6UmFsh2|_1ymIi+ldY$
z+XaNF-MpP-%8`$uX%ARndx4UseLPUt)P78g>;O<<%8BQpz(KIU4gsYChk1aT>RK?|
z%yOa|6$h_ObiV3FeYbv7&UL>-KQ$7XKQ{L!!|qnU-?s^%Ik(1oH!cq$?RiD5cZ^M-
z<&Qlw7o@>dYIFYq`2}=Fb$^#6+ah|~Hm*3|rHI%ObP(B5AWR+O?Icsj`52l`fCY9E
zC}}#y1JdsI??}N;V^U;ifC^J*c^(R!0}JdtP%7Ze1KiZdAwgls3?fK(-}H#J%_3;r
zQLCr=ufynRho5Jz+Jw@_`Gs$r`h}2OEk3qcvk=<dE6`0i_?mputS@vPpGb3~eHxE7
z%%-^BeQ#XskWYfi!4}snaSnU|9Yl5!2ve7MJIT~#K8B_%V1Zo)N}8_mz+Y3fZ>id^
z(Jq)2*>#}8)D5180<K_zxdEjDH+g`Y^13#z$>3vh>gd>g_+Dr5(_f$Fev#9+HVy8T
zTNTlFn_uyM%Zq5|mnJ*qgNrCiXO@=MixL{V@$EHhto-=peZ;j9zDH(vn%lC8dj<7b
z6tFSm5e^i$&_QIkfiQK4x06iW<zr~N2Nsw+P}1bV1JagK94PK%Qe>V$g{cQT4+Xry
z0(%IQ3V8DXH`V><E|b>yuD>wA;JC9-CTXjeKAUZqOZ!*$tvG>Y0uOzrCXA97(I~w>
zXYF>DP|@@Gf6A*$X?D^{r|h$3ba?w@&#<a8YI@`I+=4%46w_Telu-{?aUP+A$b5h>
z^_aJlO!@LLH2HxA_5>(t^5=oNrUEc2vZp|Wsb@S71p>hWdk&Neyx;+DN=H8<x!15F
z+HYxAaD6}#x!0_2x;|e{x<xO?mwzjyXF==Mlr$@(%O_?lpRdWM<m`1XvQi4j=Uuv{
z%dJ9^w|W{i&QVU|?yOk*F+om){0+az(r}06B|3;K2nbWbpugW%_yrGxK?7L`W-M7K
zP;sgeh8J9HI9Om2JoNWeBa#Q|x_yO7k-Y{g+(z*{6o>{3><v&V5W@r9?fE{9!kse0
zDWd*6kJ)w+r2aXscv=4#YV)JZr=Lctw64?0Ta|lqsY!oJotW`*QZMS}InS$v7QYI;
zW%;9w())kho)c3+WfOKEXtJb|KA&u$SF{$N%CYDmvN#~z#`AX4CQRUCXi5YN>@85z
zl*9wIrvCG7B^i?<O93iOrSd!!NCOKj9VivZ-~n#xzIN4$RFxoFU{WpHe<+AneR^oK
z;a&taTA0<>0(bhlYgpVKX_QT`3<nJM^eUuJ2ReB)+)zq8o;fzG?p;AHCu1V+*Hlv6
z>t@cU#VY#Tu!mjF4qWfZL<f;&0bweex06ie@G&&yf(4cblr-h@KwVP>m=swdP+>~W
z^H87&EU;prRG@?hxT)g#&0iVb_o2q;HpczQ^r4n5>nZDY@}nKg+(SS622y+LrE_)K
zg;4!nH{6#Feoe!SdN{9*PNf-FN3Dx=$fxaTHy0~!E1_KV$-^2umQmO~zts8exME+5
z4k9Z9!c;kLCz-0?V`!=b3#<w#X{zRdx~AS?Qe-tig{k*E4+TDe1@;jr75Kyh+|;*T
z_Rl@~ds6FCEx)qX4`|qwFJDKiKc*$iE*Ty?8A9XDcdtEvJC;mSx5PL(rBTuEHw`Mr
z<1m?ccJtO3<z!%062I?w6=^4be4N#znv!Sq_UYFjXS$!!L1bTmF!hzUlT3Z%V`%yg
z7T6D<r0FLQ)HU@BlOp>KRG9k1^H88(eO^EbC>2n~u-24ja$m!g=aFPTXLaEQhgY=L
zqg8pImhrT&up)lW?M%9=mNECRZvnM--m~N|&X*<|*}S;&tAdUms5$C$s){1kFV9?T
zQ%zUy9A52krkV`qOdiqC8K&x^16DHtVX7f-Cz)!*$I#RmEU+d(NmElEs5OO)8vKu6
zLsY?}$eIBarke9S6i@|Y0Rm7e(2@uEmI|2Za=X{Y7z%#mymjCH7+Uu0M%wZTu~cQA
zv-j<jINE>jhZvO|Pu7M5lh4daq+MaoUemLZXtcY=Yc=;&T6Di!=h)3mYMZ!yqq9#o
z#m_sFYgdHz1Fg^js~LbW)rPl|Ots}>Xle%*SbLzPsRIwxHKm40k*Na}rZjjS3UmYu
zOcN*-=)?ou)ZM1hxuabQ>89!6%E)Jhbfc#HlT)^w7Ed#6y(FoaZrA9yS%TY}X+`-r
zojaA$x0b!ck=;tk_H1#c|HooF5Ta{Vi3K$$_KzHNy`G#lk6Th;BF7D*&gdYrE<l*-
z%G*h%wD=gBx`73z4U{x>=YhYbYByBvsYVY>imWG4VM>SRp+GONz;uC90X-hzrUciJ
z4|j)`lltYYr~F5i(=u<LUOlIm)1897w~Mi!Ve^@#J-l)6#C^hDlPxKwl%za$d;8`U
zlwLg}v+ev!+TE!6W?$T(?z7dp`p=>&YB+lLw;MJ%`t?Q!k@W$>RA1gsGS!cdp-CSs
zu>L?vlK~IZ-BN~_6qylFVQK)+LxF){feivm1qSl~H?`?-c%hAcIE9Qe`MA|4oKzw;
zCf!Piq&9lZt;5}7>6$~7>!PcPq;+ZawQJW?$z@}Y@t1Tn=xnD4ljkI6(XEt*12;X)
zrN$HWRh*9J(Qn&U7yr2C5gURIA{z>XDPz!ouK*YhWWz9J$%X?JAIl@~f~y?~7T72r
z`unj=JWy-*Kc8TuF)1<=pu+ALo`(XaV5~j>N(IL80JnSL$2oRUH;YVv+{zA^oJI4;
z)F1UUAe-JjSZ=x&&pKGOwOFu08+STiI!@iZyMRnJR;o?*mQ(*Vwpu5wi|Lx@lrA@$
zmr`T>2WDLZO38SBlE$J}xMVsW9k7%E2)h$`J82tE;$vu<3>MfFprpx+2kM%dib;`8
z11e0J^E?!o4#w&Opj1HM0d8vByBxjr2mbVD#<PTG>Hd^x{xr7g;ZtfGdEL{!=}S5{
zCy2!w1=HiG9v(f^!boFCTCnZBBuajMqkiw7xwK0^u(WwWG0h0LWbHq#jDBl0E*rlK
zpI|f50ZSQxFg2UElT6LwV`!QS7T7$Xq-j16{M|AC|J}9#lOnSKDoicpc_?5B#_9v0
zRKSV{xT%lR^zUxW%BCW{IVthE*%Wg(CrPVG4mlleF|lBA4xN0kMJvoThkp9$+}oRx
zLw3&J9LHSDp<SvIRCiCxp;<vScN-;V(+kb_^!<7^>2(ZTzU>N*i`M9Xr3^rrvf=F{
zQ%m?5nwEkEwhSm~TFwJ?O|8JB$ZUZMQ!9BM3akPP%nm3OSj_|6)R8{vjgEATCGUbs
zT|4%OrP*djji2Dkr}g;Ws(-e`)8v8I+d1N)fTFD*W~4PuqV`V*Y+s9;gU3G|S80QH
z-^SAeCiiZfNqq-SP~Y}Bi_{LkjvJr@Q)|#cWNU#iwT`!wOxg1>G|9jM!>z=Bn%46`
zT~i_^@wgLEVakE$p}<D4z%~ITQ=55!n>rcM$Kh;XF}2BOzfV0crYVUwt1tK#lfx?u
z?~~(-Y2Zy|>q)p_RC@jTNV`Hg9eZrp@#i=!^Jrf@Q>AYi9ho|9S)ZtK`muULV&&ut
za$g&l9x@AG!M31-$hHDuY8!7SncB|B(6j?Au$@3j(=HyUYic(pMdk=pnA*ehP+%`u
zVEcelf&DzdO|@jXbMzwfNZZ6l?N34;H7PM#Gmr9VgjzfEOW$#==H1t0O?MU0o_Aw3
zqFNNvy_qxLE%{!A{O7HGrj?Rj+p{6LOUtR>hN4G0>J{|+{yTLq9hf?R4kB{`!qh>~
zf4;2@He!b`W62H!6^F$mc)`UU1q<vL5B)tX9_N9&ZckuRWG8_Nx2Jd>3Y-RG-2qT4
zaFz$S+cTERm-8=XQqvhKp60%pH2B>N@osSzdAKG%9B-UUoqPW1(zZ_?wO8sqc6pP0
zQnmfFwfJcsg$}yf^0Rdwjkta0;N6s5GRuF_VbHN$^0?LJYySgSlza{yu!sQ&x6Zts
zv<WZpF*IES3+xh5(sY>zYEAv;SbqhRBD)Gyn7YREP{0L@bq7GHzzrVYrow)IP&>5f
zH7R>+h;n-Un(P`5jvea}OD&D>1UVHXQDmNOx0|Ok$jIF+;7yx6+OcU2Ju=Cs12w0U
z$CwsS;Nj7#`Idzg;5xE#*W^M<D7o15MJ`T4T+soG7=SQ!led#h-Qr_tx(ycC9iXJ?
zE)Ue*Qui<^GIyZDln2j4f%{;rI{-=r9`FD+)yKGTax;^?q`&!OY|W3o)a=jvm6ych
zbmnDCql?3D((dL#@`btnRQ)l{*!A8UIzKB<&DJoHdN+)+FUPkQAv;H-ARSvu=B+=t
z2i8h#jY^(2zkrw*I$#k45T?9&JIT}|K87YAu)rPzB~89OP<KoDVNzsIfC^LoJP!o|
zzyf;;lnOlK0dDHgD4}P%$wRvF$vp1f=ZCcKXn@^C*;7*4+d^ZIMJ(Mtznqr8&7hD`
zD|FV!%BjWs4z`J=rR1=DTGwNE5Wyn-?yhZJD`_VgM2~+|NuKB97i@fiFUx`GAhPE`
zn0f)Kvxe5dkiEpLB?|&7-f+Qq!3Bqa1s2Lff8TInJW$thI3`6F0aQ4S<asFY3XG)(
zK&e0!4{*mrPmR85p}d!ZB4eBPGT%q~zJ}3O>1Sx^uSzO@<wMDdUW4RoW9X4%JKLdw
z`Q)wldCCT_V*1+lrjk;48Lg`?^n0&SNjC}}UAmrBNryGst9li}aWpz$9Rm=KV|Y7h
zGsf~UG{u1h77vs(CGbF9Q;C=q*;}B(R1(iafn+e29ss2RsXV|<U3la2C#s1H&Hujk
zYA=5m>gk)%GgIpxo%t1;^ilIAdF=Li<DnBr4~iNtn3|AHLz>sqzVITK+ONHMszFX3
z`AqrjI!32}^8G}w3xNgX+Ocx@`6%q0Y3P7;3_zI5;O!(+nS2aQSzv)>10_v4JW%%y
zmy1b}<pC9@@_8N#6o9ex04NoZ^8h#XVVR1Mw6%~7qI&Hcw6Bo-9NrBZonA<1FAeVM
zI!I2IKEm+lCUWxb_oYQzU)(aQIC^+&7diFoY&$)typX)V?TdVfqsX-9BkvukDIk5P
zHNUlrV5$fmu#N!;Qzg8eWU7>pp{WckuyUZJse%XUnySR4$f|$}Q`I~V1>S+N^Z+Oo
zc+Uge6uY*m#rGBYv@*Mk=@)CPt~qaNZQL@Sc2_S~y*49{f`$5--cPYmc-oy4i~8o$
zT>Akd>^9<2;Hu93FNfw*U;if8R`txIwT1J$FEr1im8o6VeG_oj{{bDajsXZ$pLjdT
z)Mq}1rY~TDeFaLIzVSd^Q{OQuvL8T&sh>O#1%81A_8TY__`?I-lv0OzVa^MSNMn@a
zw53anDB^M70k7SP=*XL%#(rT%^ef;%!^0PfNYk{H+l?_rR2(@>dFJn8%AYc8_4<{i
zl-Hr*+^6@-sN3)+6&fn#bnIdL&z7xm)wEs%u0{z6Q_7%a^BUpocG~>e{R{@P`k1g}
z4S<RRNJG5fS{s4!C<G7vJ%BXffx7RtrkKRj3qXb0W;_oCnuD?O04Nn`!2{gv@<US=
z={8NJcLwgU`?{vm)BbCZ4Ao4d$jZFum*!^B&}LJ|Cl1ddA>p=1_lY^w`pe_+3C9Y^
zc(~4}K>TB@ylUv%thki4>~%UcZdp#&zlWQp^n%%z=zwJmK$vaK+eurn4Ie{OTQDAl
z07{zL^FUow9WaTf7k~;=>O2nxG{9JS0F(-7@&Grrb?k)xzjWWy)%?wue=dDXB}2}Z
zd$&!Y0@IRhpW0-RnZvx%QzCO|RQRK7;aH-n-mSW5#@u4+IX3&;unlEoV3u35yQqTX
zraL$J-@!L-+r5J|p5Un12_3ME0SHrFcst2dS3ZU&EwI440VPe^JW$tEcT9?`2T)<E
zC(lCx9WYiN0Hp%DJitx8-uCiMae^DEeH=1oRMT7Z-eo}bY9k*ywEehu_rN&1IJ1+w
z_$ZGCeF*Q@*Qbo;>`ANGisvJ*7AyDm8CyxCa(n!8D5@k6zr;xf9jnM_jn}U=U2xez
z4;`?K0SHrlcst2dUp|JWeqe#=10_xUd7!Q-15Ap{5U4O^#Pd*K02nI|fKq`$Jitvw
zol6*`v*jw?5MuUB^Swq(nuIuZ9CU}q6ju)Yo$!ndm7G3pKN?M^T#g<i=M3_$3~v1P
zSP2<E_LMip<?M^g7k^wdwVFN#HEGbvrJ7R9K9{Y(iA!;V(E-aCfG{<bx06g6^D#6H
z0}E_8P|`Gl2kM#{iAj--0xC=q&qIOHV1bzcr2=DkfSYRKy=vF7l1K_^-s*Dyp0DU*
z^Cuma7DrKg_qRu~x5ZQ1tn_Q=jFKq$l~HxIeL7VJ_VT)jqi5SAQT_Tk710wDI$AKZ
zjEci1DdpWPr|ss67tRIYJ!OgxA{z^Ysd2oWWNJJgL(>GXz$OAEO_O+_uBpkG6xkG@
z!ju`$LxHJaflUKS1<ZMXoAPcH(kpEHTQUuP-`LvsEosGQ`#c?(LgiogH3>eJN!#9z
zX|q8;m*fTG3Vo~cXrfb(jXm)=Z&cWr$%AdmsFm$pA-HuVogC#K**v9^hNkv!uZG=`
zO-Bcj%>crbz}rctX7VvK%>oN-Hc-+uhX-m+{kxHgi>7liDYAJ$g{k>G4+R#01!e)1
z3M}LSZfd#dn)O?D-N9MyX1#Vncj&qJ`%#Za?i8Kgs=ZFqLwc9=;9BK=Uuv(N+@$T4
z=cL(PQ`sZ#6<vBF4^DSVq+o-D4%Mr(D9``!=08{?Hv0L{@s8W@`DKX?B3lH6DJ$Mi
zGG)!j(6ksVFdLwxX$cSfHTC~jJeFcoWXpgGQ_FcC3akJN%oZpWSjhw2ltICES))c}
zRC7Z<i~T93eYi!@>rN?czH+Lm=9Ln<s6R>h)E+!05w5wTp;i%<_}{x)RZ>jTUf8Ku
z&n~5m4Bf$pu@+3vdixr$PUSQ|W{}B*p4k0Xp@YcmfH1Y1x06h*;bUl83l`Wqprpy3
z2kM%VVNzras4%sj=b?ZI7T5-$RKS4;xT)yN8?xW{U!f1(`u$K@a+MmJUNJMi<VyVp
ztneOk+neSa_wBiBLJ-b0#_ehx8b>z2gY(X|P9dF=UVYnQ*|YC(v$i#p3Teij<PAaR
zi|CQtv;|SGaVLEvI*4o&5T-VR6p!N?8nZ2UZOOI*6`#r5@PbR-4i?xB9{T&4ypsp&
zzSVYNQe?Y<3a^el4+ZvsvFrdS71+lE+-r!V*72WJJ~V4+HyfSyzSL#bx%s*MLa1AI
z=XPDNxb&B+-x9Ndg%l!>+CQv8IXyo#q1_?eF<RsGu=BjaQnKizd-zd82|Zr3&ihAC
zJVC$T!*xYFtnb*54p_qggjXltPTGJ6`52lGfdzILC}}#v1JYZqcDLmBCyruLWXFIC
zQ^$E83Y-9A*#S^0aEb@GsZ{mPzT<16$ynyCR6Qq#md;(#*?UkbB`j8Xd^Wg<I(W35
z{jp62ZN0TTv1VHp-5R+&^g_ETO5b16J<P6>%r~yo%==zWy#uz)iYvy|@i-H!=C5&a
z@iaPM4FeFS&hmDWsdIb`P3OS^a|TM9F7QCzEp-u-BD(}sn7YjKP~Zv}%MO52fonX#
zO&xCa{#NX3Kg!x!Hlm5k6RL^RqLQ}3G-S@o77ZpOlG#|j?n523sKjyq5kEX{FW-A`
z&*qN#l+exU#WGwb9{t!hD<mU}M(m`qffX6#oa`GWdznGZ1s$-40SHq!cst3ID<4CX
z8(3gBfs&?MJWzK_-NvNI?f@00?(#eoxCh3v1E5sEg9o^&H~lv@iC1|=<`ySE)o=fr
z=GquJ%Ci#b#R<J`g?sbqp6ECrXhAV`Xk4MG@eB7SVwc(mdzVr~P0V-K2Bnlf?#E0=
zy<%FeKfo;(ivrJY*?VG{6K=5IM+dB70K(J*-cB;*#mCU}5G*inprq*$57af~gGrG+
z1}aSX@;nsq17q0%P%7Zh1KiY)-1IpQW?(h%>(<@;2jTp7d4u4mSMo?>bJ&JidS&!P
z)1lMV3zhWa?c>iox>nN!S=t`I>?+#a&-rtpLlt#ce5={)6P2W_vflDYP6fRh_vC%o
zXk3R2KnJX00K(KW-cB+V$j8w194xRGKuOa}9;j<72$LcU1}aR2@H`X<1q&<;C>03j
z0dDHl#p(om!vvaqu7!Niw*(rq^VIEkD>A9@bPKncvBk7}+oYScpaS>TFP(P7UD%I{
zmYzwOS4D}N&YGxxtE6+5ntq4z!wD05od0R<SV0zN_qKX%4^t87AhJjxOugdmBvY^X
z7@DHM0*eMpn%?lhUsM18@fCwfk;MWPrs8-W3dDm2mH?CrB=P_^)kP`v_A9J94_T8>
z-MSXjq|b%w>suC+Z6}ju?gxwL*4fLYom1sx^)9=@5vzO6-iOZ9&@ZMDd)AMxkKJ$m
zlHc27@$uC*q}~TF!*bGn&~|7nZnv?w=peErAWS9mc9N+SK8B`Lu)xxQlBRSXsB0<%
zlOoFmDokbZJQT<V3oHjH70BfQZfdO7d>6YZDP*^7!@Q2$QmAxd;iCJlY1HU$bY8Eq
z*))FF5{uR`xpcd~vPG*&`D7F@-k`&zT-s}Ot@5c_Cf(PM@K~#pOi3yFi=P}$q?exQ
zN0$nT#PZNVWcfgtD&XxTQ-ypCO>(fnihz=)Vjifw`;}l)WTilbsWP620_9+VRRE;|
zl{~;r*=VlluX{F>{8|rhaO42)2Wju<IyNPacHHuwxi}$*dhPS9`i;kZ`;u&le0&L=
z?7g{R2spaXt?}F*C8Vai{YdKNVsgJXbLzUz#q>}m>rP?|EZM6<2a#0+Vd@=kCz-0@
zV`zF07T5=%r0F9M)ZJ2_Fe$RnK!vF<JP!rFf(7;sC>8k51KiYzyNB4vTPf6{FuVLs
zlT^|x4BzG~W>OnAZBmVUKE3qKaB9{AKLNDu=g}qsML4=OEK9+&iPLsB4GqV2n$Y8K
zMm@QePi7z6SbR0kqYhTh>=$*wDd7)v5ZO;4O#K4==Yst}WA+;}mh2BuajH=dKXt^#
zDuD&243xyy=Yd+c|GBo^0FxqX2voRj#Pd*~F&N7ZfKq{`Jiy&<Us64D=HwvS(Yg2W
zSM`Fa{=V3DCdN@D|KmF*;%f%kjhnkur@fqNOe(c|;9m95EHmYCc;cv}WAf=rJRIHU
zuG`39jVjvx<VkWA?nQT*)8U!_QhcdYK?kg10K#o^-cH(ts(cJhEx-b636wOo;(@xR
zT4PdVZGZ|>ZFwFFv;$+=0Z=N?fd{y$Vezp`RKDM#CfaFJ)1Kd>SFbbw#7}=pAJrFX
ze;J!fD=Lnrf5bUpbob-s1IJX7d|O^=)Ts&z*Btt1b7VQWztic|<4_qbvgkUh`VF2N
z&M{ot_$hAAs-XkcFaTjngSV4Rb>w4c(gX{v6HwCBnFs2c>Vip;bp<L+Y4JQ1=my5J
z1E5r(I}dPEUtYMhjNBMcJ8TN(S7V9js*^kBUHz3o1GOJV447I%Z(G=2Qyo%C`x=aJ
zR2^4MxjQ{fZA_|Z-~9eRS`VqFE@t-%dTLkGp!LgZ#CBliUk6p%;+t|0bif)0AWZ4-
zc9N-Hd<;#xV1elYB~86~psuMtm=sxGpu$u?o`(YZU@SWTN(BsffScN<G&-jK$6Ivt
z#<kNe=G>-OH{JKi6Yi7dk&df=S_D#w`@^gKonMjMQ~&Jfxo@ey*Ds@M?~|!*@oVeo
zk7=ZlQle(@G>hsD=xyC&VGe1VT9;hjm4hF#K?kg10K(J&-cB+#kdL8h5LjS?fs&>n
zJWy-u??Rl~^&LYoDKcZA!qhOHhXTXF0viF83XJ3dZpt}leRiY47fJQesjTpJm&jAp
zFxz_Dl_s8gr#9Wu57#(uZ}Jc$$nN;8Z{B5Z>8}5mPr52e6x4A1)U!L1sGiorKCdE@
zX|S$a&g#$<EQJ`EtBEH@*eG-m83AExG;b%FGT~!r8Uw~74?s!NSRSZrY8)o<$OBMe
zY68zgfr(%|@&J?yOy&V@s*$=`<%p#o^u+E-kW+#O?W^D4;f?u2I@mX*f#-r}6eEtx
z8D$ktBOBP-=dX*W{)2L5VS8if#L$bndo-h{&xszHgZD<#`!^5O+Nwp8f54!bbJ|7X
zkq30ZBM(5Bn#$Wrrl#>RG?|0($OBN)G=m4~ni4RHM;?F*Q?qy;3d{!Mkq4kuU@i}E
zQx}%@5ifVlrH1y#PgFbSQkQk#HOz<PQuxl|JrCFB(yyuW@48&br3Zh4{NzEols9hn
z;cx5nXz?B8L+br-!KU^0vQ=sY6uNS|o#D0u(h0L49=;dX<oS=(;gJU*Of3NY?+*(M
zWEPmWWD9|cBfceGaLJ3n0<+?wzh5q`dEl>S#WfGM7?UEi0V+H%;dv;q6pXb8K&ilT
z9^jtM=Wf{edV4&*?`o%h^H)5T>{{7vadrw_AMUf!{!KpBJ3nCgdi?I?&1+{mUc85I
z|J$353Mr_h_Jv`cW(QW$mdVpX`i-ilu%9=_IxWHDnHznT9GBxbvH~5jkO2tKD|tI<
zL$2auXtDzfY&B5Qw1x-jnp%rVk*xzNOxg216p(?j_5dgqSkD98)NAwKFI@DV;0GoX
z_Un)Kr=QadTF>YgN~g6C-1zOAN{`~&EgN~ekY4MmzYCvNPHLT+u4+7{f>x|F7_s7Q
z1vO~%+R9;CCB3<_`0;|xm6TSnGq_?KE>4N)fQ1Y|m~!ClBvTvt7@9VL1-2O|Y1+aA
zb@$G#m=xJIpu*I4o`(WEz*u_#lnU(P0dC5qAVcZYhyv2s+HQ>g_yW>2@brpwDWC$?
zBig4r6q3!--&THvCmsfMoa%y08?DY-&ry3*NNLXdD#`}S$?e_3ys22|yt^=a@cG?x
z+7fa1@}Hgf9=01Du#f==Q+s$j$<$sxhNgXBf$axMnhx-Qw51f60h}-?vV%Z{sY5&u
z1rCD+b_6IDILZUul+)gI2P>Q0qKWtN9zJ+;i=1hllKvwvIyFmeOUZx`8ZtR|TE>xh
zTAS&vG{iKM-tBN+d+u@;UEcD0d1>z)I=X59^a%mElw%%ndoJ!6HHwH>-PRJPSI5vn
zWXFLpbprJN*Kr22lXz{(P5~8fwbOXPrJex`>?{xceXE`0fq#1?%<}geJC8|`IRh15
zFYr7RxCq9|1E5skG7oUClaBtLZZyY-x)@)f23LHr%sx3pUF``4nlH;6c_N6+RCg?1
zv>!iS7M3?qWp+Fbb@6ZW=~)`-$%l=oewatMG`k%w-cU?&_xsNt=3Gi@4i1?Xp*Rn@
zf(}^50EE|Tyq&ZGUHBN9u7d@311M>7<$+pL|2boH!=%V=0u`oi@jMi`4aUj?pj6;4
z4{%czZZB2yzeLj9s@onrOkPn(q1EnJccZDEedx_o1qr0{v&q$lH&W=5mTLR@=d-9_
z^yT8kHu-ewo&49dJ#s86i^$K#uZA7VPd$7;y_9;_{5Ib90jBPu1C}uWVakKIlT6*`
zV`%aO3+w?<(&WVhbxl3Qq{zI13R91G9t!w?vGM>Y74YQ&Zfbwp9b1og*Xf!jOK{wO
zgF0{db8>Ijd${g<A>ojB0Qoh!dg4=21QoJpj$ewBXyJLQil5Ep)X3FJxi6KIlkdUd
z{3^kM!=1+W>0C|64021i>f)g2hYncA0E8)j-cB+Vz{k+^6fCf3KuJ>|57gcLo?}vE
zFMtYDFL@pc1c3z>43r9l@BlY;;;Zd~d&k{rziXS&qmk}(XG5~`?{H6Yv-gPeYWbLs
zoGDwh`DY;c`Ww7I&^dxWyzg0=P!mG8d^OiZ9SNq6_cw1%IUh_~_50Yn>V}ZU(ea%o
z4+|j{iVh+R1Hx1|Zzq|G;A3ct1PkmHP}2082P9L9<3|)GMHUTIn0mwWP#^{@uvnl}
zAdUyPspDf#YRS{{XztXZ#l2$jW6}r56g`}lN2=ZAzNI5`>B@}byUk*<DaY`g?g6bV
zy6LRh%xzjBU5Tl<+pDyM%+g+NTzsg4hL5?qZuhB5GCpy9sq-xy>EqEsWC=i+O62V%
zQ*Zegnv%c*O9o1sQh4C6saj99XN;+s6j>TjVJe;Hp+E*$V3|OvKo$>hQxP%?^<T#;
z=+UcJS+U?duTE<3+)+WUXM=mYDOFIhb^cw`N#zvXRlm{8k!7T0z4!fe{5~?R-RLoD
zQx)w=cgt;LT1_wKS}D1#ucl_jgDx6u#O{}k4kF6|!c;DACz;CQV`$0;3#<SrX)5G_
zx~Alq6j>2aVXBztp+E^(V5LB*Kp78kQ}xbQd_KP<odRN%1B_Ou)Ai^5G=8L|)2HVV
zzrHKul7?CT1|u3}QJU2V=P^xj3FPLe_9yygkjsGGZ-!4wr2$uWuI-rjmb&anZFzQY
zB8|NoGyn7oT=FPK2a#0(VXBh1lT20dF*H?!1@;aoX{zCYx~AS^Qe+>13R53>9twN{
z3+yvcD)5B|xG8hFl2>E35?bAI>9IWg7Gl(gy{DezLFR@vv(IZ=l+ef-!;|+lE}_Ub
zA3rJW#XXj5t_C~t=y}uRJ|8;k;>Uk(20LrkFCiVvq*+TJ;~MEl_joyevWk602a$aP
z!qj)rf36M;GG;$8W66F36(3-~@Pdo|4HnoR9{T$MRu78`{<{6ov0e$2B2xw`+}7uL
zD9`{bu!caXKqDUDZWsEuTNE=tf?O^n-Whv2f_BHu8dz}qHEp|HGyH@}JYDdrbg<P(
zq6Iy=^!V^JjnWM}nx5U0PG1w#vzL#`pvU=bV(Smdq%ns+zTW*Tlk(I`Ob)!lv+9k}
zL1ayUaNCr(lQy9WA45|!u)vxFB~7Y4P`e5LW2yxvMb;9iFx86Zp+IY}z}f(%0&RJK
zn;KfWEyQJV0<FB*DcjpMfwV39XvCdNCd)%jO%5;1p+k-n5-hQz?W?0@Tkmm2bk=m1
z#e2gd8nR<zSg;0`pSFAc(5gmGBTi-<y@u<MWo;&n=#-4ps&?ogvi3lj>cHDcrquWt
zn$*Dp(*R1EI`Tl>{ZbQ?BI^WHnCi^)P@oG~U|oSy0WBWjrj|~XMUD5wn%PT7e&LQk
z?Q3-+aa%?asi;=$9$rvLn-4jRThb3}5kr(tTyn^xhRqL^_QP_$+wKc~Gd!4|6moLF
z9ISxwUi0T@oK6YVoYeV|!XF9lh7Kar2EtT#-cB;rgO8!9Cs<%QKuJ?C9;myebTKJ1
zJ)pu=Z=Qz&eZT_i3zQ1<;{k4J@%T=Q7OqJqx%>If(eh-9xzhd03q7nZJvH&$0sI<w
zi~VCJuKrd`@iQkU7Hlo26FM$W&o(Ki9M9DmyUvx-Ano<1`~NJZS1r3Y9&oRe4!R9&
zGx7#bHT2OzWc`6KWdQo`GX(}im?5StnGsO&ZX19XT<t(Go`B$?zwfreJn(lL;zaT9
z&dG*g5)Ul^6?Tnz9tsQtW8DEzDlmcvxZSg}UBn-Ohshzpq^RY)!(=^Hw}ZyPW8{9t
z*{F2c8H#JwN4<LPRSJKd<-b$7PlJE`4p!IqqZUV0RQ<|=DC$zbIW_%a@FUNfo;TVh
zQs9l_rOUVBfH@K!u!sQ&yTsc`+i)}=Lz4+uU}J!iCQ}}$yA8);Qe@+R3RB~G9tunV
zW8DEzDlmx$xT)Wo_qv9(KT1)j)pPeBIZ9dMADyl+Jw?NnTpooyyiP+m2Q3>?{*aES
zUAj5rS18?`;BwmAEP+m}GPm=4kxrKFc6u56=hH<N>TH0^uKpRTIxoa8u(QeNfJF>I
zm@?z-BvVuQ7@DSm1!fMEG)?D$x~67eQe*;9VQMDNLxEXftUCZo1?KPoH+5s;_)a52
z57Xf)t%S1DBc$Cu>uOK?a}?q~$nV)t4;uY-gYo{cFDW_TSgTIB4Y;9`%IW4QIn?ab
zs_v_%%IUzC%LnFbmC^l-cN5wMlv7wr(7K=LxJ5A+9k7T22vhTUJIT}nK87X>u)r1q
zB~6w*P<QuRgh`QE0TrgKc^(QZ24md;P%5y52e_%Md3hJTOycOkpeY&kN5@g)XpI@o
z2E>v5gpWp#!FjJ5<p08di|~4^<-yQ{adhIU@*W5L$nn?LJ{Jw<=h2)bwfGl!Mt$v^
zl+$^u@$2%57BdR*dv<IoI$#k45T=&%c9N+Td<;#tV1ca!N}5*jz+Y4Ue`aNeNs+Av
zDom~6c_^?JjCBV<senBXa8m_7`Thsm<Nv{MJ~qs}1wMDX%y_JF7l(}SkKN9%Dx&&6
zK1scB{F}03mzT=>0@Avf6uD_vA+>oo;aMOa7+co4^;Q*IIb9!>aclA_Ikj$8QoR?i
znG7AUhye&w>v=oLl*q@>v;izI2cV>BBM;Qw{Wf7zWSfBsQ(Jf*3Ty=nY#UH2u$>3E
zsoP5LUv{aNNk+zpgoL~dde-iA+7177YT%gc7ui3Jdd_)%shM9Y%`x8L+x>nj&0T!{
zM#%Rp+L?BL_?HR!^m=;Ja=*MHx-+4}h1CHiG%>xe?d4o-sU7GbvYkMf+6Ah!%-PV8
z?Z&Jna|9|5jeGEd3*HOHLl8Xl_t3bX2kJULfJr>J08}_W$n#L(5E!cufKq`YJir}0
zT+ht88Rbt0Pk73$ngx(v3w2tZ`;3A{`yM*-Gl+tecl2N{!l=XdL*X@jqey>Cw>@e<
z<n(CGUT^Cgl{B}bgJaLz)zsA^qxy1YHD#PkYP~N9C&)+90ZSQxaD1G%lQ!cCK8B`~
zV1b<iN}5jdKwVR3Fe$RLK!vGuJP!rVgR%MmC>6NC1KiY>FvEe9uE&zE+lJ<{W3jYg
zf6cs^uVd(zNAv2i8E?qMb-2=s;n6q|*`V?rrmjxtuQpGwfDR5_YxKmXoWzEk8lH)-
zBA5KeO;Z+E)4CZ?G+Wxi)J1f_QU)MQUFPj1Q&;#Hny!Kcb`2<La^ZoxZ@BB26xj`+
z!jvn|LjgB1Rv!SR0=Ia8o4RT{U8(u?66*8*dHwQhB{YKesOcQQZ;ae|eeexVjMwe@
z_T|u`Vw!mNoALK8MHFZ5qn=e(N-yU1GwBvxP95z$W9F(?Qt7kqMng7M(xR?DIW{}+
z{p>b6U?~F-rtb1~lBs)q3{CD}fq4KWP4{`A?w0bzq{tos6{fs+9tu1JWAy=0D)5L0
zxT!9CI%*%CQ%cQsHM+l?TS}kazOQyzR7#7L?@XRyQc7aTBQpG4LdyN}A7t8<(EpWn
z<#9Ew-@7Xc;ff5Il9@~;GPVqvZpLGbGK7#zY{#XNsWQuyA|z8OB&8HdN<_)fJkRrF
zC^P+@b-mwnefs_LtoN+-e%4ygXP>jr+54=$R#eW0GM^WC5=qHm<#@|V(r7WSw$!4E
z8X7cJ-*uyk%njXEq~FHBWOon&r3^rj3gU2r)LmX9r(j5#6ao}E-DAVABlX{n;QLU@
zrBI-P)C1lJ2Vsz?J^%^_;cQ@1IwhwLDi27crJB*<57ZK=m5uLRmyme0*xl61$^12?
zuhYN!()lG#F{~&*d?1!$h6HLwuE?T2ElyNz$S<M`!=GhOx?MqiyG`*_ZdpZ6rT#%C
zeL(6FBA}E32vU(8PLO)cYvdFKDU+T6MNUuIP@hyZlyd1AP(kWB?}LLFNK_vHg@YGt
zU{b?kG#eUbq|>ss25<E4q*JO{)g`UgY2@wwWx<<|Nwl<4%oL3)2^3YMb}@hA8&co@
zer1z_9Ln06WnOXyeM?PWH$LTwUgBRjy}e#uNePc)jlMSksW?PHDFYCsUU4`<>NT&C
z(;G;c6b}?RC9pw|Qrz=Mgi<ag0TrZ@c^@34K%)8pC>*4*fk~Yoa-zpTc_cNqc{RVH
zGJ-;0hVBTO9ZHi`G@6fY<VQZUr>ei~dzE%+o_OFAeU9#*f4HSiexEk1KHSo=Jf7yN
z8_Cy}6j9^Ddhs7+sEFu1EZ~F(er}{A0!kTxAeG7C1gR`uBd2UgnUn(*Ipwn97pecg
zne(8OOZh+rsRG^y2ZfNRJ^%^_#cW_wjVetJ8>A=F^rI%*D$smCQ`I=SbWRefzx(I!
zF44&}KQduS*>H52JwNw(+vqg<xou_PvuCN4+v#h;!m<=fnCe{js8<RVMvoc$Ryl>f
z4T@9m0x6Xs0!kTxAXUcU1gUafBc}>TnN$fBIaRTtKB;Ob<<eWAg48?S2M6yVQGEau
z4r<uIr1pKz*zxC&LTan|>Sp)iLR$HCMUzoQ1+?PSQk`CV@w}VLDYq4FIV69eH6VXS
zCSADelGl6z8h!S<W7YgyKFzM%s&-thm`1vuYGX<zwEp8P)9QuzX0Am9lrjK8>Jx<G
z?zWz;^ck<sq%T0lck)-f;8eds%B1gX`t_asgAM<3{p~vXCzNuj4yfR&gyDAV;15XD
z9sq>{Wi~L^+xs_LA8L|BTK-3GS}seX=}wmWN1aWhGg^ioZB_DVd)t!S`?v+7cQE<%
z?F(hJKW)>uQR~Yov}R<`KbSV;J?mk^MkW;$?zA(v<JJmVIk0coTx%2`HADmyG62C<
zg~N#hsLE^P)EH7GH35p8nzErjDK#kNQZt}}RCC@32Q45`djJ#;TC#yj9Xi&vtNcbF
zHP-&?z=ymbTH5E+PWxt$X-bE<rqZc2YFyr;Q{i8D_G70sa@bFlf#t<siPx{BtvW&Z
z)(a}>p3c{muQub^tDL?5&Fm^k-YRWnG^C`C2q<I#f>dh`CrD}V8acIrlu2!YBByq2
zs86arlyXTEs36sW_rXC&NYowxg@ev)U{X=d)87<DU!=f2vz@GeUZREirjGR4pXMfy
z{v2rgnk+}Vsk)Tq(Hv*D@Eo*PvHdYWpy!eba$Yg~fD6jPPAgw)L06!kcB<8NzY3bO
zKXTPwocX0Lh=4){AV_uPaDtRJuaQ$XNDO}fiky0|p+2deP-0#KP(i9U?}LLrkf=QX
z3J3kzz@*B9C#zV!zCf1IM-IMU<UwnKoo>&Z;YC5m)?InjJ)CBirVQIKEdhV!-9BAs
zUPxnF^&gz(R!9qL<Abx%L-&t2!_=PE7E=E%<Ms5L7Ew^V&b4OEP!Ozx2q<I#g46&G
zCrAzCHF6pRiQx}Gk&`YP>XXug67w2>3R3#K4-SSxqV@nN92l^HNlhBqF=<wyEA{NZ
z@Y$z9XK2H=&~-B^yy#@%p@tKdMNuoIx82Vdq>{@&W$WH&6q8GuX@iJYxZRWQ>lGVb
zOkrQ@KCbn`<;}*-g^{N*yl~Hs^ck+Wy*L~ZP{;rTsgWE`kQ&8n<TM%*!ykYmCqp*W
zCq+<VUIS1;Y8>x_gYl54Jpc*^Mr>eG7WejA1;%I4`ot-VLLpPSHj2HKmO<n0&u#WW
zmP!66oUW<rX47Kb>H6h9xwP+u^4^Ip^GW$eJ6-Le1+;j=a;5Cmg|yy0rg11vsyl4-
zJ2d-%Yp{ukfI<c!NKNK&g47gVBd4j582$hhIT^E|KB;L?VqODKL23r?gM*ons67A*
z2eaA0q?8Zp&9h4Nq-j1!g9|WI>waIqLiH*?I+rAid3ZX4g4^y|KhW|u1<wsE^H$C#
z1HC@VJ^N%(X5jaWL9bKEJ$uv3H*=F|sBKbW(WE4D%M4Ra$6!ZEh6pHR0D{z94kt*>
z<27=c4~gLqK#>#o)cr;3x5^b0C^4@As35hN_rZZFB-$+kg@Yw*U{VuCjP8<QmqUlz
zAK5z`eW<Iwf=AoC=1@qt36WPra%gC6>{8#SIplNYfx8_B%4uW`w^6^CLva@ZGQBNx
zXzHfx8MicYXwLVJ>9=CD$=YJ8&bvsI<19r4^lk)#)N&3dNSX5*Ijw-i@CTsCX%!p(
zJyO4svVaov8h{E?Yj__Vtc66oMWAr7jtxx8F*D4!*Nk}5i`GALY(+fTxlgc(+K@o&
z=f2yZW0Xv(%ciHdOG~ATCSia2A4Nk|OOHW5L(*wiub>R~52@rhblLiDxQgC)$*scx
z^aHMH^Xgo;armL}HzJ^SBM_uEK>U6eQFo|h1(lh!5vcfD-h>yNtptg=4{ZANwJc{t
z{Ws<oC^4!5sGz%*_rbw7NR%A_g##NlFx_bL1>c+_(6i0$^q||ZWt24Uk!#JMavHID
z^}H>Z>y?_O_vF@q3bL8*Fm9JrL7$UthNR(IWI#>a=fPA>2gep_Y_qGT(%RH6s)wt|
zZ}s%B)<<wvX^RM`VE}@z9fuRga0joE(@sdteE^D_cC(>AsXb6)R0B{!%7OR6fg>c!
z4uHbJJ~l9^ao0=kpT8SO52M{n43&cDN_c~JKc|P#<yD$X*S35_jSI%59dGjlx4q3a
zSU-9}?#BuO4!_38PS?|?l3j|(JN)ei{myt6?pUowmy?zB@=Ld`_P+R9-j4{VVE}^E
zK@KNK9pW`|It+=q4?vNV6C3K2ItnF5H2@W)j`Kb^H~}e>P6CAkXErdYCf$;JJ9;P5
zg11e_d-!9FRjoroWe3bIvp3CJ(g+o(<4X?pM`h0JRTFK7SEtbE;QJQ?`W8^|wX+Wn
z+LqFo(cOGo^}sW58wUTetaBAEq2HbRV=%5&P9cI^ash&rD~A)LPV*W$oq@#M2cXF5
zEF0>7+Rj0VQ4K%^DR<rn2Nxh^k_S*YxX1=3^)0Kv*7ld#6z91(Kk!L5*`&47TDLTt
z23~T%Ib>}XwQSNv)%s#4wdmM+^6VUpA((vensIgxrE9Ip?bRfocJ8{R)kPnV=6jiL
z*RLxi(}6R#fA|w$%a;&AE?owK)D;dVNL}SMa`J@4+y|h@=^7jAlk$cVqZ)t;QrCGO
z9Nd7ENjHJQfiD}FRJ+3_$;*4^)0B@_W_WhaC;7$&I*o7Sk!rzeXBE#}S~IB5ZFD-?
zROCO{b_a8ibk?_R{{+vMj&`gv>S&is+6jxbyQ2mov><3(NhHeJmmQl*n1dqSLIk<w
z2L!3x98Qq(=QVN)fW+JfpvWnZ4fRO{L5WcfKn1B_-UkOEkTU5WP&l~H1}4?%<?Nc0
zDM|EnYCD4t3zFz<m2abFO_RxcQ=9!Kzot-@Nkh$dxbPYlzU0TE+ZhxP+5E$&fOKj<
zD$nIiY#NQ|H&=Z$=A9_T-r9WHH<exwZGG*TPbx{Fh#;3906{8@!wFIkd5xUHA!X7d
zpvWnL4fRPyLMfLX0~MsAcpn@*fkeMWpl}e)1}1geHp|YoeFDwC*=|{6+XU)2x%ZXA
zUI`SW|9bnQngmj7lX}dkc@kCXWgls5olNcpLvAb8#*^*ieTDZkUr|7Ef7RguG1R?f
zkGx$`&*;Qv>z`pB&q#WP2x#621gRJfCrHKe8acgylu2<wk<&{y)F<@{O1bnJs37%*
z_rXCtB>F7^g@Z&kFsYG;#|@jayPQ5hjLLX&pqzGJ3GbzpP)^;oC*(M`uAoNp^Jil|
zm(!Ux5j(cHV~%T+#81oiR?@&*p{0vfSCPx1+XZe#RTMm@g;B3Y)if|Are>BZPU%UA
zfaZ-rkV@fjf>bK6ky9F^OiBleoHE!@pHwE4aw!X_AeGJg;2;MQ{T6}3K^_~Jl#yeT
zhugasP|f|x*48=&WRzTbp`ip>(D2iZfC6g&G`lqUb^$H#>(OHFDZE}%xpbmcA<2F0
zx6SWYMBzi<U0MCSh{jH-zJCTcH!9{FwRA4Ui6S2n(7X`{QiTxzzxHyNR0O@5R18#n
zH<#cAr(6n&;SX&3_1#>~hJQK#cF9};CFV5%6`ZSh9~@LeqWl0T9K2%#bG{KjW2C=w
zI&C+f*)Y60rb67d92~kWof@}a)zQ~5gVd}RPyP^=L3eyl+fN*vNm?aW=E=}CEY#)Y
z)b<u-6na{>Fm+=kz1CKVSUR+d=A~I~h?#&h{(D3~EdvmoYdD-ZkhQ!<P9GsL`~fI(
z`pky<q`p9jc@01Xsc*aw4!%Rmq#r=x;3pfHRB_jyA<E0oQOg<fko-gE&@4qMbNvB#
zGApS&wzvHy3OTXY(%|@2S~tl$V|@SXq<7WkQAMZcbo|<aS$6nwIDT_P%&5a?wSBUC
z`-DN2WN@hJ@5iS2ZmvTFxuk@lGa&T`#Q#V&;5BkmhLlMSfg-0yY^eY4RDn`1sR9+G
z8uLClXaXsdngWFbH8wCQP3=z=FIH#KD&MbuEf#0d{_kf`c&BEPVTZY4ZT`rlP9sMf
zIXh&~(vi#e4Z~m3C%em^Z^jKBr<#ow8V$3lLAG_X;Ilboqq=wfkiNNeL}$0bpkXLe
zY=#JOsW}j&T5vdVr2gbJa%u@FlhlDCr&etEMe4shzOA8@OBz50sW!Y14%$M>q;^2z
zpgkLy)abz)9`^Ci==s#IiOJobQ=hRm)>FM-km|By$x40V>9b++=lkoEX-Pv(bqS*e
zE-iDlJ%2loJT<yCxYQBTt|z{<J7ii)*>5H#e{(J)Q@2*tdv4(oOie_POC5k9)se#q
zQk{5>oH|3wq%J^_lNKB5AE~ZT$|Y@}f>bx&2M66DWl|5IaL|(tOse~prc0V_&Y=gp
zN)1OGz%|O0Ny*xqa>=XM-1bPbd|DN@di<nD1#~4PaL{PX96ZqRLatL-F-@MCAo;d0
z!<}E_o_*R?(7lX{6~=ceNd3}|xjVvfbE6j`$fe#ukm|$X1gXBfMo#@8Ws(k1<kX)H
z^+^qYQZ5YyDo73DeQ+=sQYPsFg#$e{FsZN0`xx9$D58-&rCD=c7ts`hKZ?hkEuulq
zb!2l<YtrM@)sd+xxQP+dFn%lMJj?UqeDBUGB!el7N(ZA$^+&_d2K%jx=|P|3?n7{Y
zAhRIAWLpyM0}Md~xug#Ssi7QBkQ&Bo<YWLTlZFFDP9xY*pVUYw<<cmig4Af<2M1#y
z(R2|g91t6rRFT!C#i4la-1teI>4aUyq+`7G$Pc4pDm2c!)fjzNbb9OAw6HBC_sww$
zA8ZO}*n$L)9vB4P&B6Rg=LtoWqm!`G9*>7rSzUa(jf&}9qm8CY<3VaHBA|mK5TwRK
zC~l78Dr^E?n@L7M#hHF0UT~_DAZ5~IHvKx&PhrEqTz@;$PlZx0{RLETHRgSAFbz^B
zO$Q1GGuXgfbDs1qtzCPM9M0xh?H+lb));rK)x8@*gN}E*v(XrTt{A%Q--hX~N45mV
zEjy7<OIKGPEW}5vXw$Cb#O3(^!6kHG>&X~OdG&>*a!eH!53^s=J`rd7nTQ~lW&y!<
zHir`jP{wQIGzU^9%>{~_=CPqZsrgXKr3F9*DW2}d4oo0r(juU6u$T=@YP7u5Kd1Hm
z$+m;0)lzSN+Tx{Ml(;^a<U^X81l^0G0=0wt#-NTlLmK_WdPh3F)YGpu^+T^aZE1U7
zG|WlREgjW-Tp7uK#&otiSWXA@ryltS{VgO@M375nK#*F(;RLCryhcvTAZ5~WpvcLb
z4Zq%%|MksX0i|482~?0;#rxpE0#YWe1_}pj*ubPr<D!FhtjQyjk2kt(a>}E5KH8n-
zngui}%DJmUY!SUbkYGG*MG4)XIsZq%iBg<5TGm~%EhEpu4NkkzV7udIuSjz=ytzK1
z$%dz`%PIPJY?g0JkXnlfa>)`1QtLRJAoVw|k<)renX~~Ya<XEBAjPcsbT4g$QZ8)*
zDo9DZ4-Pg%$|N~ZIM~7lCbciu>WGV}7Y&bGAN4xTi{6*-wAj|tmxi_|%ukOB!4W(=
z<M`Al(rN3vU3K&uay+_U%dkTlUHBO<-})+-u6@tzYyYH}JmtMVZTf^qyaPX6)#!mU
z{Xd8xm$m{yY8!_Wq^x<3oNORvk}Xi=w4Dw0N!dXumv#UZq;~Q?IM@X#lXe4zgFS3u
zQbQk_Em`g8K)#06dsF5(Qk&VvVMgPRQAEh($?vMJQ0ug{f41BeOh%Q*`ghnJM^7%;
z=dEx}q|KU(tCB`!(&?KcRuAe?Ktq+U&N+yFPWw8kyxoEtImsRo<dOprq#QY%Ahnm*
z$Y~#>Oxh0=IUQg_{kPviDCN>2pn}w4-UkOqAklIWC>$JR1C#PMn6usB&LMK!)H`4M
z`eAysyvLWv`ll$j`Rm;-E3VVgpUDTFREAN?roH{0+q@#zJ<mG?%y>hU^NtkijZ7p@
zy`9hZOu!!#_dfTp(#xcrrRR4Xw#_8z7$Ts5BM_ucK-8<{&>JS5gw{-Q1}aXCr|^Ok
zc7eq72R8jWHJ)a}zl?u7HJ*VI0~>$}#%FmS9GrucN#}vWfjb+R@%+s$ch;PIOdkDL
zMcwp$Oja)2x7%h!kwf8+>?+r2Dsu{JqBJUoKHK&>v>+&slIPjPS;?{~-&o$TZ*UQv
z+1P2SUzZBX)#-O-t5YR)K^Mcl-uS(E0TJYq2M~-eayW4qFYy{VU53Q;2cXF5DjVvP
z@`Mrt8-NN@*LWWsctgr0AE0n>oefMX`b=Q)mhQ3SvCz~hM;=R|X2(WmHF-rVd#Bww
zyDfnhMC%?jd7eTkQ->%SX=Tyl`-{ugDi_j?DwDuDozUj<!|}xxqblgb)BgJ|=vI<)
zW`RY!i6C_Y5#-WMAV~RgI6>+buaT1<B&I(AMNa-~s81>YN(^iODo6$LJ~#-1lu37i
z!a*<_n3R)siQ)O(K6JbJg4Ct+u2Wpv<nYY{11RxCuk6UKkEzG8_#5+=yrJ-pZ94n*
z&7lKk5j|F+@@w%N2lJylDrm{X*{hzOsiFlMs@A#c)pX#}pl6#p;JPye5#-W6AV}Tk
zaDr4QuaVOONKAhKiku#@p+2c_C^4`Bs2~-=``{oF5-k^j!a)=pnAGf!8EMaKV`=%|
zXYMa=#ZvNl!<qSS<EW%%7ulIMZzyWP%Zlg|2~;~~dcV`*$>ei6W6~y6c)yvW5x(L$
z#y`h()_Xa>nu5IVwO{X!Qlc>$g{?1x)DuKN|3)B4MRPbo>KU(*({o6f6ay4F#j>G3
zsTWYnr8uC1)Jxt62d^N}auFyTykP^A+FBlK*y7v~8W7=8nQQDs#>T$qJ^wyQX^Y!+
zxgFp^8~z&LyP?^2y4Q8kqvhR$$>B$Br^0?OXzeB!rLl+8$h5yvvgN=+Don|mHCt9f
zKI#6JH=Oamb37uTe<Kj25;>e8mBef0lng18Qh*|-R5tue>bD#BX;8|gbfAJ%2JeG|
zOh~j`1PTY)Y+zEu!vhO0Xy(#@{87duo8^+L`S6DyY;wq5tyO&2`&pFX{xrmGI3}ZY
z8e^wknNC?>jJLX)W|E%e=uEFh*>uXQclo|x)OJi+UH%uQk*Dmow{|thHCPTJpnoF}
zr1BsXzZdm}NcnheCKUh`C&fa%;8cqsG5di{zfOuJY^Z+#OQFQL2B3m#Iq!pm3P@BQ
z0EL4pHZa$tb<aJk+nu416szG+JGhZm^kNhBu2*P!LCE>jcK4{wyKO7J1;&v^`Q))B
z{d39fqHf&wA4N1lCA#wxr!w+UcKB;%=SmuHa4AVEuad&*&NuYIgbb+~5m3ee1lM;Q
zP8`7ZyhctRATj#^C~~T0Lw!;op~Scbpn}wA-UkO?AZ5~5pm6Yw4NR)plXiNim4nD_
z{OsPE*Mi9R?9AFl!C_Q({jVPvD3<Pg`*baCZ4$Yk8rATLel8jA@YrE<q=bSthj;z5
zxPoM7Z+>aJql&_3s4mJ>s-{0r%S_fb#YgQsBFLp5K#=;$;RLBVULz-^CY;kBK#@}e
zEXBL>zl&&PC^4=9s36sd_rZY*Bswkvg@eXyU{Wzd9j0~K>qoX*t}MRX={8-}&a)`^
z8BE>ZX4M!Tj-i`zE$@0brI1m#FwI?N`E>NEan-tu1vKlH*BOh$h2)^oPbtH!h%%!d
zFPVI;h_d@1yt2s`ze1WI0@^nMK}wCo2~y2?jhvc8V)g@2<n$*S>VMQ)LWyw=Kn1B*
zyblgqL!#p%P&jDA1}3#}o_d^~`vWqMZ!>J<pfKu^aW3q0cm%cAzx-~vMI5zMap-+r
zN+hdi#>Sa5GHAHthDn)sv+3Aa{Xfm0<&xiu=Wi`6Fi>dO>}{Qo6_Um2`we2f@fU4d
zL_qsSAV{_6aDtR3uaQ#+NX&i!ikv#Jp+2e3P-0vIP(ezI_rXC|NOW8T3J2ZTz@&PF
zR3)8!lR=hgN$uqY8C0xUTd@3QCM}=8`pxQ@Xk|IJTrF}~4o&DC^!kKv4!yU@cR2qj
zm)=-??K%rxA&-YEJ-XVyfLu)7+_TRV&}m)A*3MVZx4AnapnW3{q<TX9_N>k@U8xr|
zW>RmU;tbLUFF4u0kTR(sn|_@^bl6b;!`2^4xikQ%U^|fa!NDL%)E)qZ16?*S+ie40
z+p3K`K_6NgHtE^!Bo*}7>v6>P0&Pt<Drp%ML^;LnOgtl>Q^QH!D<6ExB(wOfM=qol
zQclP}`y<Ad(WE(h4ot%c#!H_n2Da{9MUP3nUyuQ=IQ0+#g$zKj)#q^H5Dw)vavBCH
zlMH|&r{Qd<Pih2|a%m(`L24B5gM-nKs67A*2Zn6;_rrGZqMEAebz14_;4o_H4RUyT
ztaCwD5Ivs0Qm@3|1&uhf`Nv{ZKu<LtAgOxe7QtHg_pR_}r|izsVdpF=>Cdmt?rA^5
z<nEVC{<_npn)c3&8k*S?XAnXJ6fyuoY8;0Xq{j0aIZc4X=m((4X(AiylbQr2rZoT+
zq^9sbIG75F+5@0)V9W+4Wv9`?FKt->`hSn9wK^R@FZ+$?_h5Ap4bC=iWz{wu&BS*6
zb5ipuY1<cgr7wR?E1q?oHmhS6JusY9Jdg^>?&zJe(f7+}_2-J4c{okBc38jK#};?2
zrXd0f8Gs-)gTo0@GkJ}iW<g@~15o57W5cfx+kdCUIZ$F+15iO~9`A#L`H-kR015{S
z*}$X*?RY-;lX)a%dxfm->kvt;>@%k|>4Kqfq0Z|&KYm6R|8cZWuZg4LZJ$bSH;kvt
zHS3<u8JkPDhCdAX_^p_%%l+QD7gmscQgGZzG?aPpE%tWbL3sO_AOZ>*fFQM)!wFKR
zyhcuDkQn^{6ge$rgCM0Syj%t)rZoT+q|A9A9ISvu?Ez3YSj7e=C6CbCGZ24q(PRCh
zzP<2VTW#vQ^4%F^zG{R~9_FYu@!xW1$opJ6c=Jb~(-Pdpu?jN(nO{Ql?=3jg77v*w
zwVv2G@o5zu8`kt^3VO0=d~CJt&<wo&ED!;O3_y@t!{G#}wY)}7mXH|z02Dd>&4yp3
z{wus(4<)8G02QRHcpn^Wgp^5}fWm>q1}2r;)h_>5NFf~@KhgLnX6Y1o?o6AD+i>2)
z$7d%GMKSZp-QSmMl#xkBozh$kfm1QtuhA`|l01XyV^Z^KYBkOOWB5@#<E$@h>gZEV
z{?k5BUV8(iHY0*uk^@0%3x^Y={^2!p+6sx$4?vNVH5=-avVjuQ8h{E?+j$=x*g>M_
zB2YNk$p$8MeDkuUJ;x?d)VAXV_9K(%(5Jfijh&MyqUu2JBfFAl=R&Q{6NV;J8~wPh
zE4!!AqPm8MI`vJaJ%2A$H?~Nk#SVXW4X#L}u)q_?PJB<KZ3E9vs#Zcd&Mrhi14kf8
z?cs2Ols&JJlLI84dH{-?_OhWqseMr5sRy8f)B)ZH2L~bX)B{jBILrnn<umHgGyAea
zS~OsE+hus-@Q~}yMk@I0%Ufr?=XMO)vTmhX@vIa63R6+CaEvabjaTERwS;}fC7sK6
z&5LN=@@2Dq78jAdO+fVO0r;uml7HH+H-2gyK?FSY00gO{98Qor#%tts91>4G07XtG
z*-)R9Gn9Dh0jMD5!u#OB6%tQ90EL4yY+zDNk}NKHVW?L4k6s6zKbKHmRd7N;9iF7E
zyr?=6LrLYH`#Vg{ET-}}Pw&c_BGS)mH_-1~DfLlyb1IORlX7X(<Ck!vID4ToVA?D+
z><f!jowX3ZblngEPdxxZ>Kw#x-^xRErSs63N$x<!xAFzN;AA}@G5Uc`zrK|(vEg5~
zzuoM(3?-&D02OSn@;*55ghc59P&l~824<UY*r><mlIyhT@(-WkIydR6m463s-Mf@B
z;zG^IPO&soC9`V4%~bl-t#9<{)B<YS`-a7!#--HrLgUUC?90i+IBCis%PYy^ZrEb8
zuV_E`<&x@i6`U5m5dn1!K(M{e;lv@l!E5An6B45zfFh?`Y^YDl4@yjH04hlN^FBBT
zfJEs5P&f!=1Cv_V@6Ne?q3-0OJKUtm{{oTakg)MbJV|fqncjE%-J=?<mM`seVyJ~x
zyVMS-2u{?J2X(uXMs9Nrk`o%Gll>aM1TS>CtQ~y4`TUU?G_}W>lSd}uvEU#?Kpg`R
zq=Gq|AQi%E<a7@bqaT1Gr%*Q3C-neIOltrtNIm3za1aiO(gUDy5Wxl}75n1G;f9{L
zC8~4mf|L6bsvBrJ<-+OLq^5NKnkH_E+G>8e;~b3#SiU^XUp}sqCSLVA>xL%<hTLp$
z$o^|N?eyvE(dt7Pb$Zhw^}9(K>Ey^${r<u&(?~=>9Rm=gqBxu&^@P{R=_#a4iUx|D
zp0VNAhwZ;_<>ye!r5K=sR4ng<gBOq}Jpc*^FWJDP8aia&(TjaTtt|$|hqjNVw&BZS
zPKPDXwCkVl7b+K#t<{Y2{kK(6uct<D=Ad6i+p!;h#N<}d!Qx&fYG@d#?7Y^XOOq<H
zRB1KzV^Ae+?Bz7#DA!HBLIl(?072>vhZCgYd5xSBAZ1b_P~?=vhF_%q`vOabQZA(c
z6{J#m9~`7XqVxbL9AvP8N!^}4<DRx<Dh=C`dBxr$m6jU?o^LxZl|1y*>T=IOKjG^G
zt!8PIZ?V@l#U+jEtcMSLJvEi?&U^ROuT=^)wVc<;Cpn4ECygx2-;SY3I*069ZpFDi
z6A@6y00gOQ4kt+E@ESSgLdv8(pvWno4dO^CergmzDVGX?3Q|S94-SeUWl{-HI4ETU
zlhW>zQ@p?MCe4kt%Q4^LOPANmlB5KG^8Q+I`0nV3<fT6-Zl!rN8I|^QaaN9}(Fc#H
zht<VXna_c7R<cC$x)rLa?TT^o`{x+XdW-G}!>`_xTA>=J3=!l~IS`~OAby<|7foL}
zWboi2QYAEIQWa2fTCBzkPWCM%#y+s=*J<%R8|r`9K0t{{4L}9kTHXf-A0bh602B^B
zvw_)m$+!I!sp&%@?oRfmD|~3}yFcfR&AUm?&If~(mG4lmv58hQKZa06=WQJ|mOi2l
zSG0akb&jLcIn!Eof0Rm-rRp<A4>1g*B0jkno&cYv?<wsbjr&z!5CIhoK(PJB;lv^Q
z&THiK0}^8&fFh?lHq<Aj)Rf))0aTD`fEVmQ84^VYK;fVf8<^BJ4WC3u)z_rHMklhx
z%-7@;u~9O!NubW*4aR7r;l+-XvcA*S=hD?3&E^+j=6Fl>&&|iX;X`<Ivg@Q<#nkId
zgBEY&OK9(mE84lUOG$hA*pd2XcwAit5m3PZ1gXXxPLOKCYvj}vQYNVZMNZAwP@hzD
zDCJTMpn}w&yblgqLZavZC>*q61C#PjE}8o<z?057P4wK3f0Qa5BTBB{y+y9aECzYk
zJRrUNiUk{uo>9YIqX*}!BvE7CWs5VM3n*;fy;(t_s1X=pf4?234Loh8@go|4V7@nX
z`MWL;Utp~f0Tm2DkZQx>1gW;XMo#S@Wm0>f$Vrn8^+|PrQZ97_DoAzWeQ?kj5=93<
z;XsQGOloe|_TTm|xJ8q^Yk$7$<VP-J)t>C~2_&<k>D_BHBgwzX%T<FqyrR%BHS4#-
z(?}`K=7whbd@}q{wLI}?F*OQGY&d9QIo&M!bEDj*f@)%?Tv>b)S72Qc0Tm2Dkm|<a
z1gY-4Mov8-Wl~R|$f*|_>XYgXrCjO*RFLY+`{1A-B#I7z!a;vFFe&>RZP%<@h1Rf9
z`Nge1#nW{46wTm+X*6kmWjCcm1vJ9){WRANn2q^|Q}+Q`6%@EecV5%EmE?2a@T0w6
zmGpOJr!8SMl~iN<&}Vk5DjK?B-lA+Y;Fbm;0xB4QAT@}?2~vZ3jhu8LWs)9H<TQj0
z|C0Lc9<V-?a%m_~L24N9g98Id6deGCgAr_CQla66b8S=bLD$%OI1IhbJZjgwTK-i)
zQ3ob=y!#`M9A1C$U09M!9u5<hPf5!myXTFT&PSb1*uy=VM^ENbZPee9&A;Yi>fZ~s
z10D0IW5BxN%Fd{~9Ek|1U;u*DXbvYxjo~$NGK7>#1Qa=qWkY>Z<DisF<ADlN6L=pS
z7(vRUi9q3C5*wIQp2lw1Ms+W!-Norg9S!5i{q4^APp-eB7~fv=(()3>?Ztw2H!!8I
zMF*`7Zl)PzJS0Y$zNV2?pSCxvT~q0P;i3GZ=BPw3x;5nLp=3I>`K0yyJ!lX)84={t
z6d*`Vg{c3ZpPuv=v}TeqP;tJWh8LXhbV$s8VAHSj{Y*C0e{ar$65|?x3dS<t2M2Q?
zQF;Ir4(73e88@msddNlBf#QdHZ1tMqNXBO-<(A}~Ahkl9u)U?eG`}RxBd6D6nw&O4
z%h*4aY6{%`krfotn>JSp)FxGsc4#w~(-~DXb5cy>2WHju#n@@nMs7Md9}!T;00iTO
z98MfY6J8^yMUa^N02DcyvY|dHGbk~x0jMCgl=s2GGDwsj0EGi{HZZA2{&7~dIsr7I
z#q-Fw`vRy@Oma<4lVGaxC|~|XC6XR=+)?5`F_yY`4E~}Vm_%1^N3PLw$)PdfK8r4*
zSf}ALkK2w3r8H=^dheJO<)|nq>_6B6l>;jf0d)*OkXpsz1Stz%Bd67nnEe10Ijv>G
zulMG^pE&;6C|N>@aScEPslRz29IS^#=>bqUuwny~nsy;`?(_aP=}UE;U+xBUVhs+{
z7$)<h++hy3*A@iQ%X6a|+<$kUBL2u;`Oqek796daX<GXV{Xo0;MNk@Pod|JP3&^MI
zvKB{tPZg1N|D{vjqLNA4hzO`-0D_dn;RLD8yhcuPNX&i!ik$vo!*8VcT1?ssCB`)X
z6{M_r9~{^~qVxbL9BgL;lN!<VQ=)FlPD+t$%e<!TBFp496mY<iypI)yev5G-17};m
zKXtt+iVD~14|zo9+U3*6zlo)mZ@q@n!X%p7{6m8=$8*W~z=8lx)TUZJ$PajkQS*`=
zBA|`|2vR#aoFKJ}*T`u%q)gfa6gk<mp*|@GDCLqPP(f-h?}LMVkSIL>3I_++z@#F^
zJwB*{ep+RIo6mJ`olf5mrUcg=%_QeP%8g!)$4te(B?BA2&!cl4zFFF=D4;R7_T3u2
zpoF^ob$R4alXBYgy7Kc3bU^Es{a1iZU!3p5y&@j!<9vS*5m3hf1gXOuPLMjnYvkkv
zDU*%@MNY@q@Qc)czwM7hDVI(F6{Jq`J~(iOM9)Q_aNxoQCN(~4_~o^haio@`-Db+_
zI6Ci|Hm9>=9Bp}icxyt7S0o>wv$w&7*JS+lhJWJMH)OW$(e2O$FX`xmi1+?R&uOdD
z@rCKxQ50;Ur4)DUF?D%m+?KXICdm~M(7+K0QfD}vAmzqu<a8ENCY=L{oX)eMJ}Gx7
z<<bS9f|LjEgM*8Z=(z|K4lc8SNwt+l^ov=XNwYdDpGcjTNsg|ER`k`%q@nAQRQKm*
zP~RafF8L*8kj2&aNyd2@<eB*6g{4O(-7f2;ddV(}B3|l+NqC$pb7{ASWrjI4c-FC~
zF_ZDQ!4*V614kf8c|s`u$ka8Eyztsgx&~C7>AmrSQ}uzw><2deI@8}^!@pdAyWw#Y
zN{nj&D!AU_eQ@9hiJ}9baNy4d=IZ&a$)MwPq2xFIp0TRO18VSbNx_q&5wvvYBiSge
z7c^}4gsd+Y5~!zj`#r0MWzyG`e+?d^T1YR|&d)T%kKS#@o;1<Ef);#_tj$Fy*6Nc6
z*Crpq6=wh<pn?Giu7Mm*9Kax6Bd5ENnEe10Ifbx6eANDT9WC9165|?x3R0oG4-Otc
zqUZo996V$LlhS#y@$6vN7;5ih@~&n3SgI~xKW%gK*EG=6wdLxUDRivAk;a**Z0i3*
zb;pXe`Sf?A5nYDAFGS~#VS)Wy<6c*8)p;8fz6~+ikXhWKoMJ3%EEBYFrVmF1R4@QR
zDuTlaQjxqyPLClm`vE9&dcp=lN^!;c6iSS104hj5<9%@O91=wbK;a;k4NS_hM``xr
z08Dd!Y<QvdlL*RQF=fp6@3?cX^00kyTr{nITD|ERuDL7+zFV%>DV|c*(}qmClS2E?
zpHaStLCN9jzKOke71G?ZQzxGNR7^>|qP_BUa1HhX5m3PZ1gV!CPLO)VYvlAA60;wG
zBByvZ)PMUWK#6e;Kn1BJ-UkQCkSIC;3J0ldU{d{zmriv1gv;n7Hci!U;)lgI4KF?a
z0BU=r_<Q{EaC|}KT+gz2MoacJ`LVZI0(H2h^0det-Thm=v%B3OkFtiX95-bYYQr8*
z?eH9=u7+4^nV`>~l!geXU;u(t28R=*GI@=hvLG@00Vs0HVZ*N@_1{S`7fOt404hl3
z^FBBzfRssvK;fW>4NR)NZo=e^C>Bjnb9Z%pkw(iB<e|}@(`b+Gmlf|V)9J*wMu8y_
z>7@RT$~AM147%Q<r%5EFNrze0J-??@Q$I&{vq$N)hEzU}^hu}p0s9i)T*H)^VnmQj
zB|wlW<#2*j8LyF3IV5I307XufY^eX4s)7>Z8h{E?Z+Raayn~cU?}5U>2R1ONgEQ>X
zKeo%IF_+R^%`gJ8{F!E<y<aZ<u)mhH;btB+4%u#TYIZ)&Tl`ZeF*Tp;)Ku@-;V-@_
zgPOGd6(v+wb@9%)pJg;{&7Izj?v&H5v8#eJBJi22K?J!}3k0c;98Qq>#B1dA84|M}
zfFh@_Y^Z;vzCnp`4L}8{AG{9^enO(-B2YL`QiFj>>3XOyn6VT8J0`~aIB!F-P;k|f
z?OiI#ZcI`7hU{{hT1C}gjmoHe=G0^F(U$6*&GtVJp!w~~H7b5*7L}8=e8`C1=!MWw
z+3?7is0!*Bwdz*vb6i;dfe2{d2m~o*4kt)8<TY|?1SykLfFdVVHq<B87)rU+1gId@
zl=s1b8YDU{0)>O-Y+zD9y@%TV^P-HZ19Kl5L}6xZj&F%HuZ-H%ev2qCK+nswllt_V
zS4;`HZXLHf;)LM(=Jd>lm?zeG@72;7C8T|^B3f&EDH+(d>)P{mDXmyre(q-y?rXF_
Q1hj7if>cY0Ra*xBAA2*GbN~PV

literal 0
HcmV?d00001

diff --git a/behaviour_overview.py b/behaviour_overview.py
index 48418b9c..24180ed2 100644
--- a/behaviour_overview.py
+++ b/behaviour_overview.py
@@ -2,6 +2,7 @@ import matplotlib.pyplot as plt
 import numpy as np
 import pandas as pd
 from one.api import ONE
+import pickle
 
 
 show = 1
@@ -21,27 +22,39 @@ def progression(data, contrasts, progression_variable='feedback', windowsize=6,
     plt.ylim(0, upper_bound)
     plt.gca().spines['top'].set_visible(False)
     plt.gca().spines['right'].set_visible(False)
-    # if title == "7 rt" or title == "7 fb":
-    #     plt.plot([271, 286], [1, 1], color='r', lw=3)
-    #     plt.plot([429, 439], [1, 1], color='r', lw=3)
-    #     plt.plot([632, 636], [1, 1], color='r', lw=3)
-    if title == "KS014, RT session 12 / 15" or title == "KS014, performance session 12 / 15":
-        plt.plot([599, 603], [1, 1], color='m', lw=6)
-        plt.plot([658, 662], [1, 1], color='m', lw=6)
-        plt.plot([703, 708], [1, 1], color='m', lw=6)
+    # if title == "KS014, RT session 12 / 15" or title == "KS014, performance session 12 / 15":
+    #     plt.plot([599, 603], [1, 1], color='m', lw=6)
+    #     plt.plot([658, 662], [1, 1], color='m', lw=6)
+    #     plt.plot([703, 708], [1, 1], color='m', lw=6)
 
     plt.title(title, size=22)
     plt.legend(fontsize=13)
     plt.ylabel(progression_variable, size=20)
     plt.xlabel('trial', size=20)
-    if title:
-        plt.savefig('./overview_figures/' + title.replace('/', '_') + '.png')
+    # if title:
+    #     plt.savefig('./overview_figures/' + title.replace('/', '_') + '.png')
 
     if show:
         plt.show()
     else:
         plt.close()
 
+    # if progression_variable == 'feedback':
+    #     means = data.groupby('signed_contrast').mean()['response']
+    #     stds = data.groupby('signed_contrast').std()['response']
+    #     plt.errorbar(means.index, means.values, stds.values)
+    #     plt.ylim(-0.2, 1.2)
+    #     plt.title(title, size=22)
+    #     plt.savefig("temp {}".format(title).replace('/', '_'))
+    #     plt.close()
+    # if progression_variable == 'rt':
+    #     means = data.groupby('signed_contrast').mean()['rt']
+    #     stds = data.groupby('signed_contrast').std()['rt']
+    #     plt.errorbar(means.index, means.values, stds.values)
+    #     plt.title(title, size=22)
+    #     plt.savefig("temp {}".format(title).replace('/', '_'))
+    #     plt.show()
+
 
 dataset_types = ['choice', 'contrastLeft', 'contrastRight', \
                  'feedbackType', 'probabilityLeft', 'response_times', \
@@ -49,7 +62,7 @@ dataset_types = ['choice', 'contrastLeft', 'contrastRight', \
 def get_df(eid):
     # dangerously many changes, check again
     try:
-        data_dict = one.load_object(eid, 'trials', attribute=dataset_types) # download_only=True messes with my code, returns path for whatever reason
+        data_dict = one.load_object(eid, 'trials') # download_only=True messes with my code, returns path for whatever reason
     except:
         print('lol')
         return None, None
@@ -82,25 +95,34 @@ exclude_eids = ['a66f1593-dafd-4982-9b66-f9554b6c86b5', 'ee40aece-cffd-4edb-a4b6
 #                      project='ibl_neuropixel_brainwide_01')
 # traj.reverse()
 
-for subject in ['KS014']:
-    eids, sess_info = one.search(subject=subject, date_range=['2015-01-01', '2022-01-01'], details=True)
+subject = 'KS014'
+eids, sess_info = one.search(subject=subject, date_range=['2015-01-01', '2022-01-01'], details=True)
 
-    start_times = [sess['date'] for sess in sess_info]
-    protocols = [sess['task_protocol'] for sess in sess_info]
+start_times = [sess['date'] for sess in sess_info]
+protocols = [sess['task_protocol'] for sess in sess_info]
 
-    eids = [x for _, x in sorted(zip(start_times, eids))]
+eids = [x for _, x in sorted(zip(start_times, eids))]
+protocols = [x for _, x in sorted(zip(start_times, protocols))]
 
 # for t in traj:
 counti = 0
-for i, eid in enumerate(eids):
+for i, (prot, eid) in enumerate(zip(protocols, eids)):
     print(i)
+    if not prot.startswith('_iblrig_tasks_trainingChoiceWorld'):
+        continue
     # eid = t['session']['id']
     df, _ = get_df(eid)
+
     if df is None:
         continue
 
     counti += 1
-    if counti != 12:
-        continue
-    progression(df, df['signed_contrast'].unique(), progression_variable='feedback', upper_bound=2, title="{}, performance session {} / 15".format(subject, counti))
-    progression(df, df['signed_contrast'].unique(), progression_variable='rt', upper_bound=4, title="{}, RT session {} / 15".format(subject, counti))
+
+    rt_data = np.zeros((len(df), 3))
+    rt_data[:, 0] = df['signed_contrast']
+    rt_data[:, 1] = df['rt']
+    rt_data[:, 2] = df['response']
+    pickle.dump(rt_data, open("./session_data/{} rt info {}".format(subject, counti), 'wb'))
+
+    progression(df, df['signed_contrast'].unique(), progression_variable='feedback', upper_bound=2, title="{} PMF {} / 15".format(subject, counti))
+    progression(df, df['signed_contrast'].unique(), progression_variable='rt', upper_bound=4, title="{} CMF {} / 15".format(subject, counti))
diff --git a/canonical_infos.json b/canonical_infos.json
index c2ca57ff..5f8da998 100644
--- a/canonical_infos.json
+++ b/canonical_infos.json
@@ -1 +1 @@
-{"SWC_023": {"seeds": ["302", "312", "304", "300", "315", "311", "308", "305", "303", "309", "306", "313", "307", "314", "301", "310"], "fit_nums": ["994", "913", "681", "816", "972", "790", "142", "230", "696", "537", "975", "773", "918", "677", "742", "745"], "chain_num": 4}, "SWC_021": {"seeds": ["415", "403", "412", "407", "409", "408", "405", "404", "410", "414", "401", "413", "402", "400", "406", "411"], "fit_nums": ["773", "615", "107", "583", "564", "354", "142", "184", "549", "185", "924", "907", "105", "531", "9", "812"], "chain_num": 9}, "ibl_witten_15": {"seeds": ["409", "410", "401", "415", "414", "403", "411", "404", "402", "405", "400", "412", "408", "407", "406", "413"], "fit_nums": ["411", "344", "496", "600", "716", "18", "527", "467", "898", "334", "309", "326", "133", "823", "740", "253"], "chain_num": 9}, "ibl_witten_13": {"seeds": ["302", "312", "313", "306", "315", "307", "311", "314", "309", "301", "308", "300", "304", "310", "303", "305"], "fit_nums": ["897", "765", "433", "641", "967", "599", "984", "259", "853", "385", "887", "619", "434", "964", "483", "891"], "chain_num": 4}, "KS016": {"seeds": ["315", "301", "309", "313", "302", "307", "303", "308", "311", "312", "314", "306", "310", "300", "305", "304"], "fit_nums": ["99", "57", "585", "32", "501", "558", "243", "413", "59", "757", "463", "172", "524", "957", "909", "292"], "chain_num": 4}, "KS003": {"seeds": ["404", "407", "413", "403", "414", "405", "400", "401", "402", "410", "415", "408", "411", "409", "406", "412"], "fit_nums": ["846", "256", "845", "945", "293", "406", "420", "109", "690", "421", "54", "866", "784", "81", "997", "665"], "chain_num": 9}, "ibl_witten_19": {"seeds": ["315", "311", "307", "314", "308", "300", "305", "301", "313", "304", "302", "310", "306", "312", "309", "303"], "fit_nums": ["179", "951", "613", "6", "623", "382", "458", "504", "406", "554", "5", "631", "746", "817", "265", "328"], "chain_num": 4}, "SWC_022": {"seeds": ["411", "403", "414", "409", "407", "412", "410", "413", "415", "404", "405", "400", "402", "401", "408", "406"], "fit_nums": ["408", "884", "62", "962", "744", "854", "635", "70", "320", "952", "8", "67", "231", "381", "536", "962"], "chain_num": 9}, "KS022": {"seeds": ["315", "300", "314", "301", "303", "302", "306", "308", "305", "310", "313", "312", "304", "307", "311", "309"], "fit_nums": ["899", "681", "37", "957", "629", "637", "375", "980", "810", "51", "759", "664", "420", "127", "259", "555"], "chain_num": 4}, "CSH_ZAD_017": {"seeds": ["401", "409", "405", "403", "415", "404", "402", "411", "410", "414", "408", "406", "413", "412", "400", "407"], "fit_nums": ["883", "803", "637", "806", "356", "804", "662", "654", "684", "350", "947", "460", "569", "976", "103", "713"], "chain_num": 9}, "CSH_ZAD_025": {"seeds": ["303", "311", "307", "312", "313", "314", "308", "315", "305", "306", "304", "302", "309", "310", "301", "300"], "fit_nums": ["581", "148", "252", "236", "581", "838", "206", "756", "449", "288", "756", "593", "733", "633", "418", "563"], "chain_num": 4}, "ibl_witten_17": {"seeds": ["406", "415", "408", "413", "402", "405", "409", "400", "414", "401", "412", "407", "404", "410", "403", "411"], "fit_nums": ["827", "797", "496", "6", "444", "823", "384", "873", "634", "27", "811", "142", "207", "322", "756", "275"], "chain_num": 9}, "ibl_witten_18": {"seeds": ["311", "310", "303", "314", "302", "309", "305", "307", "312", "300", "308", "306", "315", "313", "304", "301"], "fit_nums": ["236", "26", "838", "762", "826", "409", "496", "944", "280", "704", "930", "419", "637", "896", "876", "297"], "chain_num": 4}, "CSHL_018": {"seeds": ["302", "310", "306", "300", "314", "307", "309", "313", "311", "308", "304", "301", "312", "303", "305", "315"], "fit_nums": ["843", "817", "920", "900", "226", "36", "472", "676", "933", "453", "116", "263", "269", "897", "568", "438"], "chain_num": 4}, "GLM_Sim_06": {"seeds": ["313", "309", "302", "303", "305", "314", "300", "315", "311", "306", "304", "310", "301", "312", "308", "307"], "fit_nums": ["9", "786", "286", "280", "72", "587", "619", "708", "360", "619", "311", "189", "60", "708", "939", "733"], "chain_num": 2}, "ZM_1897": {"seeds": ["304", "308", "305", "311", "315", "314", "307", "306", "300", "303", "313", "310", "301", "312", "302", "309"], "fit_nums": ["549", "96", "368", "509", "424", "897", "287", "426", "968", "93", "725", "513", "837", "581", "989", "374"], "chain_num": 4}, "CSHL_020": {"seeds": ["305", "309", "313", "302", "314", "310", "300", "307", "315", "306", "312", "304", "311", "301", "303", "308"], "fit_nums": ["222", "306", "243", "229", "584", "471", "894", "238", "986", "660", "494", "657", "896", "459", "100", "283"], "chain_num": 4}, "CSHL054": {"seeds": ["401", "415", "409", "410", "414", "413", "407", "405", "406", "408", "411", "400", "412", "402", "403", "404"], "fit_nums": ["901", "734", "609", "459", "574", "793", "978", "66", "954", "906", "954", "111", "292", "850", "266", "967"], "chain_num": 9}, "CSHL_014": {"seeds": ["305", "311", "309", "300", "313", "310", "307", "306", "304", "312", "308", "302", "314", "303", "301", "315"], "fit_nums": ["371", "550", "166", "24", "705", "385", "870", "884", "831", "546", "404", "722", "287", "564", "613", "783"], "chain_num": 4}, "CSHL062": {"seeds": ["307", "313", "310", "303", "306", "312", "308", "305", "311", "314", "304", "302", "300", "301", "315", "309"], "fit_nums": ["846", "371", "94", "888", "499", "229", "546", "432", "71", "989", "986", "91", "935", "314", "975", "481"], "chain_num": 4}, "CSH_ZAD_001": {"seeds": ["313", "309", "311", "312", "305", "310", "315", "300", "314", "304", "301", "302", "308", "303", "306", "307"], "fit_nums": ["468", "343", "314", "544", "38", "120", "916", "170", "305", "569", "502", "496", "452", "336", "559", "572"], "chain_num": 4}, "NYU-06": {"seeds": ["314", "309", "306", "305", "312", "303", "307", "304", "300", "302", "310", "301", "315", "308", "313", "311"], "fit_nums": ["950", "862", "782", "718", "427", "645", "827", "612", "821", "834", "595", "929", "679", "668", "648", "869"], "chain_num": 4}, "KS019": {"seeds": ["404", "401", "411", "408", "400", "403", "410", "413", "402", "407", "415", "409", "406", "414", "412", "405"], "fit_nums": ["682", "4", "264", "200", "250", "267", "737", "703", "132", "855", "922", "686", "85", "176", "54", "366"], "chain_num": 9}, "CSHL049": {"seeds": ["411", "402", "414", "408", "409", "410", "413", "407", "406", "401", "404", "405", "403", "415", "400", "412"], "fit_nums": ["104", "553", "360", "824", "749", "519", "347", "228", "863", "671", "140", "883", "701", "445", "627", "898"], "chain_num": 9}, "ibl_witten_14": {"seeds": ["310", "311", "304", "306", "300", "302", "314", "313", "303", "308", "301", "309", "305", "315", "312", "307"], "fit_nums": ["563", "120", "85", "712", "277", "871", "183", "661", "505", "598", "210", "89", "310", "638", "564", "998"], "chain_num": 4}, "KS014": {"seeds": ["301", "310", "302", "312", "313", "308", "307", "303", "305", "300", "314", "306", "311", "309", "304", "315"], "fit_nums": ["668", "32", "801", "193", "269", "296", "74", "24", "270", "916", "21", "250", "342", "451", "517", "293"], "chain_num": 4, "ignore": [9, 11, 0, 1, 14, 2, 12, 13]}, "CSHL059": {"seeds": ["306", "309", "300", "304", "314", "303", "315", "311", "313", "305", "301", "307", "302", "312", "310", "308"], "fit_nums": ["821", "963", "481", "999", "986", "45", "551", "605", "701", "201", "629", "261", "972", "407", "165", "9"], "chain_num": 4}, "GLM_Sim_13": {"seeds": ["310", "303", "308", "306", "300", "312", "301", "313", "305", "311", "315", "304", "314", "309", "307", "302"], "fit_nums": ["982", "103", "742", "524", "614", "370", "926", "456", "133", "143", "302", "80", "395", "549", "579", "944"], "chain_num": 2}, "CSHL_007": {"seeds": ["314", "303", "308", "313", "301", "300", "302", "305", "315", "306", "310", "309", "311", "304", "307", "312"], "fit_nums": ["462", "703", "345", "286", "480", "313", "986", "165", "201", "102", "322", "894", "960", "438", "330", "169"], "chain_num": 4}, "CSH_ZAD_011": {"seeds": ["314", "311", "303", "300", "305", "310", "306", "301", "302", "315", "304", "309", "308", "312", "313", "307"], "fit_nums": ["320", "385", "984", "897", "315", "120", "320", "945", "475", "403", "210", "412", "695", "564", "664", "411"], "chain_num": 4}, "KS021": {"seeds": ["309", "312", "304", "310", "303", "311", "314", "302", "305", "301", "306", "300", "308", "315", "313", "307"], "fit_nums": ["874", "943", "925", "587", "55", "136", "549", "528", "349", "211", "401", "84", "225", "545", "153", "382"], "chain_num": 4}, "GLM_Sim_15": {"seeds": ["303", "312", "305", "308", "309", "302", "301", "310", "313", "315", "311", "314", "307", "306", "304", "300"], "fit_nums": ["769", "930", "328", "847", "899", "714", "144", "518", "521", "873", "914", "359", "242", "343", "45", "364"], "chain_num": 2}, "CSHL_015": {"seeds": ["301", "302", "307", "310", "309", "311", "304", "312", "300", "308", "313", "305", "314", "315", "306", "303"], "fit_nums": ["717", "705", "357", "539", "604", "971", "669", "76", "45", "413", "510", "122", "190", "821", "368", "472"], "chain_num": 4}, "ibl_witten_16": {"seeds": ["304", "313", "309", "314", "312", "307", "305", "301", "306", "310", "300", "315", "308", "311", "303", "302"], "fit_nums": ["392", "515", "696", "270", "7", "583", "880", "674", "23", "576", "579", "695", "149", "854", "184", "875"], "chain_num": 4}, "KS015": {"seeds": ["315", "305", "309", "303", "314", "310", "311", "312", "313", "300", "307", "308", "304", "301", "302", "306"], "fit_nums": ["257", "396", "387", "435", "133", "164", "403", "8", "891", "650", "111", "557", "473", "229", "842", "196"], "chain_num": 4}, "GLM_Sim_12": {"seeds": ["304", "312", "306", "303", "310", "302", "300", "305", "308", "313", "307", "311", "315", "301", "314", "309"], "fit_nums": ["971", "550", "255", "195", "952", "486", "841", "535", "559", "37", "654", "213", "864", "506", "732", "550"], "chain_num": 2}, "GLM_Sim_11": {"seeds": ["300", "312", "310", "315", "302", "313", "314", "311", "308", "303", "309", "307", "306", "304", "301", "305"], "fit_nums": ["477", "411", "34", "893", "195", "293", "603", "5", "887", "281", "956", "73", "346", "640", "532", "688"], "chain_num": 2}, "GLM_Sim_10": {"seeds": ["301", "300", "306", "305", "307", "309", "312", "314", "311", "315", "304", "313", "303", "308", "302", "310"], "fit_nums": ["391", "97", "897", "631", "239", "652", "19", "448", "807", "35", "972", "469", "280", "562", "42", "706"], "chain_num": 2}, "CSH_ZAD_026": {"seeds": ["312", "313", "308", "310", "303", "307", "302", "305", "300", "315", "306", "301", "311", "304", "314", "309"], "fit_nums": ["699", "87", "537", "628", "797", "511", "459", "770", "969", "240", "504", "948", "295", "506", "25", "378"], "chain_num": 4}, "KS023": {"seeds": ["304", "313", "306", "309", "300", "314", "302", "310", "303", "315", "307", "308", "301", "311", "305", "312"], "fit_nums": ["698", "845", "319", "734", "908", "507", "45", "499", "175", "108", "419", "443", "116", "779", "159", "231"], "chain_num": 4}, "GLM_Sim_05": {"seeds": ["301", "315", "300", "302", "305", "304", "313", "314", "311", "309", "306", "307", "308", "310", "303", "312"], "fit_nums": ["425", "231", "701", "375", "343", "902", "623", "125", "921", "637", "393", "964", "678", "930", "796", "42"], "chain_num": 2}, "CSHL061": {"seeds": ["305", "315", "304", "303", "309", "310", "302", "300", "314", "306", "311", "313", "301", "308", "307", "312"], "fit_nums": ["396", "397", "594", "911", "308", "453", "686", "552", "103", "209", "128", "892", "345", "925", "777", "396"], "chain_num": 4}, "CSHL051": {"seeds": ["303", "310", "306", "302", "309", "305", "313", "308", "300", "314", "311", "307", "312", "304", "315", "301"], "fit_nums": ["69", "186", "49", "435", "103", "910", "705", "367", "303", "474", "596", "334", "929", "796", "616", "790"], "chain_num": 4}, "GLM_Sim_14": {"seeds": ["310", "311", "309", "313", "314", "300", "302", "304", "305", "306", "307", "312", "303", "301", "315", "308"], "fit_nums": ["616", "872", "419", "106", "940", "986", "599", "704", "218", "808", "244", "825", "448", "397", "552", "316"], "chain_num": 2}, "GLM_Sim_11_trick": {"seeds": ["411", "400", "408", "409", "415", "413", "410", "412", "406", "414", "403", "404", "401", "405", "407", "402"], "fit_nums": ["95", "508", "886", "384", "822", "969", "525", "382", "489", "436", "344", "537", "251", "223", "458", "401"], "chain_num": 2, "ignore": [10, 12, 4, 1, 0, 3, 2, 6]}, "GLM_Sim_16": {"seeds": ["302", "311", "303", "307", "313", "308", "309", "300", "305", "315", "304", "310", "312", "301", "314", "306"], "fit_nums": ["914", "377", "173", "583", "870", "456", "611", "697", "13", "713", "159", "248", "617", "37", "770", "780"], "chain_num": 2}, "ZM_3003": {"seeds": ["300", "304", "307", "312", "305", "310", "311", "314", "303", "308", "313", "301", "315", "309", "306", "302"], "fit_nums": ["603", "620", "657", "735", "357", "390", "119", "33", "62", "617", "209", "810", "688", "21", "744", "426"], "chain_num": 4}, "CSH_ZAD_022": {"seeds": ["305", "310", "311", "315", "303", "312", "314", "313", "307", "302", "300", "304", "301", "308", "306", "309"], "fit_nums": ["143", "946", "596", "203", "576", "403", "900", "65", "478", "325", "282", "513", "460", "42", "161", "970"], "chain_num": 4}, "GLM_Sim_07": {"seeds": ["300", "309", "302", "304", "305", "312", "301", "311", "315", "314", "308", "307", "303", "310", "306", "313"], "fit_nums": ["724", "701", "118", "230", "648", "426", "689", "114", "832", "731", "592", "519", "559", "938", "672", "144"], "chain_num": 1}, "KS017": {"seeds": ["311", "310", "306", "309", "303", "302", "308", "300", "313", "301", "314", "307", "315", "304", "312", "305"], "fit_nums": ["97", "281", "808", "443", "352", "890", "703", "468", "780", "708", "674", "27", "345", "23", "939", "457"], "chain_num": 4}, "GLM_Sim_11_sub": {"seeds": ["410", "414", "413", "404", "409", "415", "406", "408", "402", "411", "400", "405", "403", "407", "412", "401"], "fit_nums": ["830", "577", "701", "468", "929", "374", "954", "749", "937", "488", "873", "416", "612", "792", "461", "488"], "chain_num": 2}}
\ No newline at end of file
+{"SWC_023": {"seeds": ["302", "312", "304", "300", "315", "311", "308", "305", "303", "309", "306", "313", "307", "314", "301", "310"], "fit_nums": ["994", "913", "681", "816", "972", "790", "142", "230", "696", "537", "975", "773", "918", "677", "742", "745"], "chain_num": 4, "ignore": [12, 1, 15, 14, 8, 6, 4, 10]}, "SWC_021": {"seeds": ["415", "403", "412", "407", "409", "408", "405", "404", "410", "414", "401", "413", "402", "400", "406", "411"], "fit_nums": ["773", "615", "107", "583", "564", "354", "142", "184", "549", "185", "924", "907", "105", "531", "9", "812"], "chain_num": 9, "ignore": [14, 12, 0, 10, 9, 4, 5, 1]}, "ibl_witten_15": {"seeds": ["409", "410", "401", "415", "414", "403", "411", "404", "402", "405", "400", "412", "408", "407", "406", "413"], "fit_nums": ["411", "344", "496", "600", "716", "18", "527", "467", "898", "334", "309", "326", "133", "823", "740", "253"], "chain_num": 9, "ignore": [14, 13, 8, 4, 5, 12, 11, 9]}, "ibl_witten_13": {"seeds": ["302", "312", "313", "306", "315", "307", "311", "314", "309", "301", "308", "300", "304", "310", "303", "305"], "fit_nums": ["897", "765", "433", "641", "967", "599", "984", "259", "853", "385", "887", "619", "434", "964", "483", "891"], "chain_num": 4, "ignore": [3, 5, 15, 0, 2, 12, 11, 10]}, "KS016": {"seeds": ["315", "301", "309", "313", "302", "307", "303", "308", "311", "312", "314", "306", "310", "300", "305", "304"], "fit_nums": ["99", "57", "585", "32", "501", "558", "243", "413", "59", "757", "463", "172", "524", "957", "909", "292"], "chain_num": 4, "ignore": [0, 2, 14, 12, 1, 7, 11, 6]}, "KS003": {"seeds": ["404", "407", "413", "403", "414", "405", "400", "401", "402", "410", "415", "408", "411", "409", "406", "412"], "fit_nums": ["846", "256", "845", "945", "293", "406", "420", "109", "690", "421", "54", "866", "784", "81", "997", "665"], "chain_num": 9, "ignore": [8, 15, 0, 13, 7, 12, 11, 1]}, "ibl_witten_19": {"seeds": ["315", "311", "307", "314", "308", "300", "305", "301", "313", "304", "302", "310", "306", "312", "309", "303"], "fit_nums": ["179", "951", "613", "6", "623", "382", "458", "504", "406", "554", "5", "631", "746", "817", "265", "328"], "chain_num": 4, "ignore": [13, 4, 10, 9, 2, 1, 3, 6]}, "SWC_022": {"seeds": ["411", "403", "414", "409", "407", "412", "410", "413", "415", "404", "405", "400", "402", "401", "408", "406"], "fit_nums": ["408", "884", "62", "962", "744", "854", "635", "70", "320", "952", "8", "67", "231", "381", "536", "962"], "chain_num": 9, "ignore": [4, 8, 7, 9, 1, 2, 10, 6]}, "KS022": {"seeds": ["315", "300", "314", "301", "303", "302", "306", "308", "305", "310", "313", "312", "304", "307", "311", "309"], "fit_nums": ["899", "681", "37", "957", "629", "637", "375", "980", "810", "51", "759", "664", "420", "127", "259", "555"], "chain_num": 4, "ignore": [10, 1, 0, 13, 5, 9, 12, 3]}, "CSH_ZAD_017": {"seeds": ["401", "409", "405", "403", "415", "404", "402", "411", "410", "414", "408", "406", "413", "412", "400", "407"], "fit_nums": ["883", "803", "637", "806", "356", "804", "662", "654", "684", "350", "947", "460", "569", "976", "103", "713"], "chain_num": 9, "ignore": [3, 4, 6, 7, 5, 0, 15, 12]}, "CSH_ZAD_025": {"seeds": ["303", "311", "307", "312", "313", "314", "308", "315", "305", "306", "304", "302", "309", "310", "301", "300"], "fit_nums": ["581", "148", "252", "236", "581", "838", "206", "756", "449", "288", "756", "593", "733", "633", "418", "563"], "chain_num": 4, "ignore": [8, 10, 13, 5, 12, 9, 7, 1]}, "ibl_witten_17": {"seeds": ["406", "415", "408", "413", "402", "405", "409", "400", "414", "401", "412", "407", "404", "410", "403", "411"], "fit_nums": ["827", "797", "496", "6", "444", "823", "384", "873", "634", "27", "811", "142", "207", "322", "756", "275"], "chain_num": 9, "ignore": [9, 0, 1, 7, 11, 3, 10, 8]}, "ibl_witten_18": {"seeds": ["311", "310", "303", "314", "302", "309", "305", "307", "312", "300", "308", "306", "315", "313", "304", "301"], "fit_nums": ["236", "26", "838", "762", "826", "409", "496", "944", "280", "704", "930", "419", "637", "896", "876", "297"], "chain_num": 4, "ignore": [11, 0, 4, 2, 12, 13, 8, 3]}, "CSHL_018": {"seeds": ["302", "310", "306", "300", "314", "307", "309", "313", "311", "308", "304", "301", "312", "303", "305", "315"], "fit_nums": ["843", "817", "920", "900", "226", "36", "472", "676", "933", "453", "116", "263", "269", "897", "568", "438"], "chain_num": 4, "ignore": [15, 4, 8, 0, 5, 10, 12, 11]}, "GLM_Sim_06": {"seeds": ["313", "309", "302", "303", "305", "314", "300", "315", "311", "306", "304", "310", "301", "312", "308", "307"], "fit_nums": ["9", "786", "286", "280", "72", "587", "619", "708", "360", "619", "311", "189", "60", "708", "939", "733"], "chain_num": 2, "ignore": [15, 9, 8, 14, 1, 12, 10, 3]}, "ZM_1897": {"seeds": ["304", "308", "305", "311", "315", "314", "307", "306", "300", "303", "313", "310", "301", "312", "302", "309"], "fit_nums": ["549", "96", "368", "509", "424", "897", "287", "426", "968", "93", "725", "513", "837", "581", "989", "374"], "chain_num": 4, "ignore": [0, 14, 5, 8, 7, 11, 13, 10]}, "CSHL_020": {"seeds": ["305", "309", "313", "302", "314", "310", "300", "307", "315", "306", "312", "304", "311", "301", "303", "308"], "fit_nums": ["222", "306", "243", "229", "584", "471", "894", "238", "986", "660", "494", "657", "896", "459", "100", "283"], "chain_num": 4, "ignore": [6, 5, 9, 15, 0, 8, 4, 13]}, "CSHL054": {"seeds": ["401", "415", "409", "410", "414", "413", "407", "405", "406", "408", "411", "400", "412", "402", "403", "404"], "fit_nums": ["901", "734", "609", "459", "574", "793", "978", "66", "954", "906", "954", "111", "292", "850", "266", "967"], "chain_num": 9, "ignore": [5, 12, 7, 10, 11, 2, 6, 4]}, "CSHL_014": {"seeds": ["305", "311", "309", "300", "313", "310", "307", "306", "304", "312", "308", "302", "314", "303", "301", "315"], "fit_nums": ["371", "550", "166", "24", "705", "385", "870", "884", "831", "546", "404", "722", "287", "564", "613", "783"], "chain_num": 4, "ignore": [15, 0, 3, 4, 7, 6, 1, 11]}, "CSHL062": {"seeds": ["307", "313", "310", "303", "306", "312", "308", "305", "311", "314", "304", "302", "300", "301", "315", "309"], "fit_nums": ["846", "371", "94", "888", "499", "229", "546", "432", "71", "989", "986", "91", "935", "314", "975", "481"], "chain_num": 4, "ignore": [14, 6, 3, 11, 15, 13, 4, 12]}, "CSH_ZAD_001": {"seeds": ["313", "309", "311", "312", "305", "310", "315", "300", "314", "304", "301", "302", "308", "303", "306", "307"], "fit_nums": ["468", "343", "314", "544", "38", "120", "916", "170", "305", "569", "502", "496", "452", "336", "559", "572"], "chain_num": 4, "ignore": [12, 8, 5, 1, 9, 3, 13, 15]}, "NYU-06": {"seeds": ["314", "309", "306", "305", "312", "303", "307", "304", "300", "302", "310", "301", "315", "308", "313", "311"], "fit_nums": ["950", "862", "782", "718", "427", "645", "827", "612", "821", "834", "595", "929", "679", "668", "648", "869"], "chain_num": 4, "ignore": [8, 2, 7, 12, 3, 4, 13, 11]}, "KS019": {"seeds": ["404", "401", "411", "408", "400", "403", "410", "413", "402", "407", "415", "409", "406", "414", "412", "405"], "fit_nums": ["682", "4", "264", "200", "250", "267", "737", "703", "132", "855", "922", "686", "85", "176", "54", "366"], "chain_num": 9, "ignore": [12, 14, 1, 2, 4, 7, 10, 15]}, "CSHL049": {"seeds": ["411", "402", "414", "408", "409", "410", "413", "407", "406", "401", "404", "405", "403", "415", "400", "412"], "fit_nums": ["104", "553", "360", "824", "749", "519", "347", "228", "863", "671", "140", "883", "701", "445", "627", "898"], "chain_num": 9, "ignore": [10, 11, 6, 7, 12, 13, 1, 8]}, "ibl_witten_14": {"seeds": ["310", "311", "304", "306", "300", "302", "314", "313", "303", "308", "301", "309", "305", "315", "312", "307"], "fit_nums": ["563", "120", "85", "712", "277", "871", "183", "661", "505", "598", "210", "89", "310", "638", "564", "998"], "chain_num": 4, "ignore": [11, 14, 6, 13, 5, 12, 15, 8]}, "KS014": {"seeds": ["301", "310", "302", "312", "313", "308", "307", "303", "305", "300", "314", "306", "311", "309", "304", "315"], "fit_nums": ["668", "32", "801", "193", "269", "296", "74", "24", "270", "916", "21", "250", "342", "451", "517", "293"], "chain_num": 4, "ignore": [9, 11, 0, 1, 14, 2, 12, 13]}, "CSHL059": {"seeds": ["306", "309", "300", "304", "314", "303", "315", "311", "313", "305", "301", "307", "302", "312", "310", "308"], "fit_nums": ["821", "963", "481", "999", "986", "45", "551", "605", "701", "201", "629", "261", "972", "407", "165", "9"], "chain_num": 4, "ignore": [9, 3, 5, 15, 6, 10, 2, 1]}, "GLM_Sim_13": {"seeds": ["310", "303", "308", "306", "300", "312", "301", "313", "305", "311", "315", "304", "314", "309", "307", "302"], "fit_nums": ["982", "103", "742", "524", "614", "370", "926", "456", "133", "143", "302", "80", "395", "549", "579", "944"], "chain_num": 2, "ignore": [12, 4, 11, 6, 7, 14, 0, 1]}, "CSHL_007": {"seeds": ["314", "303", "308", "313", "301", "300", "302", "305", "315", "306", "310", "309", "311", "304", "307", "312"], "fit_nums": ["462", "703", "345", "286", "480", "313", "986", "165", "201", "102", "322", "894", "960", "438", "330", "169"], "chain_num": 4, "ignore": [3, 12, 4, 5, 2, 0, 13, 1]}, "CSH_ZAD_011": {"seeds": ["314", "311", "303", "300", "305", "310", "306", "301", "302", "315", "304", "309", "308", "312", "313", "307"], "fit_nums": ["320", "385", "984", "897", "315", "120", "320", "945", "475", "403", "210", "412", "695", "564", "664", "411"], "chain_num": 4, "ignore": [0, 2, 14, 11, 7, 10, 13, 9]}, "KS021": {"seeds": ["309", "312", "304", "310", "303", "311", "314", "302", "305", "301", "306", "300", "308", "315", "313", "307"], "fit_nums": ["874", "943", "925", "587", "55", "136", "549", "528", "349", "211", "401", "84", "225", "545", "153", "382"], "chain_num": 4, "ignore": [11, 12, 0, 8, 2, 14, 5, 1]}, "GLM_Sim_15": {"seeds": ["303", "312", "305", "308", "309", "302", "301", "310", "313", "315", "311", "314", "307", "306", "304", "300"], "fit_nums": ["769", "930", "328", "847", "899", "714", "144", "518", "521", "873", "914", "359", "242", "343", "45", "364"], "chain_num": 2, "ignore": [8, 1, 0, 3, 2, 5, 10, 4]}, "CSHL_015": {"seeds": ["301", "302", "307", "310", "309", "311", "304", "312", "300", "308", "313", "305", "314", "315", "306", "303"], "fit_nums": ["717", "705", "357", "539", "604", "971", "669", "76", "45", "413", "510", "122", "190", "821", "368", "472"], "chain_num": 4, "ignore": [7, 6, 10, 2, 15, 13, 1, 3]}, "ibl_witten_16": {"seeds": ["304", "313", "309", "314", "312", "307", "305", "301", "306", "310", "300", "315", "308", "311", "303", "302"], "fit_nums": ["392", "515", "696", "270", "7", "583", "880", "674", "23", "576", "579", "695", "149", "854", "184", "875"], "chain_num": 4, "ignore": [3, 12, 2, 6, 10, 14, 4, 1]}, "KS015": {"seeds": ["315", "305", "309", "303", "314", "310", "311", "312", "313", "300", "307", "308", "304", "301", "302", "306"], "fit_nums": ["257", "396", "387", "435", "133", "164", "403", "8", "891", "650", "111", "557", "473", "229", "842", "196"], "chain_num": 4, "ignore": [7, 8, 0, 10, 2, 3, 12, 9]}, "GLM_Sim_12": {"seeds": ["304", "312", "306", "303", "310", "302", "300", "305", "308", "313", "307", "311", "315", "301", "314", "309"], "fit_nums": ["971", "550", "255", "195", "952", "486", "841", "535", "559", "37", "654", "213", "864", "506", "732", "550"], "chain_num": 2, "ignore": [0, 7, 15, 14, 3, 10, 11, 13]}, "GLM_Sim_11": {"seeds": ["300", "312", "310", "315", "302", "313", "314", "311", "308", "303", "309", "307", "306", "304", "301", "305"], "fit_nums": ["477", "411", "34", "893", "195", "293", "603", "5", "887", "281", "956", "73", "346", "640", "532", "688"], "chain_num": 2}, "GLM_Sim_10": {"seeds": ["301", "300", "306", "305", "307", "309", "312", "314", "311", "315", "304", "313", "303", "308", "302", "310"], "fit_nums": ["391", "97", "897", "631", "239", "652", "19", "448", "807", "35", "972", "469", "280", "562", "42", "706"], "chain_num": 2, "ignore": [1, 9, 15, 3, 13, 12, 7, 11]}, "CSH_ZAD_026": {"seeds": ["312", "313", "308", "310", "303", "307", "302", "305", "300", "315", "306", "301", "311", "304", "314", "309"], "fit_nums": ["699", "87", "537", "628", "797", "511", "459", "770", "969", "240", "504", "948", "295", "506", "25", "378"], "chain_num": 4, "ignore": [12, 13, 4, 11, 8, 3, 15, 0]}, "KS023": {"seeds": ["304", "313", "306", "309", "300", "314", "302", "310", "303", "315", "307", "308", "301", "311", "305", "312"], "fit_nums": ["698", "845", "319", "734", "908", "507", "45", "499", "175", "108", "419", "443", "116", "779", "159", "231"], "chain_num": 4, "ignore": [8, 10, 1, 13, 4, 15, 14, 7]}, "GLM_Sim_05": {"seeds": ["301", "315", "300", "302", "305", "304", "313", "314", "311", "309", "306", "307", "308", "310", "303", "312"], "fit_nums": ["425", "231", "701", "375", "343", "902", "623", "125", "921", "637", "393", "964", "678", "930", "796", "42"], "chain_num": 2, "ignore": [11, 2, 5, 1, 4, 9, 15, 12]}, "CSHL061": {"seeds": ["305", "315", "304", "303", "309", "310", "302", "300", "314", "306", "311", "313", "301", "308", "307", "312"], "fit_nums": ["396", "397", "594", "911", "308", "453", "686", "552", "103", "209", "128", "892", "345", "925", "777", "396"], "chain_num": 4, "ignore": [11, 13, 7, 15, 14, 3, 0, 4]}, "CSHL051": {"seeds": ["303", "310", "306", "302", "309", "305", "313", "308", "300", "314", "311", "307", "312", "304", "315", "301"], "fit_nums": ["69", "186", "49", "435", "103", "910", "705", "367", "303", "474", "596", "334", "929", "796", "616", "790"], "chain_num": 4, "ignore": [15, 12, 8, 13, 0, 2, 4, 5]}, "GLM_Sim_14": {"seeds": ["310", "311", "309", "313", "314", "300", "302", "304", "305", "306", "307", "312", "303", "301", "315", "308"], "fit_nums": ["616", "872", "419", "106", "940", "986", "599", "704", "218", "808", "244", "825", "448", "397", "552", "316"], "chain_num": 2, "ignore": [7, 11, 2, 15, 0, 13, 5, 10]}, "GLM_Sim_11_trick": {"seeds": ["411", "400", "408", "409", "415", "413", "410", "412", "406", "414", "403", "404", "401", "405", "407", "402"], "fit_nums": ["95", "508", "886", "384", "822", "969", "525", "382", "489", "436", "344", "537", "251", "223", "458", "401"], "chain_num": 2, "ignore": [10, 12, 4, 1, 0, 3, 2, 6]}, "GLM_Sim_16": {"seeds": ["302", "311", "303", "307", "313", "308", "309", "300", "305", "315", "304", "310", "312", "301", "314", "306"], "fit_nums": ["914", "377", "173", "583", "870", "456", "611", "697", "13", "713", "159", "248", "617", "37", "770", "780"], "chain_num": 2, "ignore": [4, 10, 5, 0, 13, 8, 6, 7]}, "ZM_3003": {"seeds": ["300", "304", "307", "312", "305", "310", "311", "314", "303", "308", "313", "301", "315", "309", "306", "302"], "fit_nums": ["603", "620", "657", "735", "357", "390", "119", "33", "62", "617", "209", "810", "688", "21", "744", "426"], "chain_num": 4, "ignore": [14, 7, 12, 1, 3, 4, 11, 8]}, "CSH_ZAD_022": {"seeds": ["305", "310", "311", "315", "303", "312", "314", "313", "307", "302", "300", "304", "301", "308", "306", "309"], "fit_nums": ["143", "946", "596", "203", "576", "403", "900", "65", "478", "325", "282", "513", "460", "42", "161", "970"], "chain_num": 4, "ignore": [9, 12, 4, 8, 3, 7, 0, 1]}, "GLM_Sim_07": {"seeds": ["300", "309", "302", "304", "305", "312", "301", "311", "315", "314", "308", "307", "303", "310", "306", "313"], "fit_nums": ["724", "701", "118", "230", "648", "426", "689", "114", "832", "731", "592", "519", "559", "938", "672", "144"], "chain_num": 1}, "KS017": {"seeds": ["311", "310", "306", "309", "303", "302", "308", "300", "313", "301", "314", "307", "315", "304", "312", "305"], "fit_nums": ["97", "281", "808", "443", "352", "890", "703", "468", "780", "708", "674", "27", "345", "23", "939", "457"], "chain_num": 4, "ignore": [0, 13, 8, 1, 12, 5, 10, 9]}, "GLM_Sim_11_sub": {"seeds": ["410", "414", "413", "404", "409", "415", "406", "408", "402", "411", "400", "405", "403", "407", "412", "401"], "fit_nums": ["830", "577", "701", "468", "929", "374", "954", "749", "937", "488", "873", "416", "612", "792", "461", "488"], "chain_num": 2}}
\ No newline at end of file
diff --git a/dyn_glm_chain_analysis.py b/dyn_glm_chain_analysis.py
index 04907a5a..c44c7897 100644
--- a/dyn_glm_chain_analysis.py
+++ b/dyn_glm_chain_analysis.py
@@ -11,7 +11,7 @@ import pyhsmm
 import pickle
 import seaborn as sns
 import sys
-from scipy.stats import zscore
+from scipy.stats import zscore, norm
 from scipy.optimize import minimize
 from itertools import combinations, product
 import matplotlib.gridspec as gridspec
@@ -20,6 +20,8 @@ import json
 import time
 import multiprocessing as mp
 from mcmc_chain_analysis import state_size_helper, state_num_helper, gamma_func, alpha_func, ll_func, r_hat_array_comp, rank_inv_normal_transform, eval_r_hat, eval_simple_r_hat
+import pandas as pd
+from pmf_analysis import pmf_type, type2color
 
 
 colors = np.genfromtxt('colors.csv', delimiter=',')
@@ -38,12 +40,20 @@ bias_cont_ticks = (np.arange(9), [-1, -.25, -.125, -.062, 0, .062, .125, .25, 1]
 contrasts_L = np.array([1., 0.987, 0.848, 0.555, 0.302, 0, 0, 0, 0, 0, 0])
 contrasts_R = np.array([1., 0.987, 0.848, 0.555, 0.302, 0, 0, 0, 0, 0, 0])[::-1]
 
-type_colours = ['g', 'b', 'r']
-
 
 def weights_to_pmf(weights, with_bias=1):
-    psi = weights[0] * contrasts_L + weights[1] * contrasts_R + with_bias * weights[-1]
-    return 1 / (1 + np.exp(-psi))
+    psi = weights[0] * contrasts_R + weights[1] * contrasts_L + with_bias * weights[-1]
+    return 1 / (1 + np.exp(psi))  # we somehow got the answers twisted, so we drop the minus here to get the opposite response probability for plotting
+
+performance_points = np.array([-1, -1, -1, -1, -1, 0, 0, 0, 0, 0, 0])
+def pmf_to_perf(pmf, def_points):
+    # determine performance of a pmf
+    # we use this to determine regressions in the behaviour of animals
+    # therefore, we exclude 0 as performance on it is 0.5 regardless of PMF, but it might
+    # overall lower performance. The removal of 0.5 later might also be a problem, let's see
+    relevant_points = def_points
+    relevant_points[5] = False
+    return np.mean(np.abs(performance_points[relevant_points] + pmf[relevant_points]))
 
 
 class MCMC_result_list:
@@ -477,9 +487,11 @@ def return_ascending_shuffled():
     return temp
 
 
-def contrasts_plot(test, state_sets, subject, save=False, show=False, dpi='figure', save_append='', consistencies=None):
+def contrasts_plot(test, state_sets, subject, save=False, show=False, dpi='figure', save_append='', consistencies=None, CMF=False):
     n = test.results[0].n_sessions
     trial_counter = 0
+    cnas = [] # contrasts aNd actions
+
     for seq_num in range(n):
         # if seq_num + 1 != 12:
         #     trial_counter += len(test.results[0].models[0].stateseqs[seq_num])
@@ -518,9 +530,13 @@ def contrasts_plot(test, state_sets, subject, save=False, show=False, dpi='figur
 
             plt.plot(active_trials, color=cmap(0.2 + 0.8 * seq_num / test.results[0].n_sessions), lw=4, label=label, alpha=0.7)
 
+        trial_counter += len(test.results[0].models[0].stateseqs[seq_num])
+
         ms = 6
         noise = np.zeros(len(c_n_a))# np.random.rand(len(c_n_a)) * 0.4 - 0.2
 
+        cnas.append(c_n_a)
+
         mask = c_n_a[:, -1] == 0
         plt.plot(np.where(mask)[0], 0.5 + 0.25 * (noise[mask] - c_n_a[mask, 0] + c_n_a[mask, 1]), 'o', c='b', ms=ms, label='Leftward', alpha=0.6)
 
@@ -557,8 +573,51 @@ def contrasts_plot(test, state_sets, subject, save=False, show=False, dpi='figur
             plt.show()
         else:
             plt.close()
-        trial_counter += len(test.results[0].models[0].stateseqs[seq_num])
 
+        if CMF and not test.results[0].name.startswith('Sim_'):
+            rt_data = pickle.load(open("./session_data/{} rt info {}".format(subject, seq_num + 1), 'rb'))
+            rt_data = rt_data[1:]
+            assert c_n_a.shape[0] == rt_data.shape[0]
+            df = pd.DataFrame(rt_data, columns = ['Signed contrast', 'RT', 'Responses'])
+            means = df.groupby('Signed contrast').mean()['RT']
+            stds = df.groupby('Signed contrast').sem()['RT']
+            plt.errorbar(means.index, means.values, stds.values)
+            plt.title(seq_num + 1, size=22)
+            plt.show()
+
+    return cnas
+
+def pmf_regressions(states_by_session, pmfs, durs):
+    # find out whether pmfs regressed
+    state_perfs = {}
+    state_counter = {}
+    current_best_state = -1
+    counter = 0
+    types = [0, 0, 0]
+
+    for sess in range(states_by_session.shape[1]):
+        states = np.where(states_by_session[:, sess])[0]
+
+        for s in states:
+            if s not in state_counter:
+                state_counter[s] = -1
+            state_counter[s] += 1
+            state_perfs[s] = pmf_to_perf(pmfs[s][1][state_counter[s]], pmfs[s][0])
+            if current_best_state == -1 or state_perfs[current_best_state] < state_perfs[s]:
+                current_best_state = s
+
+        if state_perfs[np.argmax(states_by_session[:, sess])] + 0.05 < state_perfs[current_best_state] and state_counter[np.argmax(states_by_session[:, sess])] > 1:
+            counter += 1
+            if sess < durs[0]:
+                types[0] += 1
+            elif sess < durs[0] + durs[1]:
+                types[1] += 1
+            else:
+                types[2] += 1
+            a, b = state_perfs[np.argmax(states_by_session[:, sess])], state_perfs[current_best_state]
+            print("Regression in session {} during {:.2f}% of session ({:.2f} instead of {:.2f})".format(sess + 1, np.max(states_by_session[:, sess]) * 100, a, b))
+
+    return [counter, states_by_session.shape[1], *types]
 
 def control_flow(test, indices, trials, func_init, first_for, second_for, end_first_for):
     # Generalised control flow for iterating over samples of mode across individually across sessions
@@ -602,6 +661,25 @@ def state_pmfs(test, trials, indices):
     return results['session_js'], results['pmfs']
 
 
+def state_weights(test, trials, indices):
+    def func_init(): return {'weightss': [], 'session_js': []}
+
+    def first_for(test, results):
+        results['weights'] = np.zeros(test.results[0].models[0].obs_distns[0].weights.shape[1])
+
+    def second_for(m, j, session_trials, trial_counter, results):
+        states, counts = np.unique(m.stateseqs[j][session_trials - trial_counter], return_counts=True)
+        for sub_state, c in zip(states, counts):
+            results['weights'] += m.obs_distns[sub_state].weights[j] * c / session_trials.shape[0]
+
+    def end_first_for(results, indices, j, **kwargs):
+        results['weightss'].append(results['weights'] / len(indices))
+        results['session_js'].append(j)
+
+    results = control_flow(test, indices, trials, func_init, first_for, second_for, end_first_for)
+    return results['session_js'], results['weightss']
+
+
 def lapse_sides(test, state_sets, indices):
     """Compute and plot a lapse differential across sessions.
 
@@ -875,7 +953,8 @@ def state_development(test, state_sets, indices, save=True, save_append='', show
                     trial_counter += len(state_seq)
                 counter += 1
         pmfs_to_score.append(np.mean(pmfs))
-    state_mapping = dict(zip(range(len(state_sets)), np.argsort(np.argsort(pmfs_to_score))))  # double argsort for ranks
+    # test.state_mapping = dict(zip(range(len(state_sets)), np.argsort(np.argsort(pmfs_to_score))))  # double argsort for ranks
+    test.state_mapping = dict(zip(range(len(state_sets)), np.flip(np.argsort((states_by_session != 0).argmax(axis=1)))))
 
     for state, trials in enumerate(state_sets):
         cmap = matplotlib.cm.get_cmap(cmaps[state]) if state < len(cmaps) else matplotlib.cm.get_cmap('Greys')
@@ -914,8 +993,8 @@ def state_development(test, state_sets, indices, save=True, save_append='', show
             temp = np.sum(pmfs[:, defined_points]) / (np.sum(defined_points))
             state_color = colors[int(temp * 101 - 1)]
 
-            ax1.fill_between(range(1, 1 + test.results[0].n_sessions), state_mapping[state] - 0.5,
-                             state_mapping[state] + states_by_session[state] - 0.5, color=state_color)
+            ax1.fill_between(range(1, 1 + test.results[0].n_sessions), test.state_mapping[state] - 0.5,
+                             test.state_mapping[state] + states_by_session[state] - 0.5, color=state_color)
 
         else:
             n_points = 150
@@ -924,11 +1003,11 @@ def state_development(test, state_sets, indices, save=True, save_append='', show
 
             for k in range(n_points-1):
                 ax1.fill_between([points[k], points[k+1]],
-                                 state_mapping[state] - 0.5, [state_mapping[state] + interpolation[k] - 0.5, state_mapping[state] + interpolation[k+1] - 0.5], color=cmap(0.3 + 0.7 * k / n_points))
-        ax1.annotate(state_mapping[state] + 1, (test.results[0].n_sessions + 0.1, state_mapping[state] - 0.15), fontsize=22, annotation_clip=False)
+                                 test.state_mapping[state] - 0.5, [test.state_mapping[state] + interpolation[k] - 0.5, test.state_mapping[state] + interpolation[k+1] - 0.5], color=cmap(0.3 + 0.7 * k / n_points))
+        ax1.annotate(test.state_mapping[state] + 1, (test.results[0].n_sessions + 0.1, test.state_mapping[state] - 0.15), fontsize=22, annotation_clip=False)
 
         if test.results[0].name.startswith('GLM_Sim_'):
-            ax1.plot(range(1, 1 + test.results[0].n_sessions), truth['state_map'][state_mapping[state]] + truth['state_posterior'][:, state] - 0.5, color='r')
+            ax1.plot(range(1, 1 + test.results[0].n_sessions), truth['state_map'][test.state_mapping[state]] + truth['state_posterior'][:, state] - 0.5, color='r')
 
         alpha_level = 0.3
         ax2.axvline(0.5, c='grey', alpha=alpha_level, zorder=4)
@@ -943,23 +1022,23 @@ def state_development(test, state_sets, indices, save=True, save_append='', show
         #         defined_points[[0, 1, -2, -1]] = True
         if separate_pmf:
             for j, pmf in zip(session_js, pmfs):
-                ax2.plot(np.where(defined_points)[0] / (len(defined_points)-1), pmf[defined_points] - 0.5 + state_mapping[state], color=cmap(0.2 + 0.8 * j / test.results[0].n_sessions))
-                ax2.plot(np.where(defined_points)[0] / (len(defined_points)-1), pmf[defined_points] - 0.5 + state_mapping[state], ls='', ms=7, marker='*', color=cmap(j / test.results[0].n_sessions))
+                ax2.plot(np.where(defined_points)[0] / (len(defined_points)-1), pmf[defined_points] - 0.5 + test.state_mapping[state], color=cmap(0.2 + 0.8 * j / test.results[0].n_sessions))
+                ax2.plot(np.where(defined_points)[0] / (len(defined_points)-1), pmf[defined_points] - 0.5 + test.state_mapping[state], ls='', ms=7, marker='*', color=cmap(j / test.results[0].n_sessions))
             all_pmfs.append((defined_points, pmfs))
         else:
             temp = np.percentile(pmfs, [2.5, 97.5], axis=0)
-            ax2.plot(np.where(defined_points)[0] / (len(defined_points)-1), pmfs[:, defined_points].mean(axis=0) - 0.5 + state_mapping[state], color=state_color)
-            ax2.plot(np.where(defined_points)[0] / (len(defined_points)-1), pmfs[:, defined_points].mean(axis=0) - 0.5 + state_mapping[state], ls='', ms=7, marker='*', color=state_color)
-            ax2.fill_between(np.where(defined_points)[0] / (len(defined_points)-1), temp[1, defined_points] - 0.5 + state_mapping[state], temp[0, defined_points] - 0.5 + state_mapping[state], alpha=0.2, color=state_color)
+            ax2.plot(np.where(defined_points)[0] / (len(defined_points)-1), pmfs[:, defined_points].mean(axis=0) - 0.5 + test.state_mapping[state], color=state_color)
+            ax2.plot(np.where(defined_points)[0] / (len(defined_points)-1), pmfs[:, defined_points].mean(axis=0) - 0.5 + test.state_mapping[state], ls='', ms=7, marker='*', color=state_color)
+            ax2.fill_between(np.where(defined_points)[0] / (len(defined_points)-1), temp[1, defined_points] - 0.5 + test.state_mapping[state], temp[0, defined_points] - 0.5 + test.state_mapping[state], alpha=0.2, color=state_color)
             all_pmfs.append((defined_points, pmfs[:, defined_points].mean(axis=0)))
 
         if test.results[0].name.startswith('GLM_Sim_'):
             sim_pmf = weights_to_pmf(truth['weights'][state])
-            ax2.plot(np.where(defined_points)[0] / (len(defined_points)-1), sim_pmf[defined_points] - 0.5 + truth['state_map'][state_mapping[state]], color='r')
+            ax2.plot(np.where(defined_points)[0] / (len(defined_points)-1), sim_pmf[defined_points] - 0.5 + truth['state_map'][test.state_mapping[state]], color='r')
 
-        ax2.axhline(state_mapping[state] + 0.5, c='k')
-        ax2.axhline(state_mapping[state], c='grey', alpha=alpha_level, zorder=4)
-        ax1.axhline(state_mapping[state] + 0.5, c='grey', alpha=alpha_level, zorder=4)
+        ax2.axhline(test.state_mapping[state] + 0.5, c='k')
+        ax2.axhline(test.state_mapping[state], c='grey', alpha=alpha_level, zorder=4)
+        ax1.axhline(test.state_mapping[state] + 0.5, c='grey', alpha=alpha_level, zorder=4)
 
     if not test.results[0].name.startswith('Sim_'):
         perf = np.zeros(test.results[0].n_sessions)
@@ -980,6 +1059,14 @@ def state_development(test, state_sets, indices, save=True, save_append='', show
         ax0.fill_between(range(1, 1 + test.results[0].n_sessions), perf - 0.5, -0.5, color='k')
 
     durs, state_types = state_type_durs(states_by_session, all_pmfs)
+
+    # how many states per session per type
+    states_per_sess = np.sum(states_by_session > 0.05, axis=0)
+    if durs[0] > 0 and durs[1] > 0 and durs[2] > 1:
+        states_per_type = [np.mean(states_per_sess[:durs[0]]), np.mean(states_per_sess[durs[0]:durs[0]+durs[1]]), np.mean(states_per_sess[durs[0]+durs[1]:])]
+    else:
+        states_per_type = []
+    # other statistics
     dur_counter = 1
     contrast_intro_types = [0, 0, 0, 0]
     state, when = np.where(states_by_session > 0.05)
@@ -987,7 +1074,7 @@ def state_development(test, state_sets, indices, save=True, save_append='', show
     covered_states = []
     for i, d in enumerate(durs):
         if type_coloring:
-            ax0.fill_between(range(dur_counter, 1 + dur_counter + d), 0.5, -0.5, color=type_colours[i], zorder=0, alpha=0.3)
+            ax0.fill_between(range(dur_counter, 1 + dur_counter + d), 0.5, -0.5, color=type2color[i], zorder=0, alpha=0.3)
         dur_counter += d
 
         # find out during which state type which contrast was introduced
@@ -1041,13 +1128,86 @@ def state_development(test, state_sets, indices, save=True, save_append='', show
     plt.tight_layout()
     if save:
         print("saving with {} dpi".format(dpi))
-        plt.savefig("dynamic_GLM_figures/meta state development_{}_{}{}.png".format(test.results[0].name, separate_pmf, save_append), dpi=dpi)
+        plt.savefig("dynamic_GLM_figures/meta_state_development_{}_{}{}.png".format(test.results[0].name, separate_pmf, save_append), dpi=dpi)
     if show:
         plt.show()
     else:
         plt.close()
 
-    return states_by_session, all_pmfs, durs, state_types, contrast_intro_types, smart_divide(introductions_by_stage, np.array(durs)), introductions_by_stage
+    return states_by_session, all_pmfs, durs, state_types, contrast_intro_types, smart_divide(introductions_by_stage, np.array(durs)), introductions_by_stage, states_per_type
+
+def compare_pmfs(test, state_sets, indices, states2compare, states_by_session, all_pmfs, title=""):
+    """
+       Take a set of states, and plot out their PMFs on all sessions on which they occur.
+       See how different they really are.
+
+       Takes states_by_session and all_pmfs as input from state_development
+    """
+    colors = ['blue', 'orange', 'green', 'black', 'red']
+    assert len(states2compare) <= len(colors)
+    # subtract 1 to get internal numbering
+    states2compare = [s - 1 for s in states2compare]
+    # transform desired states into the actual numbering, before ordering by bias
+    states2compare = [key for key in test.state_mapping.keys() if test.state_mapping[key] in states2compare]
+
+    sessions = np.where(states_by_session[states2compare].sum(0))[0]
+
+    for i, state in enumerate(states2compare):
+        counter = 0
+        for j, session in enumerate(sessions):
+            plt.subplot(1, len(sessions), j + 1)
+            if i == 0:
+                plt.title(session)
+            if states_by_session[state, session] > 0:
+                plt.plot(np.where(all_pmfs[state][0])[0], (all_pmfs[state][1][counter])[all_pmfs[state][0]], c=colors[i])
+                counter += 1
+            plt.ylim(0, 1)
+            if j != 0:
+                plt.gca().set_yticks([])
+    # plt.tight_layout()
+    if title != "":
+        plt.savefig(title)
+    plt.show()
+
+
+def compare_weights(test, state_sets, indices, states2compare, states_by_session, title=""):
+    """
+       Take a set of states, and plot out their weights on all sessions on which they occur.
+       See how different they really are.
+       Similar to compare_pmfs
+
+       Takes states_by_session as input from state_development
+    """
+    colors = ['blue', 'orange', 'green', 'black', 'red']
+    assert len(states2compare) <= len(colors)
+    # subtract 1 to get internal numbering
+    states2compare = [s - 1 for s in states2compare]
+    # transform desired states into the actual numbering, before ordering by bias
+    states2compare = [key for key in test.state_mapping.keys() if test.state_mapping[key] in states2compare]
+
+    sessions = np.where(states_by_session[states2compare].sum(0))[0]
+
+    state_counter = -1
+    for state, trials in enumerate(state_sets):
+        if state not in states2compare:
+            continue
+        state_counter += 1
+        _, weights = state_weights(test, trials, indices)
+        counter = 0
+        for j, session in enumerate(sessions):
+            plt.subplot(1, len(sessions), j + 1)
+            if state == 0:
+                plt.title(session)
+            if states_by_session[state, session] > 0:
+                plt.plot(weights[counter], c=colors[state_counter])
+                counter += 1
+            plt.ylim(-6, 6)
+            if j != 0:
+                plt.gca().set_yticks([])
+    # plt.tight_layout()
+    if title != "":
+        plt.savefig(title)
+    plt.show()
 
 
 def smart_divide(a, b):
@@ -1097,11 +1257,13 @@ if __name__ == "__main__":
     fit_type = ['prebias', 'bias', 'all', 'prebias_plus', 'zoe_style'][0]
     if fit_type == 'bias':
         loading_info = json.load(open("canonical_infos_bias.json", 'r'))
+        r_hats = json.load(open("canonical_info_r_hats_bias.json", 'r'))
     elif fit_type == 'prebias':
         loading_info = json.load(open("canonical_infos.json", 'r'))
+        r_hats = json.load(open("canonical_info_r_hats.json", 'r'))
     subjects = list(loading_info.keys())
 
-    r_hats = []
+    r_hats = {}
 
     # R^hat tests
     # test = MCMC_result_list([fake_result(100) for i in range(8)])
@@ -1112,10 +1274,9 @@ if __name__ == "__main__":
     check_r_hats = False
     if check_r_hats:
         subjects = list(loading_info.keys())
-        subjects = ['KS014']
         for subject in subjects:
-            # if subject.startswith('GLM'):
-            #     continue
+            if subject.startswith('GLM_Sim_07') or subject.startswith('GLM_Sim_11'):
+                continue
             print("_________________________________")
             print(subject)
             fit_num = loading_info[subject]['fit_nums'][-1]
@@ -1132,7 +1293,7 @@ if __name__ == "__main__":
 
             # mins = find_good_chains(chains[:, :-1].reshape(32, chains.shape[1] // 2))
             sol, final_r_hat = find_good_chains_unsplit_greedy(chains1, chains2, chains3, chains4, reduce_to=chains1.shape[0] // 2)
-            r_hats.append((subject, final_r_hat))
+            r_hats[subject] = final_r_hat
             loading_info[subject]['ignore'] = sol
 
         print(r_hats)
@@ -1233,171 +1394,88 @@ def plot_pmf_types(pmf_types, subject, fit_type, save=True, show=False):
         plt.close()
 
 
-def pmf_type(pmf):
-    if pmf[-1] - pmf[0] < 0.2:
-        return 0
-    elif pmf[-1] - pmf[0] < 0.6:# and np.abs(pmf[0] + pmf[-1] - 1) > 0.1:
-        return 1
-    else:
-        return 2
-
-
-type2color = {0: 'green', 1: 'blue', 2: 'red'}
+def two_sample_binom_test(s1, s2):
+    p1, p2 = np.mean(s1), np.mean(s2)
+    n1, n2 = s1.size, s2.size
+    p = (n1 * p1 + n2 * p2) / (n1 + n2)
+    if p == 1. or p == 0.:
+        return 0., 0.99
+    z = (p1 - p2) / np.sqrt(p * (1 - p) * (1 / n1 + 1 / n2))
+    p_val = (1 - norm.cdf(np.abs(z))) * 2
+    return z, p_val
+
+
+def compare_performance(cnas, contrasts=(1, 1), title=""):
+    # compare the performance on a certain contrast
+    # cnas is set of data for all sessions (from contrasts_plot)
+    # contrasts is a tuple specifying the contrast of interest (first or second column, contrast strength)
+    plt.figure(figsize=(16*0.75, 9*0.75))
+    for i in range(len(cnas) - 1):
+        if cnas[i][cnas[i][:, contrasts[0]] == contrasts[1], -1].shape[0] == 0 or cnas[i+1][cnas[i+1][:, contrasts[0]] == contrasts[1], -1].shape[0] == 0:
+            continue
+        perf1, perf2 = np.mean(cnas[i][cnas[i][:, contrasts[0]] == contrasts[1], -1]), np.mean(cnas[i+1][cnas[i+1][:, contrasts[0]] == contrasts[1], -1])
+        p = two_sample_binom_test(cnas[i][cnas[i][:, contrasts[0]] == contrasts[1], -1], cnas[i+1][cnas[i+1][:, contrasts[0]] == contrasts[1], -1])[1]
+        factor = (perf1 > perf2) * 2 - 1
+        plt.annotate(xy=(i + 0.5, (perf1 + perf2) / 2 + factor * 0.025), text=("{:.3f}".format(p))[1:])
+        color = 'red' if p < 0.05 else 'blue'
+        plt.plot([i, i+1], [perf1, perf2], color=color)
+    sns.despine()
+    cont_string = "Contrast {} ".format('right' if contrasts[0] else 'left')
+    plt.title(cont_string + str(contrasts[1]))
+    plt.ylim(-0.07, 1.07)
+    plt.tight_layout()
+    if title != "":
+        plt.savefig(title)
+    plt.show()
 
-if False:
 
-    all_changing_pmfs = pickle.load(open("changing_pmfs.p", 'rb'))
-    plt.figure(figsize=(16, 9))
-    for i, pmf in enumerate(all_changing_pmfs):
-        plt.subplot(4, 7, i + 1)
-        for p in pmf[1]:
-            plt.plot(np.where(pmf[0])[0], p[pmf[0]], color=type2color[pmf_type(p)])
-        plt.ylim(0, 1)
+def type_hist(data):
+    highest = int(data.max())
+    if (data % 1 == 0).all():
+        bins = np.arange(highest + 2) - 0.5
+    else:
+        bins = np.histogram(data)[1]
+    hist_max = 0
+    for i in range(data.shape[1]):
+        hist_max = max(hist_max, np.histogram(data[:, i], bins)[0].max())
+
+    plt.subplot(3, 1, 1)
+    assert np.histogram(data[:, 0])[0].sum() == np.histogram(data[:, 0], bins)[0].sum()
+    plt.hist(data[:, 0], alpha=1/3, label="type 1", align='mid', bins=bins)
+    plt.xlim(-0.5, highest + 1)
+    plt.ylim(0, hist_max + 1)
+
+    plt.subplot(3, 1, 2)
+    assert np.histogram(data[:, 1])[0].sum() == np.histogram(data[:, 1], bins)[0].sum()
+    plt.hist(data[:, 1], alpha=1/3, label="type 2", align='mid', bins=bins)
+    plt.xlim(-0.5, highest + 1)
+    plt.ylim(0, hist_max + 1)
+
+    plt.subplot(3, 1, 3)
+    assert np.histogram(data[:, 2])[0].sum() == np.histogram(data[:, 2], bins)[0].sum()
+    plt.hist(data[:, 2], alpha=1/3, label="type 3", align='mid', bins=bins)
+    plt.xlim(-0.5, highest + 1)
+    plt.ylim(0, hist_max + 1)
+
+    # plt.legend()
+    plt.show()
 
-        sns.despine()
-        if i+1 != 22:
-            plt.gca().set_xticks([])
-            plt.gca().set_yticks([])
-        else:
-            plt.xlabel("Contrasts", size=22)
-            plt.ylabel("P(rightwards)", size=22)
-            plt.gca().set_xticks([0, 5, 10], [-1, 0, 1], size=16)
-            plt.gca().set_yticks([0, 0.25, 0.5, 0.75, 1], [0, 0.25, 0.5, 0.75, 1], size=16)
 
-        if i + 1 == 30:
-            break
+if 0:
+    all_intros = pickle.load(open("all_intros.p", 'rb'))
+    all_intros_div = pickle.load(open("all_intros_div.p", 'rb'))
+    all_states_per_type = pickle.load(open("all_states_per_type.p", 'rb'))
 
-    plt.tight_layout()
-    plt.savefig("changing pmfs")
-    plt.show()
+    # There are 5 mice with 0 type 2 intros, but only 3 mice with no type 2 stats.
+    # That is because they come up, but don't explain the necessary 50% to start phase 2.
+    type_hist(all_intros)
+    type_hist(all_intros_div)
+    type_hist(all_states_per_type)
     quit()
 
-    type_2_assyms = []
-    tick_size = 14
-    label_size = 26
-    all_first_pmfs = pickle.load(open("pmfs_temp.p", 'rb'))
-    plt.figure(figsize=(16, 9))
-    plt.subplot(1, 3, 1)
-    counter = [[0, 0], [0, 0]]
-    save_title = "all types" if False else "KS014 types"
-    if save_title == "KS014 types":
-        all_first_pmfs = {'KS014': all_first_pmfs['KS014']}
-
-    for key in all_first_pmfs:
-        x = all_first_pmfs[key]
-        if type(x[0]) == int:
-            continue
-        linestyle = '-' if x[2] == 0 else '--'
-        plt.plot(np.where(x[1])[0], x[0][x[1]], linestyle=linestyle, c='g')
-    plt.ylim(0, 1)
-    plt.gca().set_xticks(np.arange(11), all_conts, size=tick_size)
-    plt.gca().set_yticks([0, 0.25, 0.5, 0.75, 1], [0, 0.25, 0.5, 0.75, 1], size=tick_size)
-    plt.gca().spines[['right', 'top']].set_visible(False)
-    plt.xlim(0, 10)
-    plt.xticks(rotation=45)
-    plt.gca().set_ylabel("P(rightwards)", size=label_size)
-
-    plt.subplot(1, 3, 2)
-    for key in all_first_pmfs:
-        x = all_first_pmfs[key]
-        if type(x[3]) == int:
-            continue
-        type_2_assyms.append(np.abs(x[3][0] + x[3][-1] - 1))
-        linestyle = '-' if x[5] == 0 else '--'
-        counter[0][0 if x[5] == 0 else 1] += 1
-        if linestyle == '--':
-            continue
-        plt.plot(np.where(x[4])[0], x[3][x[4]], linestyle=linestyle, c='b')
-    plt.gca().set_yticks([])
-    plt.ylim(0, 1)
-    plt.gca().set_xticks(np.arange(11), all_conts, size=tick_size)
-    plt.gca().spines[['right', 'top']].set_visible(False)
-    plt.xticks(rotation=45)
-    plt.xlim(0, 10)
-    plt.gca().set_xlabel("Contrasts", size=label_size)
-
-    plt.subplot(1, 3, 3)
-    for key in all_first_pmfs:
-        x = all_first_pmfs[key]
-        if type(x[6]) == int:
-            continue
-        linestyle = '-' if x[8] == 0 else '--'
-        counter[1][0 if x[8] == 0 else 1] += 1
-        if linestyle == '--':
-            continue
-        plt.plot(np.where(x[7])[0], x[6][x[7]], linestyle=linestyle, c='r')
-    plt.gca().set_yticks([])
-    plt.ylim(0, 1)
-    plt.gca().set_xticks(np.arange(11), all_conts, size=tick_size)
-    plt.gca().spines[['right', 'top']].set_visible(False)
-    plt.xlim(0, 10)
-    plt.xticks(rotation=45)
-
-    print(counter)
-    plt.tight_layout()
-    plt.savefig(save_title)
-    plt.show()
-    if save_title == "KS014 types":
-        quit()
-
-    counter = 0
-    fig, ax = plt.subplots(1, 3, figsize=(16, 9))
-    for key in all_first_pmfs:
-        x = all_first_pmfs[key]
-        if type(x[3]) == int:
-            continue
-        linestyle = '-' if x[5] == 0 else '--'
-        if linestyle == '--':
-            continue
-        if np.abs(x[3][0] + x[3][-1] - 1) <= 0.1:
-            counter += 1
-            use_ax = 2
-        else:
-            use_ax = int(x[3][0] > 1 - x[3][-1])
-
-        ax[use_ax].plot(np.where(x[4])[0], x[3][x[4]], linestyle=linestyle, c='b')
-    ax[0].set_ylim(0, 1)
-    ax[0].set_xlim(0, 10)
-    ax[0].spines[['right', 'top']].set_visible(False)
-    ax[0].set_yticks([0, 0.25, 0.5, 0.75, 1], [0, 0.25, 0.5, 0.75, 1], size=tick_size)
-    ax[0].set_xticks(np.arange(11), all_conts, size=tick_size, rotation=45)
-    ax[0].set_ylabel("P(rightwards)", size=label_size)
-
-    ax[1].set_ylim(0, 1)
-    ax[1].set_xlim(0, 10)
-    ax[1].set_yticks([])
-    ax[1].spines[['right', 'top']].set_visible(False)
-    ax[1].set_xticks(np.arange(11), all_conts, size=tick_size, rotation=45)
-    ax[1].set_xlabel("Contrasts", size=label_size)
-
-    ax[2].set_ylim(0, 1)
-    ax[2].set_xlim(0, 10)
-    ax[2].set_yticks([])
-    ax[2].spines[['right', 'top']].set_visible(False)
-    ax[2].set_xticks(np.arange(11), all_conts, size=tick_size, rotation=45)
-    print(counter)
-    plt.tight_layout()
-    plt.savefig("differentiate type 2")
-    plt.show()
-    quit()
 
 if __name__ == "__main__":
 
-    # visualise pmf types
-    # lapses = [0.1, 0.2, 0.25, 0.33, 0.4, 0.45, 0.5, 0.55, 0.66, 0.9]
-    # test_pmf = np.zeros(4)
-    # for i, lapse_l in enumerate(lapses):
-    #     plt.subplot(1, 10, 1+i)
-    #     if i != 0:
-    #         plt.gca().set_yticklabels([])
-    #     plt.ylim(0, 1)
-    #     test_pmf[:2] = lapse_l
-    #     for lapse_r in np.linspace(0.02, 0.98, 33):
-    #         test_pmf[2:] = lapse_r
-    #         plt.plot([0, 1, 9, 10], test_pmf, c=type_colours[pmf_type(test_pmf)])
-    # plt.show()
-    # quit()
-
     fit_type = ['prebias', 'bias', 'all', 'prebias_plus', 'zoe_style'][0]
     if fit_type == 'bias':
         loading_info = json.load(open("canonical_infos_bias.json", 'r'))
@@ -1405,8 +1483,11 @@ if __name__ == "__main__":
     elif fit_type == 'prebias':
         loading_info = json.load(open("canonical_infos.json", 'r'))
         r_hats = json.load(open("canonical_info_r_hats.json", 'r'))
-    no_good_pcas = ['NYU-06', 'SWC_023']  # no good rhat: 'ibl_witten_13'
-    subjects = ['KS014']  # list(loading_info.keys())
+    no_good_pcas = ['NYU-06', 'SWC_023']
+    subjects = list(loading_info.keys())
+    # subjects = ['ZM_1897']
+
+    # meh pmfs: KS021
     print(subjects)
     fit_variance = [0.03, 0.002, 0.0005, 'uniform', 0, 0.008][0]
     dur = 'yes'
@@ -1433,11 +1514,16 @@ if __name__ == "__main__":
 
     abs_state_durs = []
     all_first_pmfs = {}
+    all_first_pmfs_typeless = {}
     all_pmf_diffs = []
     all_pmf_asymms = []
     all_pmfs = []
     all_changing_pmfs = []
+    all_changing_pmf_names = []
     all_intros = []
+    all_intros_div = []
+    all_states_per_type = []
+    regressions = []
 
     for subject in subjects:
 
@@ -1445,48 +1531,66 @@ if __name__ == "__main__":
             continue
 
         print(subject)
-        results = []
 
         try:
+            continue
             test = pickle.load(open("multi_chain_saves/canonical_result_{}_{}.p".format(subject, fit_type), 'rb'))
             print('loaded canoncial result')
 
+            mode_specifier = ''
             mode_indices = pickle.load(open("multi_chain_saves/mode_indices_{}_{}.p".format(subject, fit_type), 'rb'))
-            quit()
             state_sets = pickle.load(open("multi_chain_saves/state_sets_{}_{}.p".format(subject, fit_type), 'rb'))
+
             # lapse differential
             # lapse_sides(test, [s for s in state_sets if len(s) > 40], mode_indices)
 
             # training overview
-            states, pmfs, durs, _, contrast_intro_type, intros_by_type, undiv_intros = state_development(test, [s for s in state_sets if len(s) > 40], mode_indices, show=0, separate_pmf=1, type_coloring=True)
+            states, pmfs, durs, _, contrast_intro_type, intros_by_type, undiv_intros, states_per_type = state_development(test, [s for s in state_sets if len(s) > 40], mode_indices, show=0, separate_pmf=1, type_coloring=True)
+            regression = pmf_regressions(states, pmfs, durs)
+            regressions.append(regression)
+            continue
+            # compare_pmfs(test, [s for s in state_sets if len(s) > 40], mode_indices, [4, 5], states, pmfs, title="{} convergence pmf".format(subject))
+            # compare_weights(test, [s for s in state_sets if len(s) > 40], mode_indices, [4, 5], states, title="{} convergence weights".format(subject))
+            # quit()
+            all_first_pmfs_typeless[subject] = []
+            for pmf in pmfs:
+                all_first_pmfs_typeless[subject].append((pmf[0], pmf[1][0]))
             all_intros.append(undiv_intros)
+            all_intros_div.append(intros_by_type)
+            if states_per_type != []:
+                all_states_per_type.append(states_per_type)
+
             intros_by_type_sum += intros_by_type
             first_pmfs, changing_pmfs = get_first_pmfs(states, pmfs)
             for pmf in changing_pmfs:
                 if type(pmf[0]) == int:
                     continue
+                all_changing_pmf_names.append(subject)
                 all_changing_pmfs.append(pmf)
             all_first_pmfs[subject] = first_pmfs
             for pmf in pmfs:
+                all_pmfs.append(pmf)
                 for p in pmf[1]:
                     all_pmf_diffs.append(p[-1] - p[0])
                     all_pmf_asymms.append(np.abs(p[0] + p[-1] - 1))
-                    all_pmfs.append(p)
             contrast_intro_types.append(contrast_intro_type)
+            continue
             # state_development_single_sample(test, [mode_indices[0]], show=True, separate_pmf=True, save=False)
 
             # session overview
-            consistencies = pickle.load(open("multi_chain_saves/consistencies_{}_{}.p".format(subject, fit_type), 'rb'))
-            consistencies /= consistencies[0, 0]
-            contrasts_plot(test, [s for s in state_sets if len(s) > 40], dpi=300, subject=subject, save=True, show=True, consistencies=consistencies)
+            # consistencies = pickle.load(open("multi_chain_saves/consistencies_{}_{}.p".format(subject, fit_type), 'rb'))
+            # consistencies /= consistencies[0, 0]
+            # c_n_a, rt_data = contrasts_plot(test, [s for s in state_sets if len(s) > 40], dpi=300, subject=subject, save=True, show=False, consistencies=consistencies, CMF=True)
+            # compare_performance(cnas, (1, 1), title="{} contrast {} performance test".format(subject, (1, 1)))
+            # compare_performance(cnas, (0, 0.848))
 
             # duration of different state types (and also percentage of type activities)
-            abs_state_durs.append(durs)
-            simplex_durs = np.array(durs).reshape(1, 3)
-            print(simplex_durs / np.sum(simplex_durs))
-            from simplex_plot import projectSimplex
-            print(projectSimplex(simplex_durs / simplex_durs.sum(1)[:, None]))
-            continue
+            # abs_state_durs.append(durs)
+            # continue
+            # simplex_durs = np.array(durs).reshape(1, 3)
+            # print(simplex_durs / np.sum(simplex_durs))
+            # from simplex_plot import projectSimplex
+            # print(projectSimplex(simplex_durs / simplex_durs.sum(1)[:, None]))
 
             # compute state type proportions and split the pmfs accordingly
             # ret, trans, pmf_types = state_cluster_interpolation(states, pmfs)
@@ -1552,47 +1656,46 @@ if __name__ == "__main__":
             # plt.savefig("temp")
             # plt.close()
 
-            quit()
-
-            dim = 3
-            ev, eig, projection_matrix, dimreduc = test.state_pca(subject, pca_vecs='dists', dim=dim)
-            xy = np.vstack([dimreduc[i] for i in range(dim)])
-            from scipy.stats import gaussian_kde
-            z = gaussian_kde(xy)(xy)
-            pickle.dump((xy, z), open("multi_chain_saves/xyz_{}_{}.p".format(subject, fit_type), 'wb'))
-
-            if 'mode prob level' not in loading_info[subject]:
-                print(subject)
-                xy, z = pickle.load(open("multi_chain_saves/xyz_{}_{}.p".format(subject, fit_type), 'rb'))
-                happy = False
-                while not happy:
-                    print()
-                    print("Pick level")
-                    prob_level = float(input())
-                    print("Level is {}".format(prob_level))
-                    print("# of samples: {}".format((z > prob_level).sum()))
-                    mode_indices = np.where(z > prob_level)[0]
-                    if (z > prob_level).sum() > 0:
-                        print(xy[0][mode_indices].min(), xy[0][mode_indices].max(), xy[1][mode_indices].min(), xy[1][mode_indices].max())
-                        print("Happy?")
-                        happy = 'yes' == input()
-                print("Subset by factor?")
-                if input() == 'yes':
-                    print("Factor?")
-                    print(mode_indices.shape)
-                    factor = int(input())
-                    mode_indices = mode_indices[::factor]
-                    print(mode_indices.shape)
-                loading_info[subject]['mode prob level'] = prob_level
-
-                pickle.dump(mode_indices, open("multi_chain_saves/mode_indices_{}_{}.p".format(subject, fit_type), 'wb'))
-                consistencies = test.consistency_rsa(indices=mode_indices)
-                pickle.dump(consistencies, open("multi_chain_saves/consistencies_{}_{}.p".format(subject, fit_type), 'wb'))
-
-            string_prefix = ''
+            # dim = 3
+            # ev, eig, projection_matrix, dimreduc = test.state_pca(subject, pca_vecs='dists', dim=dim)
+            # xy = np.vstack([dimreduc[i] for i in range(dim)])
+            # from scipy.stats import gaussian_kde
+            # z = gaussian_kde(xy)(xy)
+            # pickle.dump((xy, z), open("multi_chain_saves/xyz_{}_{}.p".format(subject, fit_type), 'wb'))
+            # #
+            # string_prefix = 'second_'
+            #
+            # print(subject)
+            # xy, z = pickle.load(open("multi_chain_saves/xyz_{}_{}.p".format(subject, fit_type), 'rb'))
+            # quit()
+            # happy = False
+            # while not happy:
+            #     print()
+            #     print("Pick level")
+            #     prob_level = float(input())
+            #     print("Level is {}".format(prob_level))
+            #     print("# of samples: {}".format((z > prob_level).sum()))
+            #     mode_indices = np.where(z > prob_level)[0]
+            #     if (z > prob_level).sum() > 0:
+            #         print(xy[0][mode_indices].min(), xy[0][mode_indices].max(), xy[1][mode_indices].min(), xy[1][mode_indices].max())
+            #         print("Happy?")
+            #         happy = 'yes' == input()
+            # print("Subset by factor?")
+            # if input() == 'yes':
+            #     print("Factor?")
+            #     print(mode_indices.shape)
+            #     factor = int(input())
+            #     mode_indices = mode_indices[::factor]
+            #     print(mode_indices.shape)
+            # loading_info[subject]['mode prob level'] = prob_level
+            #
+            # pickle.dump(mode_indices, open("multi_chain_saves/{}mode_indices_{}_{}.p".format(string_prefix, subject, fit_type), 'wb'))
+            # consistencies = test.consistency_rsa(indices=mode_indices)
+            # pickle.dump(consistencies, open("multi_chain_saves/{}consistencies_{}_{}.p".format(string_prefix, subject, fit_type), 'wb'))
 
             mode_indices = pickle.load(open("multi_chain_saves/{}mode_indices_{}_{}.p".format(string_prefix, subject, fit_type), 'rb'))
-            consistencies = pickle.load(open("multi_chain_saves/{}consistencies_{}_{}.p".format(string_prefix, subject, fit_type), 'rb'))
+            consistencies = pickle.load(open("multi_chain_saves/{}mode_consistencies_{}_{}.p".format(string_prefix, subject, fit_type), 'rb'))
+            session_bounds = list(np.cumsum([len(s) for s in test.results[0].models[-1].stateseqs]))
 
             import scipy.cluster.hierarchy as hc
             consistencies /= consistencies[0, 0]
@@ -1620,10 +1723,9 @@ if __name__ == "__main__":
                     for x, y in zip(b, c):
                         state_sets.append(np.where(a == x)[0])
                     print("dumping state set")
-                    pickle.dump(state_sets, open("multi_chain_saves/state_sets_{}_{}.p".format(subject, fit_type), 'wb'))
-                    quit()
+                    pickle.dump(state_sets, open("multi_chain_saves/{}state_sets_{}_{}.p".format(string_prefix, subject, fit_type), 'wb'))
                     states, pmfs = state_development(test, [s for s in state_sets if len(s) > 40], mode_indices, save_append='_{}{}'.format(string_prefix, criterion), show=True, separate_pmf=True)
-                    contrasts_plot(test, [s for s in state_sets if len(s) > 40], subject=subject, save_append='_{}{}'.format(string_prefix, criterion), save=True, show=True)
+                    # contrasts_plot(test, [s for s in state_sets if len(s) > 40], subject=subject, save_append='_{}{}'.format(string_prefix, criterion), save=True, show=True)
 
                     # I think this is about finding out where states start and how long they last
                     # state_id, session_appears = np.where(states)
@@ -1660,18 +1762,14 @@ if __name__ == "__main__":
 
             plt.tight_layout()
             plt.savefig("peter figures/{}clustered_trials_{}_{}".format(string_prefix, subject, 'criteria comp').replace('.', '_'))
-            plt.show()
+            plt.close()
 
         except FileNotFoundError as e:
-            continue
             print(e)
-            r_hat = 1.5
-            for r in r_hats:
-                if r[0] == subject:
-                    r_hat = r[1]
+            continue
             print('no canoncial result')
-            print(r_hat)
-            if r_hat >= 1.05:
+            print(r_hats[subject])
+            if r_hats[subject] >= 1.05:
                 print("Skipping")
                 continue
             else:
@@ -1729,6 +1827,13 @@ if __name__ == "__main__":
 
     # pickle.dump(all_first_pmfs, open("pmfs_temp.p", 'wb'))
     # pickle.dump(all_changing_pmfs, open("changing_pmfs.p", 'wb'))
+    # pickle.dump(all_changing_pmf_names, open("changing_pmf_names.p", 'wb'))
+    # pickle.dump(all_first_pmfs_typeless, open("all_first_pmfs_typeless.p", 'wb'))
+    # pickle.dump(all_pmfs, open("all_pmfs.p", 'wb'))
+    # pickle.dump(all_intros, open("all_intros.p", 'wb'))
+    # pickle.dump(all_intros_div, open("all_intros_div.p", 'wb'))
+    # pickle.dump(all_states_per_type, open("all_states_per_type.p", 'wb'))
+    # pickle.dump(regressions, open("regressions.p", 'wb'))
     #
     # a = [x for x, y in zip(all_pmf_asymms, all_pmf_diffs) if y >= 0.2]
     # b = [y for x, y in zip(all_pmf_asymms, all_pmf_diffs) if y >= 0.2]
@@ -1758,6 +1863,16 @@ if __name__ == "__main__":
 
         plotSimplex(np.array(abs_state_durs), c='k', show=True)
 
+        plt.hist(abs_state_durs.sum(1), color='grey', bins=12)
+        sns.despine()
+        plt.xticks(size=26)
+        plt.yticks(size=26)
+        plt.ylabel("# of mice", size=40)
+        plt.xlabel('# of sessions', size=40)
+        plt.tight_layout()
+        plt.savefig("session_num_hist.png", dpi=300, transparent=True)
+        plt.show()
+
     if False:
         ax[0].set_ylim(0, 1)
         ax[1].set_ylim(0, 1)
diff --git a/index_mice.py b/index_mice.py
index 59af1292..ebb58b84 100644
--- a/index_mice.py
+++ b/index_mice.py
@@ -34,7 +34,7 @@ for filename in os.listdir("./dynamic_GLMiHMM_crossvals/"):
         local_dict = bias_subinfo
     if subject not in local_dict:
         local_dict[subject] = {"seeds": [], "fit_nums": [], "chain_num": 0}
-    if int(chain_num) == 0:
+    if int(chain_num) == 0:  # if this is the first file of that chain, save some info
         local_dict[subject]["seeds"].append(seed)
         local_dict[subject]["fit_nums"].append(fit_num)
     else:
@@ -62,7 +62,6 @@ for s in prebias_subinfo.keys():
         prebias_subinfo[s]["seeds"] = new_seeds
 
 
-
 big = []
 non_big = []
 sim_subjects = []
diff --git a/pmf_analysis.py b/pmf_analysis.py
new file mode 100644
index 00000000..9e7a0338
--- /dev/null
+++ b/pmf_analysis.py
@@ -0,0 +1,343 @@
+import pickle
+import matplotlib.pyplot as plt
+import numpy as np
+
+
+type2color = {0: 'green', 1: 'blue', 2: 'red'}
+all_conts = np.array([-1, -0.5, -.25, -.125, -.062, 0, .062, .125, .25, 0.5, 1])
+
+
+def pmf_type(pmf):
+    if pmf[-1] - pmf[0] < 0.2:
+        return 0
+    # elif pmf[-1] - pmf[0] < 0.4:# and np.abs(pmf[0] + pmf[-1] - 1) > 0.1:
+    elif max(pmf[0], 1 - pmf[-1]) > 0.37 or pmf[-1] - pmf[0] < 0.55:
+        if pmf[0] > 0.37 and 1 - pmf[-1] > 0.37:
+            print("___________Troublesome PMF: {}___________".format(pmf))
+            return 0
+        return 1
+    else:
+        return 2
+
+
+if __name__ == "__main__":
+    all_first_pmfs_typeless = pickle.load(open("all_first_pmfs_typeless", 'rb'))
+    all_pmfs = pickle.load(open("all_pmfs.p", 'rb'))
+
+    fewer_states_side = []
+    for key in all_first_pmfs_typeless:
+        animal_biases = np.zeros(2)
+        for defined_points, pmf in all_first_pmfs_typeless[key]:
+            bias = np.mean(pmf[defined_points])
+            if bias > 0.55:
+                animal_biases[0] += 1
+            elif bias < 0.45:
+                animal_biases[1] += 1
+        fewer_states_side.append(np.min(animal_biases / animal_biases.sum()))
+    plt.hist(fewer_states_side)
+    sns.despine()
+    plt.tight_layout()
+    plt.savefig("./meeting_figures/proportion_other_bias")
+    plt.show()
+
+    total_counter = 0
+    bias_counter = 0
+    tendency_counter = 0
+    lapse_counter = 0
+    for pmf in all_pmfs:
+        max_b, min_b = 0, 1
+        max_tendency, min_tendency = 0, 1
+        max_lapse_diff, min_lapse_diff = 0, 1
+        for p in pmf[1]:
+            if pmf[0][5]:  # if this part of the pmf is defined, just take it
+                bias = p[5]
+            deviation = 0
+            while True:  # just take the closest thing
+                deviation += 1
+                if pmf[0][5 - deviation] and pmf[0][5 + deviation]:
+                    bias = (p[5 - deviation] + p[5 + deviation]) / 2
+                    break
+            max_b = max(max_b, bias)
+            min_b = min(min_b, bias)
+            max_tendency = max(max_tendency, np.mean(p[pmf[0]]))
+            min_tendency = min(min_tendency, np.mean(p[pmf[0]]))
+            max_lapse_diff = max(max_lapse_diff, p[0] + p[-1] - 1)
+            min_lapse_diff = min(min_lapse_diff, p[0] + p[-1] - 1)
+        bias_changed = max_b > 0.55 and min_b < 0.45
+        tendency_changed = max_tendency > 0.55 and min_tendency < 0.45
+        lapse_changed = max_lapse_diff > 0.1 and min_lapse_diff < -0.1
+        bias_counter += bias_changed
+        tendency_counter += tendency_changed
+        lapse_counter += lapse_changed
+        if bias_changed or tendency_changed or lapse_changed:
+            total_counter += 1
+            for p in pmf[1]:
+                plt.plot(np.where(pmf[0])[0], p[pmf[0]])
+            plt.title("bias: {}, tendency: {}, lapse: {}".format(bias_changed, tendency_changed, lapse_changed))
+            plt.ylim(0, 1)
+            plt.axvline(5, color='grey')
+            plt.axhline(0.5, color='grey')
+            plt.savefig("./meeting_figures/bias_change_{}".format(total_counter))
+            plt.close()
+    print(total_counter)
+    print(bias_counter)
+    print(tendency_counter)
+    print(lapse_counter)
+
+    pmf_ranges = []
+    for key in all_first_pmfs_typeless:
+        for defined_points, pmf in all_first_pmfs_typeless[key]:
+            if pmf_type(pmf) == 2:
+                pmf_ranges.append(pmf[-1] - pmf[0])
+                # if pmf_ranges[-1] < 0.6:
+                #     plt.plot(pmf)
+                #     plt.title(pmf_ranges[-1])
+                #     plt.ylim(0, 1)
+                #     plt.show()
+    plt.hist(pmf_ranges, bins=40)
+    plt.title("Type 2 ranges")
+    plt.show()
+
+    lapses = []
+    for key in all_first_pmfs_typeless:
+        for defined_points, pmf in all_first_pmfs_typeless[key]:
+            if pmf_type(pmf) != 0:
+                lapses.append(max(pmf[0], 1 - pmf[-1]))
+    plt.hist(lapses, bins=40)
+    plt.title("Higher lapse rate of type != 0")
+    plt.show()
+
+    n_rows, n_cols = 5, 6
+    _, axs = plt.subplots(n_rows, n_cols, figsize=(16, 9))
+    for i, key in enumerate(all_first_pmfs_typeless):
+        if i == n_rows * 2:
+            break
+        axs[(i * 3) // n_cols, (i * 3) % n_cols].set_ylabel(key, size=17)
+        for defined_points, pmf in all_first_pmfs_typeless[key]:
+            axs[(i * 3 + pmf_type(pmf)) // n_cols, (i * 3 + pmf_type(pmf)) % n_cols].plot(np.where(defined_points)[0], pmf[defined_points], c=type2color[pmf_type(pmf)])
+    for i, ax in enumerate(axs):
+        for j, a in enumerate(ax):
+            a.spines[['right', 'top']].set_visible(False)
+            a.set_ylim(0, 1)
+            if i != n_rows - 1 or j != 0:
+                a.set_xticks([])
+                a.set_yticks([])
+            else:
+                tick_size = 12
+                a.set_xticks(np.arange(11), all_conts, size=tick_size, rotation=45)
+                a.set_yticks([0, 0.25, 0.5, 0.75, 1], [0, 0.25, 0.5, 0.75, 1], size=tick_size)
+    plt.tight_layout()
+    plt.savefig("animals 1")
+    plt.show()
+
+    n_rows, n_cols = 5, 6
+    _, axs = plt.subplots(n_rows, n_cols, figsize=(16, 9))
+    for i, key in enumerate(all_first_pmfs_typeless):
+        if i < n_rows * 2:
+            continue
+        elif i < n_rows * 4:
+            i = i - n_rows * 2
+        else:
+            break
+        if i == n_rows * 2:
+            break
+        axs[(i * 3) // n_cols, (i * 3) % n_cols].set_ylabel(key, size=17)
+        for defined_points, pmf in all_first_pmfs_typeless[key]:
+            axs[(i * 3 + pmf_type(pmf)) // n_cols, (i * 3 + pmf_type(pmf)) % n_cols].plot(np.where(defined_points)[0], pmf[defined_points], c=type2color[pmf_type(pmf)])
+    for i, ax in enumerate(axs):
+        for j, a in enumerate(ax):
+            a.spines[['right', 'top']].set_visible(False)
+            a.set_ylim(0, 1)
+            if i != n_rows - 1 or j != 0:
+                a.set_xticks([])
+                a.set_yticks([])
+            else:
+                tick_size = 12
+                a.set_xticks(np.arange(11), all_conts, size=tick_size, rotation=45)
+                a.set_yticks([0, 0.25, 0.5, 0.75, 1], [0, 0.25, 0.5, 0.75, 1], size=tick_size)
+    plt.tight_layout()
+    plt.savefig("animals 2")
+    plt.show()
+
+    n_rows, n_cols = 5, 6
+    _, axs = plt.subplots(n_rows, n_cols, figsize=(16, 9))
+    for i, key in enumerate(all_first_pmfs_typeless):
+        if i < n_rows * 4:
+            continue
+        elif i < n_rows * 6:
+            i = i - n_rows * 4
+        else:
+            break
+        if i == n_rows * 4:
+            break
+        axs[(i * 3) // n_cols, (i * 3) % n_cols].set_ylabel(key, size=17)
+        for defined_points, pmf in all_first_pmfs_typeless[key]:
+            axs[(i * 3 + pmf_type(pmf)) // n_cols, (i * 3 + pmf_type(pmf)) % n_cols].plot(np.where(defined_points)[0], pmf[defined_points], c=type2color[pmf_type(pmf)])
+    for i, ax in enumerate(axs):
+        for j, a in enumerate(ax):
+            a.spines[['right', 'top']].set_visible(False)
+            a.set_ylim(0, 1)
+            if i != n_rows - 1 or j != 0:
+                a.set_xticks([])
+                a.set_yticks([])
+            else:
+                tick_size = 12
+                a.set_xticks(np.arange(11), all_conts, size=tick_size, rotation=45)
+                a.set_yticks([0, 0.25, 0.5, 0.75, 1], [0, 0.25, 0.5, 0.75, 1], size=tick_size)
+    plt.tight_layout()
+    plt.savefig("animals 3")
+    plt.show()
+
+    n_rows, n_cols = 5, 6
+    _, axs = plt.subplots(n_rows, n_cols, figsize=(16, 9))
+    for i, key in enumerate(all_first_pmfs_typeless):
+        if i < n_rows * 6:
+            continue
+        elif i < n_rows * 8:
+            i = i - n_rows * 6
+        else:
+            break
+        if i == n_rows * 6:
+            break
+        axs[(i * 3) // n_cols, (i * 3) % n_cols].set_ylabel(key, size=17)
+        for defined_points, pmf in all_first_pmfs_typeless[key]:
+            axs[(i * 3 + pmf_type(pmf)) // n_cols, (i * 3 + pmf_type(pmf)) % n_cols].plot(np.where(defined_points)[0], pmf[defined_points], c=type2color[pmf_type(pmf)])
+    for i, ax in enumerate(axs):
+        for j, a in enumerate(ax):
+            a.spines[['right', 'top']].set_visible(False)
+            a.set_ylim(0, 1)
+            if i != n_rows - 1 or j != 0:
+                a.set_xticks([])
+                a.set_yticks([])
+            else:
+                tick_size = 12
+                a.set_xticks(np.arange(11), all_conts, size=tick_size, rotation=45)
+                a.set_yticks([0, 0.25, 0.5, 0.75, 1], [0, 0.25, 0.5, 0.75, 1], size=tick_size)
+    plt.tight_layout()
+    plt.savefig("animals 4")
+    plt.show()
+
+    # Collection of PMFs which change state type
+    all_changing_pmfs = pickle.load(open("changing_pmfs.p", 'rb'))
+    all_changing_pmf_names = pickle.load(open("changing_pmf_names.p", 'rb'))
+
+    plt.figure(figsize=(16, 9))
+    for i, (pmf, name) in enumerate(zip(all_changing_pmfs, all_changing_pmf_names)):
+        plt.subplot(4, 7, i + 1)
+        plt.title(name)
+        for p in pmf[1]:
+            plt.plot(np.where(pmf[0])[0], p[pmf[0]], color=type2color[pmf_type(p)])
+        plt.ylim(0, 1)
+
+        sns.despine()
+        if i+1 != 22:
+            plt.gca().set_xticks([])
+            plt.gca().set_yticks([])
+        else:
+            plt.xlabel("Contrasts", size=22)
+            plt.ylabel("P(rightwards)", size=22)
+            plt.gca().set_xticks([0, 5, 10], [-1, 0, 1], size=16)
+            plt.gca().set_yticks([0, 0.25, 0.5, 0.75, 1], [0, 0.25, 0.5, 0.75, 1], size=16)
+
+        if i + 1 == 30:
+            break
+
+    plt.tight_layout()
+    plt.savefig("changing pmfs")
+    plt.show()
+
+    # All first PMFs
+    tick_size = 14
+    label_size = 26
+    all_first_pmfs = pickle.load(open("pmfs_temp.p", 'rb'))
+    n_rows, n_cols = 1, 3
+    _, axs = plt.subplots(n_rows, n_cols, figsize=(16, 9))
+    counter = [[0, 0], [0, 0]]
+    save_title = "all types" if True else "KS014 types"
+    if save_title == "KS014 types":
+        all_first_pmfs_typeless = {'KS014': all_first_pmfs_typeless['KS014']}
+
+    for key in all_first_pmfs_typeless:
+        x = all_first_pmfs_typeless[key]
+        for pmf in x:
+            axs[pmf_type(pmf[1])].plot(np.where(pmf[0])[0], pmf[1][pmf[0]], c=type2color[pmf_type(pmf[1])])
+    axs[0].set_ylim(0, 1)
+    axs[1].set_ylim(0, 1)
+    axs[2].set_ylim(0, 1)
+    axs[0].set_xticks(np.arange(11), all_conts, size=tick_size, rotation=45)
+    axs[1].set_xticks(np.arange(11), all_conts, size=tick_size, rotation=45)
+    axs[2].set_xticks(np.arange(11), all_conts, size=tick_size, rotation=45)
+    axs[0].set_yticks([0, 0.25, 0.5, 0.75, 1], [0, 0.25, 0.5, 0.75, 1], size=tick_size)
+    axs[1].set_yticks([0, 0.25, 0.5, 0.75, 1], [0, 0.25, 0.5, 0.75, 1], size=tick_size)
+    axs[2].set_yticks([0, 0.25, 0.5, 0.75, 1], [0, 0.25, 0.5, 0.75, 1], size=tick_size)
+    axs[0].spines[['right', 'top']].set_visible(False)
+    axs[1].spines[['right', 'top']].set_visible(False)
+    axs[2].spines[['right', 'top']].set_visible(False)
+    axs[0].set_xlim(0, 10)
+    axs[1].set_xlim(0, 10)
+    axs[2].set_xlim(0, 10)
+    axs[0].set_ylabel("P(rightwards)", size=label_size)
+    axs[0].set_xlabel("Contrasts", size=label_size)
+
+    print(counter)
+    plt.tight_layout()
+    plt.savefig(save_title)
+    plt.show()
+    if save_title == "KS014 types":
+        quit()
+
+    # Type 2 PMFs
+    counter = 0
+    fig, ax = plt.subplots(1, 3, figsize=(16, 9))
+    for key in all_first_pmfs_typeless:
+        for defined_points, pmf in all_first_pmfs_typeless[key]:
+            if pmf_type(pmf) != 1:
+                continue
+            # linestyle = '-' if x[5] == 0 else '--'
+            # if linestyle == '--':
+            #     continue
+            if np.abs(pmf[0] + pmf[-1] - 1) <= 0.1:
+                counter += 1
+                use_ax = 2
+            else:
+                use_ax = int(pmf[0] > 1 - pmf[-1])
+
+            ax[use_ax].plot(np.where(defined_points)[0], pmf[defined_points], c='b')
+    ax[0].set_ylim(0, 1)
+    ax[0].set_xlim(0, 10)
+    ax[0].spines[['right', 'top']].set_visible(False)
+    ax[0].set_yticks([0, 0.25, 0.5, 0.75, 1], [0, 0.25, 0.5, 0.75, 1], size=tick_size)
+    ax[0].set_xticks(np.arange(11), all_conts, size=tick_size, rotation=45)
+    ax[0].set_ylabel("P(rightwards)", size=label_size)
+
+    ax[1].set_ylim(0, 1)
+    ax[1].set_xlim(0, 10)
+    ax[1].set_yticks([])
+    ax[1].spines[['right', 'top']].set_visible(False)
+    ax[1].set_xticks(np.arange(11), all_conts, size=tick_size, rotation=45)
+    ax[1].set_xlabel("Contrasts", size=label_size)
+
+    ax[2].set_ylim(0, 1)
+    ax[2].set_xlim(0, 10)
+    ax[2].set_yticks([])
+    ax[2].spines[['right', 'top']].set_visible(False)
+    ax[2].set_xticks(np.arange(11), all_conts, size=tick_size, rotation=45)
+    print(counter)
+    plt.tight_layout()
+    plt.savefig("differentiate type 2")
+    plt.show()
+
+
+    # visualise pmf types
+    # lapses = [0.1, 0.2, 0.25, 0.33, 0.4, 0.45, 0.5, 0.55, 0.66, 0.9]
+    # test_pmf = np.zeros(4)
+    # for i, lapse_l in enumerate(lapses):
+    #     plt.subplot(1, 10, 1+i)
+    #     if i != 0:
+    #         plt.gca().set_yticklabels([])
+    #     plt.ylim(0, 1)
+    #     test_pmf[:2] = lapse_l
+    #     for lapse_r in np.linspace(0.02, 0.98, 33):
+    #         test_pmf[2:] = lapse_r
+    #         plt.plot([0, 1, 9, 10], test_pmf, c=type_colours[pmf_type(test_pmf)])
+    # plt.show()
diff --git a/simplex_plot.py b/simplex_plot.py
index dd8ba57d..d518ca40 100644
--- a/simplex_plot.py
+++ b/simplex_plot.py
@@ -50,7 +50,7 @@ def plotSimplex(points, fig=None,
 
     P.axis('off')
 
-    P.savefig("dur_simplex.png", bbox_inches='tight')
+    P.savefig("dur_simplex.png", bbox_inches='tight', dpi=300, transparent=True)
     if show:
         P.show()
     else:
-- 
GitLab