From 75f47b2e0b9598f8055d68359df4e74243cdb913 Mon Sep 17 00:00:00 2001 From: Michael Schellenberger Costa Date: Thu, 31 Mar 2016 17:04:17 +0200 Subject: [PATCH] Added data creation routines and plotting functions. --- Figures/Create_Data.m | 47 + Figures/Data/Orig/Experimental_Data.mat | Bin 0 -> 45552 bytes Figures/Data/Parameter_N2.mat | Bin 0 -> 352 bytes Figures/Data/Parameter_N3.mat | Bin 0 -> 337 bytes Figures/Data_ERP_N3.m | 71 + Figures/Data_SO_Average.m | 73 + Figures/Data_Time_Series.m | 40 + Figures/Import_Data.m | 44 + Figures/Pictogramm_TC.tex | 64 + Figures/Plot_Calcium_Depletion.m | 86 + Figures/Plot_Compare_ERP.m | 111 + Figures/Plot_Compare_SO_Average.m | 154 + Figures/Plot_Time_Series.m | 62 + Figures/Tools/boundedline.m | 392 +++ Figures/Tools/inpaint_nans.m | 411 +++ Figures/Tools/license_panel.txt | 24 + Figures/Tools/outlinebounds.m | 31 + Figures/Tools/panel.m | 5086 +++++++++++++++++++++++++++++++ Plots.m | 100 - 19 files changed, 6696 insertions(+), 100 deletions(-) create mode 100644 Figures/Create_Data.m create mode 100644 Figures/Data/Orig/Experimental_Data.mat create mode 100644 Figures/Data/Parameter_N2.mat create mode 100644 Figures/Data/Parameter_N3.mat create mode 100644 Figures/Data_ERP_N3.m create mode 100644 Figures/Data_SO_Average.m create mode 100644 Figures/Data_Time_Series.m create mode 100644 Figures/Import_Data.m create mode 100644 Figures/Pictogramm_TC.tex create mode 100644 Figures/Plot_Calcium_Depletion.m create mode 100644 Figures/Plot_Compare_ERP.m create mode 100644 Figures/Plot_Compare_SO_Average.m create mode 100644 Figures/Plot_Time_Series.m create mode 100644 Figures/Tools/boundedline.m create mode 100755 Figures/Tools/inpaint_nans.m create mode 100644 Figures/Tools/license_panel.txt create mode 100644 Figures/Tools/outlinebounds.m create mode 100644 Figures/Tools/panel.m delete mode 100644 Plots.m diff --git a/Figures/Create_Data.m b/Figures/Create_Data.m new file mode 100644 index 0000000..6ce7b1e --- /dev/null +++ b/Figures/Create_Data.m @@ -0,0 +1,47 @@ +function Create_Data() +% This function creates the model data depicted in Figure 4-8 of +% +% A thalamocortical neural mass model of the EEG during NREM sleep and its +% response to auditory stimulation. +% M Schellenberger Costa, A Weigenand, H-VV Ngo, L Marshall, J Born, T Martinetz, +% JC Claussen. +% PLoS Computational Biology (in review). + +% Move to source folder(assuming it contains the Figures folder +cd ..; + +% Check if the executable exists and compile if needed +if(exist('Thalamus_mex.mesa64', 'file')==0) + mex CXXFLAGS="\$CXXFLAGS -std=c++11 -O3" TC_mex.cpp Cortical_Column.cpp Thalamic_Column.cpp; +end + +% Add the path to the simulation routine +addpath(pwd); + +% Go back into figures folder +cd Figures; + +% NOTE: The data routines utilize various functions from the fieldtrip +% toolbox http://www.fieldtriptoolbox.org/ +% You have to make sure, that the path to the fieldtrip routines is known +% Edit the below statement to fit to your system. +if(isempty(strfind(path, '/Tools/fieldtrip/preproc'))) + addpath('/Tools/fieldtrip/preproc'); +end + +% Import the original data from Ngo et al 2013 +Import_Data(1); +Import_Data(2); +Import_Data(3); + +% Time series +Data_Time_Series(1); +Data_Time_Series(2); + +% Extract KC/SO averages +Data_SO_Average(1); +Data_SO_Average(2); + +% Data of rhytmic two-click stimulation +Data_ERP_N3(); +end \ No newline at end of file diff --git a/Figures/Data/Orig/Experimental_Data.mat b/Figures/Data/Orig/Experimental_Data.mat new file mode 100644 index 0000000000000000000000000000000000000000..760508e37bcc7a758c5c47374be7cb5081b67ad9 GIT binary patch literal 45552 zcmV)LK)Jt7K~zjZLLfCRFd$7qR4ry{Y-KDUP;6mzW^ZzBIv__(PFO)UG%O%Pa%Ew3 zWn>_4ZaN@WWn>^tVR9fbI3O`FIx;spGc+JFFflbCARr(hARr(hARr(hARr(hARr(h zARr(hARr(hARr(hARr(hARr(B00000000000ZB~{0001Ju>b&goLpB0R8-sdrbi4y z1f@ZcR5~TjCIkfqEJ88Bz}9DA;q#%$Lq#l}7|#|J71U?v90b7vMMV@*>5iG9B>f-e z?)Q6Z{qI_{?%bJk&;Gt0-`+<-f?dBPB;?p_5c^0;aF;ATOaAX8Be5`IeS|E#55?aM z*~i~&o=5bewGmN^BqXK>{`Y;tw4nbU*rmWeFMR&9BlrL8^stXS?kj=ckNV()L>s^8 z{2I5xMXme{-?A+uGFth&8m8XY+11Jq-~X)LB)pYxE52I!gT)J5`2(wb%S3#(h^-Xy z#mn}K_>0Ati0>v|DB_P0Zz$q74s6#TgjLhdk?$UoG;gMLxC2qZWD8 zAded4QG-0Hkw-Q1s74-D$fF8*R3VQ_QL?*A0`eT@4(zy%@E3DyJhC`CJ^$fp$Tl_IY)v|EP!%Fuo}@+?Oj<;b@J@l+u13dB{3 z{3{V(CHkQXaaN%}st|8A`lT9iSEGMw5PuE&sRnt}B9B_+QHwlkkw+c!s6!rg$fF*4 z)FY322duUo~m_ zuJU;X{hE~XlNROtqDA+})_h)>q(!q`v>xT%)TGP4#WDbzRN?9bn*c*CwEN%xX6JV3 z+l|ts-=w=9r@Ymm98V4U^y^_4c%V+@58hDzrJ+GL+#d5{`VI{$H=n`!sg8WqsaI}Q z>s)tr`XIAkLwH@C#;$kBUA0b~jt<+W^3rJ(<@~2c&vfal8y`}k{NINy|8YZ^{_7t} zANG37FlxMRwB0CAdHU$e^0fRTQuOyT$+C09+W7mF1jk4 zx%#M^;gtl($isSp@`GUAOF^<;pdvU;VK7})AZV?*QF^2r{0e7u2dRUS<52m0|(Z?(yZQSHp&ke_Me$_CC1rFyf>eidUCy7p$7-K-b>Aj+TpFYf6u#f@Lgqq^N z#}$5224|2QDr2riTLIf%M0Tw56j`AxL2gU2~76^}!%IYC41Y=QJ8F|4RGJB%RYh%xTKdKzY%K(e#{|ANyfAHM_HG($%Al{2vD|=`QP*22O8A z!acj&dnR4ghwJsTuG~sB0=1)~Ja_7wL#=)rt4B-t(ezULOt=jgzh$@=+Hv`+;{ci8 zm>gw0!u#roCe2`HP&vZ<^JzTvR4w~)Bia?(l(tx&c65h|Lme!hgp2=O0I|%*dqDnY z=3jqL;0Lj7et3enYk*;ew-==8vwTYRf)|?C#|j2}gL4@3^K5Sr$8nc8i2WPm4T+`? zYtQ+6LtXaDq;GQGu<6eQ_nhx~f%V1Ob9)24Kv{L@+01vI@T|`9+BPpw5c?z9gR`Fl zkb9BgmjSGIcbMvRcSh7}H>h`JvXbr!3NcJB_PM}}2}@NE2aku^YkMso);fV$W;Z#) zxA{!|LL6Z4<@tk##oJ^3+Hw0Ywu4nS_&&W2wj6(5TX?_wtJb4oHXzPJ8EX*t4=pQr zV> zN@4#PLi1M-#B(?mENA80riS{gL05OOwiT;G-J<5O0@0`9JQ-s^xv|xdPI}Gq@#7e( zk{!(i#h8k9uh*E)N)<94HlZHlPq4t4P|jZ_xW5U#bViaDnkn57a=KI0WlC)husFMz z(IE4WZC-h1bZ90smtsz}Pq22PJeCR-SUYJn=iXm$PU({GuWQ`QX=9spfZ>Aw?qYpW zH>c}f&4}`EH{*0cWJXWyT-JT|q8SaF*Kc3`!<5!I8n7=W^mTa;>$k>qkY8BoU%QN` z-zp~KlgCi5Ul`J~J@3s}fzb&!S$_EIW8c-K{Pio)rAaAijiXelA(863>j`&`q`(R6_gzHUm3vw@@U~hX;v_D z6k5Glf8^0uP^}QuJcwpqx973j6O z>l+<%6{tc0pB1zcO>CYX<-b9hMz;7dT*IkFG=pI?oYrnrV~LdVmI8!27&)3{C=t;`k*ZkbuznI_^ z{%@XKaxIJn;3whyBmvw!APr&CC9IHTL8@*=(zohCP(S23E3hHJ<>xSHD`#>kr-;8R z0XKiCz>qE4%$_RtPc`iGBY~S=M#8~C?C0OrfzwG%c<~Q2lcf!1MaS6ypbPa;tX$^m z!ygu2ORQuJpj?}o|7-}lCQR?tjKMm`s*_m^6#ExeaZry1s*%~Z(G5# zAp`RO_C0mm|LKERS4-?U9SC=X?GK6_uRU=Bdw(W3LTBi&8NlrVw<4lXEa-FralSus z<@)n%H<;YR+T%MnGn_pBNIeYi);mSGI7na+>0mr#vr88~ueH+wyYmn9Wixt+tcPkLn zP_*Fs#_F*k?z1DzVZUw+8z@b|(36$lX=7+|+0FcJh`Oo|Tpw45@6F@bTs8y}P7VnC zV-h7_$_9#IRH}YJhxNF=oor4UN&=fD#@JEud=lb9tIb)upZB1O32cmbGl8D+WB%!x zM8)Z)89=2~2I4lAemua!@@yK{r^nBrdDB>$t(Zy0?h6X0OEM-hK+;T`sZJaF^-NU+zINp=A} z5Wj!Q5c`TLw8XNpCdvZj`OIIFtwF4B<80x<;U>*b{dQd3uR4HO4<|TrahDqpe-yB7 zdtBgc7#r6M+_*V&B7sVkpk|2)z|F}u0I4yoP44gj@f);PI4kE}0E+^c99pXIa4d8TI4My##JGa+j zIB@5m!Ekj{?7p8x642JZoXzh|L~!Te_}N#o_?|D**_It9X>t z6GhZ#6}mfeCEG`htCu2mJU?pD>9g7ElXNIISL@JmA6Oe0uS<8%Wn+Pa9`2(@Z+tXo zjYE%$=R|XTJdf&AaXxI-r(FH(Q{~@x>f0UFqulvhk8=8`M|sk$9t65{b3U8b_UcmZ z{H=rMPaRtFgdH17wW;?q)-LR|sYk4c%?(7oYwlt1m^25BKdPE)^M;wGTST$N6!$0uFT$dkeL)S1ujIpzzjX z@15<$z~k%r@ccGnW3)g)Z%!L26S_F4WVDjJhYhnVd_`os<}z8GQUS?t__I)$Cm=8G zx5Z4@E+l{JP5Ux4Mo4bVYoGIZjDXaX==#>FHWA@=?<%K@jbw1{Ch8R1Of)}LRAsIe zkgw1EjYhZ&Ny3^l+3~?bvgN=vll?mcBtCE2&V&=qUBp{v7f^s$=Y$ndeMMKY=Zy{ex42naS3dk$iCUOKtg)ISAlFCqtT@{( zAd<3MtDy#m5}qR>=PPHWynQ4hV=m5F?A{?H$9~?}tyU%=i^?ND zUzj2wX{+Cdy=oPZyw(Q~o~ww+j(g4y{BRM;YPD@TM1{ogUGytQT_MpK;`(IW8X@_s z_!*niM8xCxf&*8ki^wW-&C6PWLbBAsVTRK;0V$Lkw|>MTArU4->K4d}2>(ufv%+%` znf&Jcuuq|_WYPZ6*?C_@Wc%Iko@bRsWS;i65f36nq~CA2x6;8@a$-~5!s&vy3>(xV6?c26TJ*|fvkV<;Izixf1KB|Pjdlmf-C4@S+3mWG0)^$y#8WZ-MUCaGiW;4-vV-$cGh8ceo* z&aQqc4R7KPKQVqF19KZ*Xa9OD3*qbgBk3VI7#4qPN=%O&9Q}25QT0PPa11M+_4uVM zX!`##Ht3fO9CF>?XQd+x*`K`+W;~XKbhp_0{zY=2IPKaZhpTcBy6~w>^FOknRPb&} zF0)r!aC6YS#Zq80(#vA}IZ2qGWOZ)cC`s`2Yr0#RBLR~I-|ZDfO29Ua&4mYhdP%y4 z!$*nwpJc1D;pKk!4g$6U(fp2f;&<6t(cp`u%IecOOv&*YmIIPMjGD3*-HNy!KUxOMw$N z*zeQ;e)*eVi+!4KZf%TQ^9C(AAAYR;(Jw7%DsN0ot<-{|Wicb)IcmWUDaHJhNKHum zw&M8FG!1wfu=udW8+EA2ZV6FM83pT8j`x3IbE5vzXo1EIHSiCPzaDL*3ai}vGOt`! zfzpJGKLcSn%(kvk+-;%^b*V!&;?0%dWaA;bQX>WU?b55RPi8}*&;P#7zO=!RH%ay{WEbObCdzirEdV?Or=?J1Ij*9i{v*F`zFP~aW4QI7}9uSqE+eCL6r z%khB7Pds=L(oi_-ryTHt)+G6@lmmyhwcDR>mII%}GtE;SXg5^THzA&K-SXf4T&2^W7n5IXeim}oVo|x%- z!AKg+YvO)ubCiY>)u*G4w4`B6-|m3uL@8JmXmr(*Ck1-8a*34@lJI_G*RIZJ3CPu} zD9DFj#Ng%JCpOmoB)@u}VTocNk>_X4G91xMvdp3$9>3c|vUe%?S>^W-@AQzKJE=Y7 z;XALvLX94B%s1uSf#hy7DE?89v}rea^SvSG^6M_5sNhvBSl2}&)R)~$+TKML=d@TH z@$M!k(q|f3P3a+#*{kOT>h%)4;`V*-_wIG9mFTrFy7O$gM|N5O=>*VLG=2yEe>ULkSVLmqH;!c5Y0Vm zt-e>vDtgyJ#?BD@{BWe5yvn~`@4u{#EbW$!y85P# zIO%o7eVp1s8l+Z#4V3OA-EB>$ht2IEGO9W|@2j+vT%urkb6Okuc+PaANU4qF#oYU_ z(Wj04qph8`-?p7xvhAL@eQXEOwL4q!qP~MT2{+&H+u2F}4*&rF|72HpJXP->HxiOi z5)q|BLs3zY-lw5mR*{k>QE4F+B`K6Nghb0qDzdV<*WM#SX54Ev%xhd(zw^zZKYr(* z^SbAr=Q+>w`MlRt-FoVByhl5X4!8;^CUoHV=hhVmRl3oT@k79RF$?PV*~-3#EIe_H zeE;xf4?5gmT))lRgFj2puGq%v#>G84hJs(a@myH7C22tq%txP;`CIiMjCt>}+NmD+ ze*Pt#Y2E|5gY}xU<{s1@^)NC##)9$n_&3rKEEq}4#L?SX*j>RN!M@MJmBr_|my+X0 zGK#nig1eEfYMVAz+6@iC`Fq}S_aJoFn*%MmU8t9xSEMi0iFf zijU5@5Y^F*1YM!&@q;W(8`RKX>-0lC?cquO+XKk$G}{(}OLXO{ex!^V|KOjAM3%5(Ve2hZGjnDyTz2_B zPmb@0g7)XI_0|K3@=LYlu^T`%r?}qZhXXj#AK4xq(vM`L{R@NId-1KX(Z6+LFNF69 z_}x6$i_!E}DaA?_t`h5(>2~#C*}K!+Ki07@=fDf);FZ1LdmvccGOHJrO*`bwPm_JI za$7ZWfrSl(@a&IVy|4`B`zZ|?d<;zahCUXaomeEu0;5d=v-gTv=!tGlu7#M}yg>9E6zdK@Fc^ z4q|?d!Fty>9K?_LmB(L`&wYGVtu3uM2(zdMT+tewM2%pke%}l(qP2y=nPALCh;-+F z*jUO%{CTEu*4%@eIGQ5xxI=CV@z{CFx#8R?gpGxooyq_Y(YR6QzN#KCVd;Oy)aE1~ zk$3dL408p3!u5IX$=s9t#1}>~Uq~e%Ve-Ry&6ieQVux0LY_1$H5w*|OcXs?#BHVlB znrp^Wi31#W`kjM$2=C2Hd{665Av*h{R0BV85i*yx+8#@A63?R8Sr_Zauvcf@SgFJq z?mKQyw>UTk%@6v{wKZc9X4`F|tN-KP;rHYd&YlIuU!wZmjiIuF6)?=*CE zAveJsXc(Tqn3D*U=lX8;cnpOBUI}KBBQR}<2@SYCgjl^BjM@8#P^jqQEAKLlrT&6P zYS)b6n{d^E>#ZCF)A6e!Er62{4n8cPdVz~b>pm{?W5EnMa#4Hk7CKx#kG5Vhv6HW$9E}f2wU7jn>0=jL4ECG z>6ht$z$abtq_t=e*1ykhGnhLFTZfuQ{+HSKv8q&0?c+(}_dA=>rf?ZM!tjXaE}?N%iG?=|^kw;dlCjeXzg5lsFmJ z3*Hr59hnv^Jok_2)1A%2VX+mW@pF3d)x=0yzpfAOT{iG5%MGAvQK0^j_yMRI2X9)Y z%*Lsj5T{=&**M^H>y`7m0i63S9;HU=w^4DJOS<9!QUzaE$X{SXr1}2iS?dSk-yQc! zqHquz%Qi*tTK5O)qs<$PI|jiU@%YPon?Zam3gGou7=+_)Z@^<|n_grYDBN~;>_vzC zgUI(LeQ+!LWaZ54L!;RKmt}^1_`x{n)qatzuXedyp>QwU9Nat-(n)_Xycm)((t|no zJ-+YO?7^t^wuRRvdJv&8W^Q`E2k}=Hn#Ma`C9G@sF?PqC5g{`7?A%)e&o zrdPLVt`-?`n_gX}SvkfXcQ7yo^&HfDIJMOO_JTFu6Z262`y=v2?h(U=Aox80?+QW2 zEPVxXTej2ej~j)8$a2p&x&nxv8sB*#~bN+@U+j&Z%r=1o}qn=qYI0`@%Jjh zhei%Tr;-wEh-i>XEc}92_r@igUX@~jdpWtzGT4fSZEw8(73E$|`pV|zh|PR>>cx%< zi0n9R`CGLTtSc8B48*=+@8I^w!BeWhuHNxdvAl}(NS9sW)$S`Yr4p+B5`>S6RKtxv(|JC?Uu-@44(K;nFwQ^~Og1YXpbRvFL$S&>p9 zm5~PQAV-lyu{k%?JNkYjweCjv_B@)e#m|6w&z?ohtqe?`<`F#Q8UwRz4cS#4)SaWSyB9 zp1UxyFoKB>iGPxw6*FO9c4x!loPYQIzh)#w(*M0TF>YeM_a~;PGLABFR5FC@FB1pf zH;4BNGV$5{6P$X+Dh>{!gTgi{V$U^x%6KfE*8Stpg!=xo;=w6sk*lEbvCLF9NBFh zlYx-0VsfD;((skx-Yra2t zoKrX1murp5!gHrLA~kK@;v78&#Rn$d9hgh2dZst3HglOyTNSjbQ(=|Y|_9G2UrSFa^~kMFJPiiBQNOYJWcQt3mt`CXrz z@A~k!9pe4i*Rs`2=}tdYkA7GSHL;904xlbkBc9u104uGGHg+}-0NQz6pO&++e%*-E zXVPcJdCSIN^yaBBE1zCRm6HV6X_x?)#N zejs4ZcE<&4+VCH%G^Vvucxr}Ulm4DVnN5@WGd_+1vBtIHn;jdWxal0Pf=vTxk9Kr$ zG}q&+#^Yk!!CI<6)ss3rq0iHs+ziC7mP3L2xjGSFhMqKAvahAkT={V77PAuElzBq> zSRowOn~*x0hipSGav%SAJ+7~*I8;~OnrxSXIx%&!FY#2pqH$HPgrj*zID#I^*bOEH zQ+#;^rJlhdzqgMjcuS*x)5d-A1#R570%@YJ^-12mr%^shGc~X#Wh#>TKbkhao>q`&7nmCH)V`%^v}NbNslQpYCzeO&({@$6fo?n}kbRA2x3`jWo}mo!;NJ(1bu z_vjd04N~vZCiSOvZ4-%`%t=2V_tzY1o?Obm^3YP>o;&YN0fv0#NW2xHXWlgu_r=)r zMu)%0qlCiQ7gU8GBIhr~^0W&|cczw6cq~I%)AJo!_ej0Z@;R{QP&vZ)=aafp0dCU< zk}nnbYtLbo6u!S<yS@lG_#EBab9VDS02}QB$iwnOnUa6!hZv`-Ue)Qx#PyJY(R#B#8Z9c z#!0+;I26o#Q`?9W4);jhGSIL>YvkBjZO*3fQ{{4SoE6)})r6#-$P~7xt zW)mWpr;z();%}QPnUlIc@%?fAVIqn91X*_`jME~Qcjf)dB?&s#B(4~f`cugip!+*NR@6Eq`ZG`kYIRFM!Tbg;5no zF;s|oWRCjh-)9$V9DW~F4#9hmZ@&9c3D4>@GS^hYdboki$+eStbX-pxCvnI9pe5O^ zo{7{Kt4RJgV?$*N=?^Wq`No6C*szV-?++^X{hIWFab4`dEVKP&Ug$y|(K183qYDla zSA9>$cH>mFJ*k5|Pzv5-9kZWB;j|a84qkj%*Fv5Pjp@~EKlZ_WUv0}bi+&34|NM8H zM{HQ;c6rR#8ie8s^Y7vPgBWp=I`hZx5BvrU$-FQ`y*~t#`Ie-Q3{(3Xq3}0C<%Utn zjuw$RH44F#_engAQsc*P?yM4t^D)XN$M9D^#W;v#DYDG1t2hWreOK~)z(JUOG$eV% zK`2crBLfKsv5if_`4|T=UMCLX$SP8X^f(B;Q4vzVI0$J58As(gh_xr4k~yA(7|+>b z6rRU$)i;mKv13#nM=?Hb6diT1NPLW7RXYjaEhAJu7{-XD#Ae3@!&H76!oYbg(l7tu zZ;aLWL+y7ELLoOu-C55puW?OEwX_!Xa4k2=Z-!~H+$ikO|z@dXHhwt1@)U< zBp-XI{M3zxH48|;=%VsT7n+ZDlRni+t+$=>&)-=6_xeTYzz?dAZ4~}maqi}7^1Rjp zJBvUvJ~U09KPL2VLN_=K`^ntZfQZo4Hh+@FTO!cIWr>N?ZzEVN)<|{Q% zDaFr{$@Ad^zbAP0M&F!lSo&mM)}FmETx2YP^2KPXUSZh!z+KsVb^vuQcAxNB%8$Ki z%l(~5zPz3Ed&(~-{hiVU+Hsbz$k^Oe+8s+TiwLv($PCPQ?wfWO&PrD|h^yM+#bF8c zn?n~-9l57QyTA}HRD=wMx{tz=o^r*n=F;SQy=G6pz7V|&mQ=(yu}l{TZI+bc&32nS z&x_O_QD**lic4L}d0FmW__{*l$$lem9Lf(D;!gH~xqy$qXypr3Da0<=Y5y9!yUSjN zKY9bs6OtJ^-F~=Y?#i(=%^$YGKbL=R3V@M`({0I?K&<*1nd!gqEkagCZz)y^#)@a* zb%MO_pc^(cqxe_|imrQDe;Rv_L!qKP3eg{+rsV1R`*0|pt<>3Y&Nd8b&a2f=eE$Sa zIUYvH{cuQz3}>|MjliXQ3L^7nM?%y~w42S01TCO)P{=C^C*<3Dw(CTr!+x4$(MU8j zhJNixY5WY06p3AC5ixMIv;TZ;cPy4(ODI{mAP$+rF3M%>I2gDzj`_OBLtuw~Rq4hA zeA*$hw(n;GmR(+}AZee7wYh%okt2x^slFm`?LiVunrWpxipgjUie`-kC!?;j_qWob z6sTpIoo_yy0w2kZpPsW*5U0>OGq3K@M9fkfFGk*f=F<_ZITX+sPRFsSOxcFZte=X!vSqR!|KX)1K+KH*_qoIyu8=h_g9W*Q_$1XSKVpn;nsPru+U z4GUOSkxAERKwf#udoCJYOpShgXc`TcuTPDwi%)^OlVe5Qi4=6S-x#H?PoQ)7q!ttQIJBKS0@gzv$PuKfI427hn7`Y|F|C)R5ySPL=SRXYz z_*EhlO_ux3jZB2a;Ch+8%%9o+r} zy`$N5%nE;+VznR@ie{;mp24Yj$&h>9E|rEd4JoCaF==3>b<#~G(-E>qMEb_Zba3P2 z2mW;#NL8q|S@kXh-J7?YIw@!3ePK`DhsaFSSN*>IMm`Jo1YDzTKK!?55wF+w%ouj_DZx?qq#0))mOxLo{I{!133P2&v`X5PV9-8u)0ZhuTSyO(H$j-JoQEH!0Hm%eA?xnEK`Ef_Fmha93`+f)jX74Sd9H~BZ^G#Vytgl z)O~edF&g(0o?0VC__l19!SjS7@K`*4CT(1V0qLiqEWskkn#9Ehq!&UeX_aTQK_P;= zxTh$y3$TCfv{&qV1qhl?R}>U4K-8Zz-WS61v0?e4dFAW#QNwG-Hf+d)$ILZ1OwQ)v zs!#aB<(;{>Yw_@B$mLvk}-k zhEe9HM`N$Epd_7GDyWo&&o5Kul1QAi=53dKPUh1&0r7ht2V}tBaOk|_k__Ct*Qikw zo({&PI=4>MbO^nWnu@G6EKzr;D%D8Cj>Za_C6p_I~Vr zf(~}aR?V{{&ZSZ{5-)e9pk~_vMXk*#kdl$f5w;?8?&38nJJOQieCpxUbM=WxSCK85 zZc6Iv=D-Zipak@L_FG2ElXuZyvWl&<=hb7&wPUd&Z;Ul}dMx^sEKEYb#$ZjJD5dN2mR)uon?9*aTLOyd}H^BCBsJh&x9#Ngr5q;qU5ay(C) zrJZ*SDsD$L@i1brQu>&LuXHS;jTf%-RgH!0qnA8oF0o{;2rRqM9t*~jmfdEWaY(b{ zKJFbDha2CG88QuVxMl3<@Ps@EtNDuEasC*OBXxIX+EVrtjqq+b2a;SvVeo z`pP#gn|~Mjx*ylLU`cqhWOBYS-}ExHGH3tcWNdfx}>dBp+zr4PX4`fQ)v z;Uwt!XcQJV{tzCTt!o}sdkmhNqra*eBtyBPtmo77DR5}i@`%^1C-9$exGihbQ_%1U zJ6hiM6nw^bGP39yXxdvpm?(S(+{52Tw0u*+|K30Sbq%TT@YkRKYvFS!b2GYs|H^Zi zeb4K0r$id8_0SVF3QL1u4HeC5pVHuD(e}|-elMW)KHQz)zJTVteFrL6zl7gkBNuIb z{1PnVLpjPC=@4@1bo7l^=}^>T$9Kc!6@>kLU`#Om6(|hw1tqAyhWM-M%M_Dd!|2~O zo8*^fKv0`ny4=SMILdqeSLTLHFgujj8jzU@sg{Ey3$3ysq*f{6L39>~iv@0+Wtt6n zUVa<PfbtyYOQHl%GjPz8cn{%LO=Uj2s z@*Ma)DB897d=6-YODk5meV8P8weGvVPxSVH@1x9&YZ&(5S{|{k0mGQTKcs2q!aP4YsoLONNHKUJ z7x^L=5)PC~Yf9(A-e5&;qkA5l7oXH#mXim(nNOk@F3pFy(niUZ7xMwq=Z~)uD1fb3 z%C*;=E&yF;aYHMCLXfx-c=NAYg|O)Q>^Exj3YWg9AafKvuzGJAKM z{3?Tkih?gPa?4@su<+9pffev_%>RzNW+i;jZS`*Lsf0$MO;yRRRbWsr|F71gDk!^P zto{2}6`Wxm_O)DH4ZM3Nw*=p*2IVw5)xX7AvIc^R3gi_Nt6_TAj|3l1 zH3Y5vd{wW#3T~>P^&X*B&{Lki*i)hkI@&V^%Umi!&TH**y;BtsY<=^C=7Vxj$~n4v zb8i_O@Ah5Lre6k)x4)$J1-}LF*juYF7L`JsQ-hDhl2TaYEaK65`3+PUB#j&qeFIHv z{Leg&DuJ4CN3?X~+ZdxZ@N675aLqhSpm+4P+0B!Rmh|eLbUE% zs;U>lMs+3E*}n@w<<+sx!vh6ycqjTttEd3nHrzg5np6NuX3iPKi3RXuc5JlD{sQni zJ_4s43n2JdaGkqy0hk%;ORe~w4}9~tsPI?h!}r%MiaQeXVP(vRQ!9?;L+BctJshWe z;ORf&28!mxf(iHSsY7`%SK-R&hKf9JdEhsbfBf-qE@&RP@|4q(3wE#NwIe^}f_Z0dXiQ};1il|Ao5Zklg>F{$^2seB}Am|2OT) z&9v{_4|n%k&i*jWO>t;CBAGJGJ=NpN)TjNWjF zt3U@H;{KRgHZsU9FPr)daucTD-5ulx(;E(QUFg7rTumC_AeZ4WNbnC3+yeyf0KwT$ z@bwd1{RB@R!O=(J_mQ~0BwjCx(@S*fA^P+XU3!Qf-9(3O!oQpF?jn4<2+uCUuaofV zBz!svk50m)gYf7eJUR%EcEY2b@MtGI+6a#}!lRAwXeB&a36EC7qlNHjAv{_Lk7mN7 zneb>PJemlPCc>kM@Mt7F8VQd^!lRMk|3+{(5WEcp=U0O73&Hh;;Q35&d?N8bk+>g8 zym}I+p6K*}=<|W-@}B7Nj_B}?@UJ7hYYE?4!n20(t0ugv37;y$ql)0KAh^p3-ZFx- zl;C?qaFq}|#RNwoiC;kC=9752Bn~F~=a7A~$bPTMKI!E5Y2W)7lKcNkC-*?U2eO{rCji-xJWl}V!SI@|EEI0wkXa?yo@do#1Z(QV)bjBZ<>U_%xDujf7VdiQ7c@ zHIewugl99s(MHiEN_=+Q>-wi8|23GQ~HPY1!@L3HXM zJUR)FPQs&;@aQBwx(JUh!lR4u=q5b636E~VqlfV5Av}5rk6yx~m+D%yyc&p+8 zA&tM~aN4b9A&(5U1Ja%&vabZcx@{iE6p=% z-ek+hUVk@U6I?{TA7NqpH#gJHyakz=&arUJbKQ(iWfs;szxKF=>Qd~uyD-76 zYzc<%_XQ3+>fj%5-{trwYvG^smJF6o*TfT2OaJ?+r;Y=!YpEqGtKv=E(=YbaE93K* z*@04P6|vpik1ltX$YLqyTrrvRQut`+npr$aG3;2{5g(~5h)Yqsh*Lk08z{Z@U^V9h z$ipb_@X|?Gf6vW-MUXJMu6*=C*BMFFaBqLFk%J5p|CFe({fGim+>@1xwUrSafB#~% zQzhu7>wPtJTdh^oxXeKfYdR@CsT0it;e8=+k{zqZ#eF-EH;cvltIn4o~8%ZegSnxg3! z9HZ{Hnj!ie_F5oXcRnk$M^)qR@cn?6@KN}@5NbXxxN+SeB5|zIz>iiDrwMCxpwyM( zXM>I=cb&P_V1sD-#V7gt8AcS`zl@NV%=*FvDMR$! z_pj(J>H5eBhvt9s&|`S>Sg3N{0xE%Zkt>HPi}~7!&i`9#$f;r#1xN`=&Z&6z`L-0A zd`;;lHcGu`?W6z;;qgjy^qOChT9*XFxaX0Sw8fE5YfmPi3D~I~sv?jd8K9?=n#ZTj#*k;B*pwHg$!T(@* z1L~1g8myzp5{Al0Tr~OK*<+*XSGArrF0W#7p{E!Yqrv}d97D4Wmd~7naiplu*Z)y` zfB-vz#;;M5c_Svfv`aO{H?ePq*)m77z}J`2fq#g@$JVhdzYc<~hcdluWlTyH4wV9j zL9~1!43%fYOPw~w+pbsmUKMA>NBvyHH4dhf8^Qu83mJ9{ij9~l=nUmQ{n9TNQv8o6 zN{+pYwOgpOnYVjnsNpcKN;(l<36!4FOyd6O(CL3Ri{|||ZnMvA6#sH@Rn&Z4EBfkC z^QU)Scsq+|BVjydeQX#a5lr0pRyeDmWixuRsgE@m7kiI(b6VVe)tvNS}VEKHcF!1P`o>2tdaF6G~Kw?*c0p) zzAR{B+De>El7uRNm7GZ+i~T(la8-qo(j$d#R9suY^cheN;eXM5KY0aBf;sq>m~S@i zVudVJoDB&m(!S{Q^4{or_%5W5NV$q1JlyD$>0TK90vowp-Ib(U`laboTvH?;pSIPC zVoJ7X!w(c%QwkTqN!TwIoKht;bYAPLtn#NMXphDJy`s=<+J49f4Sj7wJe zH{H;-H0oM5=NXKu$gkcpag|@V2MLvb*d;)ggKtd5zY923O;wB_9Q+P%2ds}}uCNWl zu(0;`D84q(M_G;u2wR8023NY${}=|&2M3yea17tLTh0kGS;TMmYQprUa*S)fr{P!4 zvyZnJB4r*-zzeaTwB;DnUUl8HOP6E2r8_?L-AOw-}BTVF#Y_p3qUQE_Kq za4hp_-ilJdJO1D628N@7%(8kLZ(*!YhU5YLk4>_W)LT<0D_3fXpR5!P@}1&x3_%;6RZ-;J@F@S4hG-s4I_{tQ zEN@XHcGoc|zjlGIlsi-Fo|BW*o9R7uHZ>upl+l@7%7~Ld?vi(TlIcq)OL&8naZb#F zqaT_~6c0gA*Dr3qlidkOFOBV)u?@)-2c}%nFZK@w09fPb%x{6XcPxsdoXw!rhA~@) z`Dc3Lu0E4A*7igs&C<~Jqh~_gXU2VrD7VFav$SudB2$^yo{VXC0@MutD|Dvtj8lQR zschTmv#D8qv#LXNcWMAfwW0SKs5<{h3XhHY+Vrg+1Ub^1=mhzs;co#E#3f7Sx7B;R z+JCGd=U)9riL+!)_*hn>*#j~tVu{`=j)T5sW;Ks>5Vu$3+%aFvFN0<|btDHsuc8Ef z013brl~xX+As-8QXh7Mz(pxl{lp~_Vv}F?n;u zak@ZrH^>%8*Fxd#SziDyW6 zWHqmMk@}>5YT3FXH$Kxa39<5fn8a}Gqv3LdQ67~OthmfaK;#`8udUlPm-@KWN@axj zqc*Z_c|mgx@X>lnkH$Q1tX=s87?YUaHH%or5!zEiPr4JNB(C~V7kEs8Q;h9yIDkXY z-_&H{x3Y$XDjlOZ)xeb4rT?R>oA=I@J@qv*S#>UTzr=3nWnIM-)ctO^28&Gz!d8B$ z>?7}bJKDOJ{tG8zC|LTElMRi%URl1~uDB8zcKYyD7INy{Bni#t`p*PZz|UDeIg#Ru zBdMUmpg2m)-T~4(W)VH!9nH)-(9X=~_8>U;a8Ua!*V^azxe@A`*fiOgkGQ<6gOR=Y zE-n9zxx1W6%}dV?nnfg9S9Kky!MFQN8A@uwBUP_BE^yQIi*y`suc+9*O~9#|JKAQu z8ZWx$xn<>49PK^`A)*UNy6&)HbPt%wJ3AS5&=OvsxhRSz(r5+wylidalKRbS?*W?d zba(SstGs8z?#A=xyhb^u@a%G(l$w?A9HN$awqhAp85vo(g-0?D%C}SNGOikRgvugP z0SCg3D_U+7hksw8H9E97XiN8N>4k}+qXDK1B>Buw3G@74E$Jzg~nc2-O8-1K(ep45JXy7ECYF9CfVZ8NR&Lw+Z{LG(P zFKp3YBGxh+1_6Jy$s0UwP`vri?xMN(6 z!%If57NCM{#%FrbuQFeUR?&c;!q(F>-@eu;-UL-P(z7y5?J~Xf2H}7Ke&gBN?GuE8 zwTt;FxCx0&F@)mkn9 zue4@6wId#f`H^pC73x4za{kf`X>u=G@TZxaYid_vU5hMcMn*0(puQ+9woU$#iPKD1 ze*Hq(d<2B(?WWLB5|9yzQ~bf_!)0>Zrre=oQ1FX>oJoL&fs&Yel4(&-4v{4o-c#x6 z#a~FMZ|^LjQ}DjS=KFci$aiyr>UyaqR$n%cbXhlD`jU9?3}q$X%E-)C-7M&jGoV~? z2Wb8l4-w=26=Q7d(}zua6GF5u{}pdyioD~wsjkOW{^67uN5p+mYfw`>4syj8tf9To zVipwqxWW5mluvUKbY;r*=Dkn2T#d!mCuLx z1oq-d^=d>n-WIXz{{0n~{=?wt zD})hy1+=(rPN;Uz{%(w)`n#8ycBv19kp2g0Rrf?vd;1%NgP~iPtk?}R>EFPRL`O8k z>4=}f=H6)d2hxGi1^d%-oNuNLVz&qHO+{>L-*CJ3k87PP-3(dWKS3}pG0UBJvsvL= z_ylmR^!G(iO^HnRmuIJ8Akx)k;fj*m4o7s3gnP<@t~0#6kaw%MN!Ass-w;*Ymx~>% zbidcOvPs8tVS#U728^SISQ78{ALk$B!B&$yMm$AKkx)}C7tEh@zumA#C|pco5@V!J zPyoV^OG;Xj(|LPE&^fSj3KH>kbVI{*`ME>GyJjZ&D*VSU~yEa2oA`-p!UlVVk zs1`Cgt2$Ro3*;1WvjA~cnY8=V^A0OJk9YyS|WSxL-)RXCf2{z)i78*XG z)4EBsw0vL9hK$)t!263(QAV+-W_)*aD@k1!pJo5PDIkzLhL!{|of1|W;D6@Bd*1aY z$NYEIU9Z>ajvxCFUggM+Cb?g#5;bgA-?<2_q6<<6+NIbh zq>Z}nyj#gd^#%#ey-$*E)BIo$8`3}|Prdjqy@*Bp%5ohY8qb`qR;gzu@{LqW@6Dq| z=&Nb@4rx;=d{T4&w$%}Kp?Ejy8=mCAZW*QAoFc99&diEKp{ zJy}S>fR6L;3Y^2>r7892D^U&x?@1M@;jwYaYumbXfNGd{|Do)EI}(r&aQomb0mr{L za-ZPwndt`Pp6yA+nI9>-Puq>l*0+BKcBrphKIYMxd`&!a9oe6Tf9?>Kb05*=d2{^g ztAX9lUjN7=*bcS&rjX8kphnN@jX#V(PhNeAls=j8MI< z1^#iM_?c6yE$65i{A8(2RBM42OlZrl?3MN7dtW;DXWkQrm19%(5M~wvMHeS^}B~%#3LIwHUT1Y8=v`dQljId#di=A8kDlav&J+t04VTYF>ra#fL ziT6?;8nsDSp=C#EnLKo&Tg!@hv(u@+d#tqj@aQ3duIyz}x`^er(whu_7zZ2F=SV>j z3_<0uXP&myvaBu+_v)hsDhL3t6*(tOB+g6W=(K8Ol%E@+v2 z~lwXgGq3PAX|m#wnN0wtcYpU4kt2zPOWN3kS;dj%(hZ`(6lF^$bJ6 z-mIS4|6Y_&Pjo5`hePiE)J5koV$g{=O#!bgX+vR&q@`LYEs5lT+d&P zdiYMW*vWr47j=p!k0I8cSkO7so{GCKpMMpLL|P63e_yA zf+wjob$_Q>*DZ;Lf3AeiRe&Y&Lu#qv9lWP$^E!+=Y14~rsq}1KI2{Cj4H7rN*HPbJQ~x&Y zVy$O!f7?NVW9jPf;koL=A+UemGv6z$)MCE8h1ey$5V89AKJVco@lz^+T2shh`ZIf% z=9f~n_!zQ*heNMQ1kxlRqu_Y8NOVd%ZHQ%=d|hDuKp*)$@l?0k0&j{hlgXB%NQUtcS&beqU|37vB{EKd5qTe`w*2T`hrl=mUGqSc-)Yj4L1R4Yi4aS4F|;~zitHsT`Dt0^Ai-HWnF@Ox+~+>6cf>vrZv z^Q>Ct@FU)*&<>ULPt~ZdCEQu0Q}u2xMIyL#^N6grgYJ=Yo_;mjNc*BhFllCQgs;pB zczcGcGiH~oXqLT*#`APj*u4JWl~w-vUWZFA4gO*0y2){PIw6c1^Ayom=-}zQs(dtC z>0@hUHL*&ww8JJYXsd>;v%hq3Otcyqz0!R|OUO~o5Cb|LO_p-&_&k%ou9L#Sa%b|S z2a`;uIY!OdiOu03r-jx!N~!rcTuT`v+q@|2T$9;~Zde-%zj{ss{f)!tZqw|odVvqb zDVkTOV#imkjdHXMG{t_{Gg^6DjZ8$G#rRkGKb&hS(x6)UeUpsuwISpAr$G^8u4F$O z%%4B(nBRXk`TS^8!|Y3YMx$8o(`L#1!{Y#e#Ky1o{yBUzBgE#nn(`B@yimRQDP5)iU51CsQsnG zfEAoa>(7Uz+^e z7Kz?>E|+F6$fI<(M;T@Y)MrwIH5`(VUjq8&?u`K?%jxFO!dzKl{_~MF0=rC;P%=;I zEV7O0;?S!E!E)5)0ljNr*57Lsztu?;XHBSlPg8`QhBl#7?6?p_HB0mJ%!TO3YG#i^ zoed+O_5q-DD6It5vkk4%e=#`yh-e`&&Q21VMH8p*_FUu-?6qR6i*M4X?l0$8nib7T zbpK!ZB?GUs7smqnaIzKDWv4Gcjn#*nvj0Uu8->}6nyPM)Ijv#4a{BlEz&Mzf`yQ;m z_b~1D3ZbZe+sgiH=euT@<*L7?qyh;Pp}Y6-SxmUd>|ow}{X75)ZF{GkgsU2=k-}fM zD32j?6Fsk+Mgq2Y2-ch_YpQ|mJNZD9KdWDxd>k*AA8RP1!VF6u2Ldxc)8^ZR8<~}W zDN_6y8brO`b$+cMwbz+xjnHQ`j1F8+9LZQDIqRa;Mh7)Y=^Pu6uA4jr8&asLBJVVB zSallMLp9Nu)39wTC-iU?T``jiOX;ejyVjc(a@#S-KictO05txOi6%&3!IQbCI}jiJ zvd#m9cV%4gmFsN32uniz?kjU;k{%Ox<4=>rWWnoQ==-q|h){XQ>#HIj1@rg3=F?pm zi`Xu-M0`_b0`zS4*%cPtLg|n_NlIVot<)+x9fY1PTQPfVaZhfUZ5nYwM(9UbXMxsI zWGr5IqBxw!0pso8FV6Q@Y0t+L1gS{~SDmSxZq3mGBkPcVU9^wiWpBy)6Jo!Jb0w_` zN~_x>O?1%=lDNYPi$0(EcSjtsfJFx8mcbX=(N74~d5RLAC~T26NI7S#D?OZPvtR&| zosiEUFm|Q2(0J-|mDL&%)ukT#OO1liWll<6gYoa)>&2_jgRb9gJ;=9lJssUSpIm#o zPq{($)kJ1K04aa1D|6I!#z8&jdj z(ozCi2foYqA1soQq=Wl2k!qnfv~g+O<%~le^=Sv>3CQqsa8y?I$52JjBX?X*`+Hzm z@k%y5Cm*&t@=qf8WLuXNeU!c5D)C|m$;Ow@di~wQLJT75u4R+59zD5cXU-y7&_Z}E zv`vR1mZ|?HdOe7}GW^Wf*Kvv0iv9eXE8yXG+APo&`e%ab`2~}o+4R;t&BTE_^31n$ zihsL@-}ra`F!X+6dnfWqsdz(2(ZM<2}n?&P?==8}eN`Hv~X((~Q?)Urm z6S`vSb_Yz)d`5TI=-pnPZ$m*H4LXV7gvMuwG$vH~)m-VhxTAu7mV!EeXR4<36G3 ziMjq~8*XoYKF5)J6mgox2#NPlMhWkD7isHyYvJS89qo_vY~_h8nCk@T*+}QV^9%e7 zvC&+!ZW=!vrFO9`xzsrC_7_u=+T18fdGLNIuoJx$>R{&JYfPV$R*fQGe0|Z8KCDJn zF!>RbsnM}r9uMm}o(=iH68d~0`87u?iqFzTL^*_5vu+guy4-P4p?^%{8JgYbro(5? zo;1Hb`n*1X5Og*A!QH6+4UcC0J1?UaEQAz7aP7G8z}s@y_m6Bk!q0C3x8UMHtBIdK zE)8DN$1WDh==#pt#K|DQ^}Ygs*r1$}%hqe#F$Y!@!7{FzPXRAaVRw^8kb3?fH8XUg z5Y%|BLl<{x%aJ#cjw=z2hqG#^As53M7F0FC&c_+6u4QoyUln-9!BG5KP^d7#SXg1u>>!JX=+v&?KxWiDqu zM_DLKOmDrea4@Q4S^9)z!H^}noW3ZM4OULEA3aADJdjhXP7k!f-mM}76}LexpBK6HS#d# z*y-+YxD6O3xtRu<_jrMsAd;`J3#mn61C`# z`Et-Y)7y1lH;Np!2^=kR5O}{o&k^O9X&UHrCc!NvPcQ34Pdq;UT6}whQv=wSj9XV2 zBsfB@zR6a)on5A(N|!bUB=I;>)%a>~KaN^prGx0y<2fnsgk*Y1rW4H3KZFjO-s$|< z6~Hu*)o~Mm6>^yT0*DwcOkGx(nY_tASS(Znwaz=IftJ{}w6*LAO8pHHxK3X@^r632 zIiN0n{}#+p#8!xhd12~z`DF?X*{GkrE(f*do1b2)2p*;}Njcp0GC zu-{B!ELYuhv*Ay*b-0AiTpgnKSu0B}Ov_TtYe>K^@7*AB)zrv$`TZyhh(ger zYjQ^xya0Q&Or0aab4}4O_e$ulRxYUOAWiUo-Yi%{m9&BCipQ4hk1L&l9w**n57O59;&S2dR{u(on_Q7 z%>+%3Y;V7w*dZB4W2wX*ynkzcv?TUAm+{$GxyDW9!V6Y%e8AQ z7}NOrMf7&pipU;!QYfL6^w?zy$4gANUNM-k8+~ztBKlR zX_fynY29PQX%X+L!11lPDN48}fW|jw#GDasE!N})S4PH!4=N#HL_Exhzim8;tZrvO zPW-qKU%|0K=mpeAs(qX%fcj5$>3Gg= z#V(bQ1i+`c`h%SGQ@_`-OT!~@8=t-G<5&{*{9=Qu`19slUyJCb9#7X^btzXyHHZePE> zP*#VayjvI066E5;J=??mc$QhpWJFn!`C7wgfoqZ?ZHq>?bTar>;|q6FFzBB!cxUM( zMDUyDK?YfhT}u;X!qq$>FA$;eGt{Yu0HAA-tTemHv4Ogy#VTgfSe6B^lKb0YGKMvqulF_G82Wi$L-1cY7W1rHp2=lY=2TgAB^zZ9IX$3r&4ivN|? z5!kbA=;<#MKXO9r{E4vTh@b}#Glv}HKKiszYWFQ0+o*Z6(FjV{B>W?8 z#uP8FfS$0xWs2AB*!~k{wHHae>esW-?qtZ{J)|cxff`>DIyg&xHp(9m1^yldrj5FY82StvdKBI7xNQu0;bHc>3;}oMudsHF z-y+9>YaqbrkWaHvm+)bkwqiBVa9~Y1*}mF7qpVS*ooLn*$;7?G*k1m2|dbYV3 zf(wvv;P;@g*GZS8ZkW{=7_bcDQ-sQz>TtOw#x=O4jMW#~))ckucAj~E7$-3}#OCs5 z4CPCxNFwOJS?sS-?^QhDdRh+LKYwAic{TL;{Plnz`0_&%NxVg)`7S+JBMraXk^XYo z@(NQi`hjlHJKboPc~D3u7+eptnnbYLOzJI?>LA?s=C;VqalVhd`@}{^_2TShp zU}w(52;m34${6)wt$AX+e9G+gZ5f{U=TaJNOrU;0eH-uV6ZtKADWHK4^Mw6z^Ze)y z^~?GNniwcrH|(Kb6?6;T#@Vgp%f0??$RbEuKTfAd8I?o%mq8^Ef;loftmn?l#bZ z@Y(~K)y&ZiH8W>jv2?{`_^Hqvz}0~um}%=LN!^O@tKE7wVTN9e>l-t7-AeP9@U>D| zZwCRaGdC)LtCipv4%$smQwu>-BR7fR1B=u~b$bG^M>!UWPs8_HKb@2euR45K?dWM_ zt)@Y&4&wc*@tQNR21YrmsK~NMZec%^<6xr!%kA&D{@3RsRw~c7oVuD}aQqDW-g;F= zaW*+hL8f}NorhpO-L_gZ9!zF^axR0UbT`yFw=jz&J|MgK$)8>__B5&{nqGW8KijX< z`P&9w8s2Gl2g+3{_1nfu} z$$MwXS=_&?r_zHzk8jRnoxH>oho*kLtGkVI;@eW|Y-uu9bGk zc~q&1q4Obm-N%x&Yb&>^{bbY6e9e+2RbJM_Va} zxUg-}%lNXkOGA-@G8jJW(tlrrh7W$+QyCd%9w%op2O*bQhz0b(oWdovO+aIv@*aMA z^u@r;ND0h!)zLz}!Z%&9cSRT4)@5tBDPU3^UiC4>(g#Ov&I{Lqn1Yp-Wo$b)B#oXr zH_Hl(FnF%#7RQ7bLph7;{S;xESF^?i3^VbEc9-C7l`{y?4Duh2Kyqvxf>o2?E=`ml z&!pp*d5q!LDsM9zNSGkq3YO_UN+{v=U9_YC3c9r;t6On>7|YhVYF}A`wU!f;=xcSj zIt}@(7jH^vKJbp;4A5vVF(sqSj}4Uxw(7sU0FF%S*Rr*j#f1|KlE#_e;UCYfoGm;Y zUYUAlm2-GvD;qWI;hLehPX1eQp8zqk672UKUIfuf&Xbw&tV{ndj`pA4(Rtgbq`4cxS#=8>i`ROTAJMnmX$Wu#zP(81buJd>+nlhB?jb+ScNB*Ojd}GT@`(4Xr2!)n97j z1u1m|ESHZ^L>~H%=U>nTOD22OmMCDo2>M8N>ZQUYLviFsBT*Ng@qUQ~DYKc?@g=f* zxUwycJryJ?ehVHnG0#%JN**#K%)iLe0 z{JEF$d+=A9VabNGlq3Pm5CYFB%=U~*RA4Y zD=4%vQHC=Tr5;~mCo9Ukesk6*S?06^mRx#VFIymbo_o6sRZ%LU1#iZF5354hRWNWx zth`WPSv5}+(R_j2X@sY~op$9k%l!_iy&83gcujr#H1Ry3-bKwytriimppSxo@bSD6 z_&uD*QFP^S!h`_7f8>2(%7is)rvWHZP4M<5K$lpbA+ilZe?m^CzGdAR{=Bir>cG6? z*&iTAgc&vKaIj$v9^i}+(=5jmJ;5$3zJruZw~BSVgAg4fL{R>4aMH)*ceSJ%rHSq$LmOd_s zpgeIM)OEVW>^0cEmH-PRl#T`ZjBHX7S_(>izfM@7HgRs)>U-YUG&93DwKCN7a&7RRUES5p0f=wF53}|n->}y+wC2a za#(ToPOi?BkYAm@vCOAfBg8JHXpJp2n8g1je_iD8Yg(>tTvrq?4>0frBN0vw_>Op^ z6PA~P7(Lb-Kk|?l12Tb5ca!TOjG(nws_-AjclcxE{t$2EGE8s2rMDbNLn{WM#s>@(7T*a8Sw1w`XlUQ2&8z0fh z+fZbU>4nwHN?(|;4>D>XDc7&;7kt#{n85KKgd$P*oBX^^XO5bj?MjYJk>Q?^;KCUk z!k?};m0rl_dhduD>q?&-30b}IG=QRFm zd-<45=EtCK=_Xm*IzU!;X$8QoHKKIb6* za1v)v5UhxMAW>X&D=|KDydQxjebT!)^}@lI2uT=DCQ?7KKK;OPaOX}6Kb`hHI71<{ zn6csoHKc}PL|~p{Z(6+=sXZ7vlcE}MH}z^LnvnVW-d;P?plG>87htxbESagX5H}}7 z@)sD~E(CNO^s{_zs>DFGLOHg& z#9n&W@NH@~fJC>FtkS^uDVcSvMlFXwRn5k;bAA~#Kn9pm^_=9I>)@H2+2?hF+=#zJ?z zR^F@-EMl^_WCiDO%x1fc#@g`ZVzTsFi-K`A`riG7vJ{ZXW;}@kw~;*W;gDip$j)T@ zaA`v*MUOAVdA#5h$<6PAyLs|3ex`VD^ejrj>FIc1z@g=~={srtSLscsfxbr*p}`o8 zHSb4mSKzVSbCHcv6(%!CeR=ce{W#_!0UF%x4&Kb6HC=ZmjByo{gJn%LNdZ4}LQpVw zmz6hv9kCV&6R9gjZC-*xz8els{ zq{)vclOi!h9%H9j(oft(E0Vru2j^R=9vqz0M{+n(d^FZVE}W|ZOPR3T4xHrN*M%e!l2^9B3S4^>$tIe1Pi#$k5K_M;{g=lm&AohPpp z)!nU>HdcBqu)%@bA<%;Ay}V*;&mHAtnw@Cs!l;GR z3u+7^VjB*JhazJsqjmh&jM@UL?g)g4zj?p7s&m6$kme}F+ZP&D$48a5M+O_g8y9$R z!G6J&9_f}TH`tt)ITd&zU1D*@-i$T)3x3Mgpl$9l;(l8HKr15aE8E|W?ZpX=@#oIp zksScA@uwU}tN_OwDP;#rF z#NEpP521bk>1Rq<24b`dU($$@0fRK2zZVPABsV4IJh_;c>Tm_6_!8&a8K!^rFr?LfXeZ&8OA%ith`LZg`xFb8t14wz$%BRpYRv}`NVvF@10f&?|tApwNjVjv_Sy$C6(z8^^G zw2<|yKGRd^I48pU*d(BUtzI#6UV7#LdTf5ILz{pF`yi3$i$9ke`ouTxp=ac#;@@-r z!3!R-fsJXe4a9Q}+8K{xuMK3;mRCoagA2e%$Nm zv}&LO_q2N8CGE6kM1TRP8}_PJ)Q$A;7XNvZx$9XY&UvgI_gFdRXWPKcndynTQi1nS zKgzdUT!IF0Kdlu1?zeqOd6lB3kZyhcW`qnDIo%EpPoudqx1PJ z07K$Q+(3mT@POXaQ%PwSJ|7@3Cf>m9U25F%e-%(s2|B?4^i*7$l<#c;7!e=f1}Z5b z1$>;IN=pm#!7Z@;j)!zZDKv!uK%WvwjFa)f0ofwsMcq(}L&Sd`I3fIr_^&x5{E7Lm zxgh*W_^-Jl{7Di2`LKtf!V-~qJWO8dfJQzm7oPZLj$j?|BeVtWM<|ElH!&UbZ(`Z? zDCD&`qXAq_6U8I`w?qQCKTVX55a2rmk|&^fyG@jc|9=wBLeXzD zSU0*taYz8tsgg(p6yN`rqED41A}IK9fJKpLl5TXx;+TN{CrQ}JS@YZ$`!Hi-S#%wW zz6tz`7vHE42xK?Onw~vqUOL7;PgK1R&R5J85-Y3WwO63Km_zGj9fH!w%WQ|eyR5ed z3diDMfcrj~UfLWXliGK#QcOOfhTh(yUjl2H5L;P)aPK!!HwI3fx*S-hk66S_En}1F zV1=hJFH&MLMg1wL-Z$teYpI)NiYl?|h+jpH>|{J95O^~e%_Y40k>SlP?Q z#&q^w8XB9+Sh4NOH*woX`-${NHsD-Q+~<0!xK;bT257HjAj4k-VIL3wy<=O z9vjS*kzqcpu4rn-r|C>~%+Lq+@hqniM_XiHiFx#-hdw8Gl528i^#OBE<>U43v)=x! zC%wB?8hN9YrW-!!YlTpY#CQLCfN-B9fmu)4w)dOi-)#6%(-*%C3m zI(Mlzc19N)GfRFS&-)NBHIw#K6!I;pv~wlU zc-qfIV*aX@ejjf{X4DQ!uULLng6RDCL720YQzdmH=*hI zJTnbkE0`Y|fo!ct_;L;|`09F80&;TyKZ)l!7}oGwDKJJhD}Ji@<=11O6l zIwBJ<{pXFJ&x`m(Lmrid>~>M3qSe?Zq^)V!Yv6oss-myivgH4lyP0XoycIXaNkCB%mRF{Q(t=4Padh0U5yZ<8jmW22z z(zdCwAde5wLT46+LU1u(d+oal2V#{t4M*Byqd)>xS)wS6G^_U6(hA?-HjZ#j75;rvl|3*~A$RpuU3W9l7-G^g z3Nw@aWcO(8KG{fgVO&?pv6fpI=xy?BzAQZeTus`_syf&sxR52Z2tKTN&nC^-e&(7J z?Mt{kX)OY+3Kj02A7 zqNC}^dan>PE~%gHzjZf`ImAXqDVxTzd6cj{?U=B1KkX<5_mY`tVw@H4UeCn5yAEvo zwmW(R*gq{-_8Q1TuxU~wn0x)ZQ5q&WC`Bc6MqM4r;heNdcZZeuFicGg zGa&*@cvVe<9lus(fu2fCZFk!|7%h$Eo+ir-39Sd;V)f}(87Hu;4477XvZ1yIys)^b zy01+3pVg(kO}$YTdwRT*XbaR z)D60EBBUs_o<%OK#xi6$x>@z?RjOQUxwcOq zY#g%qoXE<^#FjVJgIGN zV!drm{My+n%@D+YBu^CUon7t7ZVHz1oXIP!J(iJKwCYp#J5QC|ww4wC+iE@kc8T$4 zj?8xphtjruCjw8qv_Ct~MM@^SwlXnnK@G2gfJ5otq9!e-2{gRPCc8ZMWOTs~>u1m1 zji@!IKY3k}KG7q)ET0RhuTiV^=r3RQ3_n|5t%%dMnbl@a7Y__fbNGgMOOFrNXDLEm z56(?%*o#NfW;L7nx+G#2n{uHO#v>t7;Zf&r@M*)RUG(T8@jsTsKCnj}1Jk6A^<}U! zQgrkSb}1ya2ysB2Zo=wYUbhq8qErZAdgMjRqu`FEl2_gcUQ~&$t-d5t>mj|HfXI)( z9Q&r5K+Yb4Uq(c2viKnhu?9xDdCW@;nRhRDL$6v z$nKlu)qle_>OS=UhS^zY2XX9o_D&u|0j*8a?553Ah$7XS0L}41nG4>gGe^VnE$u0J zUV@UBE#6$+$2nA2oT5g;_Ho!YL{gtCvDI5=m(HfxNr>FRqky zpseqBtQ=|l4P*qTLPbjzc{bu91uk>UGDr~I_y zzC{nTLM7={Y~v2{wYvv;?7*R}<`3vE=eZRc-3C-(*N`-v998ymVQ_8hswlAy_vHuX z{n@Vfxa1t~<5tHUE)Lv67gopRlX3h5mq*b5?`AI~MWC|&qerroA`zqHwBk4S02R(m zum3v&elIkKP+O$q^p?S8pe5!j(Hutt5f0r7|0K%`zY7#asQ?d&{#w z?7el6wp+LJVU#6$vq-=zkgDYq+eVw;4rSw{N5UcvekC5i2j(#^L&=Tr`=aBzBIu$DCK$t5nXHN2JoM?7GkG)N%c&PxENfk8vY z&rF0{`5!X6cqfrAew#!(xH$>tPp(KDaU}_`pPZ64-cBYRxswcv)*shO0pziEDR9L{ z{?;M|vZZ`CGNeE|B^a5>gufZdfcxl|Ot>1C1pU>I)DI=Mk1i%aRNnULxrOn>ABK4F zWM}T}m3T&2Z;#Q|Qy;Ldvm)vGyl~QwqEL-j0_pkd>Or}LtNqvU_dDBsv#i9Ebl_5T zjW;922PvyB_Qi@1LppxHeQdH=dU(-qZskC{J33w9_RecWI@7#~zl=y#JGvf8} z8KT}Wsfgjt#Vln{&NmmGUe`Mu8rMP$k@PTNRtq83{PoIvEd=hDPfMXtnG?`bNS9E1 z(+&AK)o)=f#29Po5MrAP1=mjI4APOeO-1`ffpf>~Y9i3DeAQT-RLja~C^kws{yHdK z^j5NIl_pTPH-eZZ;jwXQszAQ|K1qaO!^De$iZ=Si2#N01jiyD3#h-thT(PX4Al|Ai ztR#G9ey*&j?p8fY007K>Sz30vX>)=g>F= zlpOA^31cM#&ncej84s$BmsO3Y`|gMX-0!Y& zn(vZ2a9t9}gjlWrW8m^7^%D(JjduKD7Y#_a>GjFi42zWc5k!AA5gZP%s~8P zp!NHEz=jD*_}s|HS{&i5-2MiStU%FukLE?pLBham_5-$TsrtFi(a%jU;;wS}{i@oC z^2ySXcsPqQ(pK*mIzj*U5ShbeVQ+g1JIQXzmEHpDmgYk~;qo{k1(^iY9n#zRydX+EsVhn&Y0@ z^LUB^NoUha^`he9+w&)_^`iI58)5yrxQkyim66j4FZD$sZeq2xd)wE!2;y5uf%bNT zy-?kKBg{_B8uF}S(@0w(QM~76p|xnJ{!zI$!$_~5YJb1fd_g`_W$?)< zWVfWN`=hvsDL0c|Uc^uzJ`^(4|JwzOo}++07^w_c0h_IC(v(8h1C4kQ8!G8?-`B+q zdE)p7hCJ1=fGJniYwRmF9u4>tmgT4f<__zo%C$94Bq;s z9;hRe&U5AF$2(u*nB%5raq`g*xKe=`nK@dW@;Kq*1)sU|@r#AKpL5igk{tf8Dr;!A zrq>s^atqI9pW`U+z8CoeDHmODoZ^aB+TNMPrp2qEd5)oN%DKwue%BcC^2$us5g+gh z`(2go&0+S;Pch>Tj`?ZbX-;*z#0h`rImXnhXLtt;K*u$`pZeeTGW02jGa31y+YEJh z+9O8&m80#?V}`hR_=sUnO?}9SAF?&yca)W?b~v?@2hANLC3+{n>bqmg-LRR-Yl0PE zmoqmsq8l0dq{suzf`C8HR=!h)?l`NDKPQ=jlm7jAb|Fe3>IU1R;@91K?BgaSd^yZz z!6xZ6o-pz~4_U8TiotGuS&t44-Xli%dBsi`YDlRqWaNuV*~U)l&J}EIYt`8m+I>nG z;ktm0Szwfs@rq4RFzl4asIE^L=C4h8Y-))@#?zr$QhJ3en2SeLyT*J>TQNp*GQ@u}f%x|Nz`8No(CmjFMzvie{Q{QuqCph|V zQNoeFqDwjDkP@y(Ls7yTVFQZog-3CI#77)+xVX#Fhupfw zF&|&d=C}=KA8^b8ELZEVr@Z{lE0GYK$2I(Lobxn(_eICbA{77%IO#+ICp_hIh64VU zQ*M6633t!94vF;<*La>&|5{VND7OYbX}2C;wbuNfHTk)0Yw}h7Y~YCO#)+}E;D!Om z7W%5Pm)ZjQ>=Zk|{CL+6l#Nc_YEQX_IcWaO0m!FU1>~6yj(|Sqs3YK9#>$Cu$DdB1 z;@@wRodNAg3m0gH2G<3U2lrP4<(@KLx&qQ%51p3N-N~1()&sSV2NVWK0zK;q$oI#+ zAVnGBKfM9VIn@BTokxx2KhsQj87X zLcO+54>+urv&ajOE{7XvZc%Ck%pH*?$WwOA)l9VB&kV$?W~$#1n^mxdY1A0UcnhPSip3 zF3JJ8fvY0!P&fMi&aUuW+WEyjssZZD02d&fJ5!D>cBH)9q$=g>`Sw75*A5VlBWwZl zScolP-a2Fh=vS?5z@P%hPu76*?Bmve{^K7Tz;VD{TR@-C-VTVb>}WpLz@Fw)6YVLN zo9sb$!Sb(G*+CT)0y~&}Pz|Q+XujIR9#AgVI?z0+xuZ7ka)Qa~gA-f;bAyX3&C@oz z0p{4=IzU`+&}n_HgN6MRBkBNsdxkrpKZ^RVo*vqMc>(4HnOjCn_OK7|0o(?Y8I()0 zDa3jmy#y!Hp2`#Jhk?)~ITe0tvuI9C=X zn)96Sq4S4-jUZF0WNd@RH{Mk5Ctg5)*$d{{E1~rS{FuZg;_r||9`uidPZU0R7fB zS4dLvq^FMNSEt;m{@u0wp6Lmg_bR-A{De2qoX-IFRD3vP0M!55M!;B>VI=eB;_yu)q(Rs$k{d2k5Kk`fAT*rsr;G2F%y} zOmt2ld0utD?Y{=?dA*^p5_qYjuDK&JC2?yxB|#IJW($qxrgro^UYO1Ju}M;Ym-RoZ|&JH#p!; z_?c}W-socllZvB(K7jWAv=5z2H8%nJtLY{nKW_p&KN@Bt{heo`>(`hF_e)H4?l;#& zdNtStuX z9OtEP^?+^CMo;~;*&Q&K-qrygo0@dwJFMLR=M;;p(H!ZbGoTH~ccgufsj7XiJ=9c! z@v9x3!*sI+oNFfAXy-iEfOElgE6tBpA%5Idh4f&K72te&uQko1KiB};g$O%seq|54 zRNS8JKy{ROCm-TO^WNnyKzYHH`gf6=#=APww_F|VXSzG_;xTtXJDRB1&K31EUwft3 zp5KFfo}Y*2FZ7`7)bcDn&F>cL>AoR)?Htmb&TFjQ3EwMqa7o2HKOLP%U37(~lKrBG zSA$Xbpf1o*5#K>4?fkha;5Z}1fzG$2Kfw9yWot;m1F0(7IzKTwxBJMDzfY7g)a{EU zOwqsFi(aifZG}ltAL>q=LD9OJJz>XZ0t}7ejLYs@wjm zQUdhxyEJ_9<;eR-p46~Oj(>Nn>KcC7<=NEZ%e~m|o3-oDUiRP{yT7T?@}F3It^UG) z=iQ3IP4ZadVzA1~@5UYo#^6VlBZ^b5y76qJu;&1847wNTHJ4m+&MY#*<#mNn zB=yhYlOcR(8@`u(D2(moKbFU1?>P%!3Uy1s6|ZzQ&wQJJ1e+5HD9b;zC;_X;KTJd_ z-=Rd@(&D9WT~n&!>%Hzw9PO{B;_6-vhYuDN?NANB{j{SiDYhE^R8VLY6u_i}Dm?^^>8$wa^HDNo1@;FK3@)M&Bz`6**m-RxNGD#N_sQRh7F<@hzM zT{yPiTP<{Kv~Px}w2t{FlEl#sj?X}NT2|iR{xa~3$0RZgq;|N^#Mu`{i54_5(0t3+ zxvNbS=(&Ga9}^o$I{3)IkN4edJm_N`@1;a*C}Gv^6AkbC@vz8 ze9FN`H9Yph-)W=Qda;g3U)S9pl<*B&9@XD7V)3H<+Spk9I7JpF1|_T1tZ-w0daxTG zrvkWfx9H}+o){JXGclOGM-o>IUWXNzcL!ond`q*Us%{L$!@-N~U5pzS%6AvHb>XbJ z;tTq?P_~PEtxKJ|wj1q?Z5v}ox>eqKy4CrsxbfUjfyf11um8C4;s2BdGq$;Kzf40{ zx5D2BH#Y4kI{mO4>yDlH=1Nh&4QJ(zZltu&bKxg%mv?;SahHn!0T<)fh08>{kKgJ- z@$R1=^l_u&RX69+g)a-({pxZhZgThHM z2$2y`EfGbML{|#L=e5eA$JICEH#}ScD&LhLqkrh9N(k(3^<~p)mC#}2gF;#=V4jqz z+La)(eRsDzDna460u-+a{gRLUp&aD9S;^&~_`3o~CTdjhIZ8of{luj!MNsS=DjsnY z1`)suK#=l}_mwgOA_!x zAYE#86I4G7ohJS+gJ#tvlaxcZ^)Fs}bygYcb}76jWWzhO1im}?sA%yLP<>tu^Vi6R zD1x_V?+}Vn2$If2`T4A;5s2+a-vGx2I$Qpd1H{uAnV`wE&Q52Xr9%=4Q97u8Oy}pO zv#&^p%~$0CGl0HkX$CYBYng%9!wN;y`v%m<|xTH5%7~y z$OZBN?pz>#AC?O>2?!BLqy~+^Lqxa{=ube*WgKP0veWHyJ3N#I)P}25V2F^h&R?Z* zK69Y<`&mzHzcUX;eE;l|eHsy2>WW!!b0|mKoClO=>#`Sw^!%xpen0F#T z>PBFaQx*o3@77o4wR3*!oK?&JW2*G1(r#fX^;%PFf3Afnelv{^L-XH;WBzxoQ3q! zGSCYkyJ1|rEByUep_P=|t;ry>F!TLns7{7G8Law^MN3};@*$BFpvN~$g=Z-dsX%(V zFqP|>3d+}~f@GC1n_dUyKdwVxYJ@as`|D=mv zLGc8ckZaukDICXjAb4y}gHE!2dZhx<)ZAp|)vMfx$*cp(T>mR@)sltz2Nr!P|12cf zzHccXInlz$B%zJ4kovBlg-;2W|2SI8?=-ToOt8vpVN3anauZeGn<%qW_?wADzgtY! zXA?yeUl{d{iPZn2OdM3FOWyO3naZ!UGm-F^XrhYCz(bN&uNZi6oN&r}45Tq=g`xcV z3Sr8%TL*ZJ-<{7RMQ=^2omKtd1f$5uvX9#uA;+qob;y&gaqjJUS{K*0EBs z<{6ic*)ly@8cJ3A_&m}Eo? zsXsgUkmz=%5BCoj1bD)SB%gKtsQi#0mG1d@t_vXbY2N@+8!rnmjslFcASyl#s(E%@ zP|YE!L98a+?YALY`N+F%Bz^8_D<3$_M$+F&Hv4cJ$=@s&bDSJww%N)bAF%PF93M{E zYRn184Hj(&pjLGwqi9+0_p4%nxDH9nn z@j+?en^?>b7Os_Jz(xz@RCKt=;y(9retgVdJ|uY>=-`ulTCD(zLw$bg!vNzTz;OrF zeELTar5dlOZ6oRT02}AZF?JsB0S?dO}PJOuEe1& zg}!Qt`!j@-WWbX{{J9WHb=q-0gaK(Z^f0!U0BjTHehhOx!;JGV(w=I6nEOA>yc)(@ z1BZ9I9E?Y@i}mBNjRa7;c<%dnl(e#_M?B{{9{Uh+#-qHgS5iExxg#F$700gA#kpVQ zAs)%5_~W@R!g#E>@3oZ|!gy5Ll3iiOPZ()lc|FX}4fA=z++QIi`*ARYUkewocVUS8 zB7~%iEkc+l$KjGV)`vJGTR0<*^B2ec5_h}5KNs$7@M0SW$o1P08(YYMw4sfcgvt** z8RYs0@h`dV%M5To29W%(5#V$CxgY(^i#}BT%7=1`(fF3d{BL0&xoGHW@$*caF;E!f z&m4YLcru-+P|{G)N3w-rz5NVa{P&Ux37{9QokTNk-8(3oE|B%S|8LlN!8 z?$b2x1C4cDP-VtQ$IJ$T6!@XZ>(=UilM)^K3Am87N&(O1^QMN7~7;w5#Jc z)Er?VjbWQi{vMNYY9iVH1dHRh`1>p*dpgqM{iMZu?dYG1*J8X{9Jht1{I9>;z4QUWD*bf-2H@ezO?$nXyxH&qGmIWEEGv0Ni zc8b?|U#fAOj$c|J(wsKg!B-TADjtpF@v6P?X)o)jm;IcV{fd|C?PcEbvfg-k-{)m~ zc-en=xorp82bF*J;Q5+UhPC_)%k>-qgUTi`O8z&D+SDpN*I(d=i z!R}sE^O=`*z>73@JniLuqL+Eki!achmv!08JmqEG_p%Rjt~6%8?L}!1s!jKD|9esG z)hzZ^&O8zA)2Odd{={ToXQG-H9Xu<3G&ydQ^W)%J^{K+^+*TsiwzvScl{CuCE^X^C53zhk~4+4zO0OKvdIuzhM2YHSOGXDp;Z-ZQ4o99@Y z{jJTu&*r+>?AL7ON1M-MvmdhYEhz|MZlS%^6dT8n`Xpz{ROfZ{eZjsTw|>J$TC=}w zBkh|;+Kdw$X`j>E=HHiW=2M$}wXOE&VVm>j#1Wk1s5+P$GpHSS4DO1J_SCXMJlzaDZ_CDvHb6xk` zltcq%O8C}#*IIWy56|~}_mBHNd!Ie5^$x$^?|t|APJ14je$HZ-Qety%{A96z*I#hk zWn1j{;U=FC&9vBw2JM{%85WDId9~r)5sRIkw*16~G-G!bzyAI2$1JulU(sUt&SGA> z{)Ky_SZvkc9lws=X|ZizpI`p$CX4Oaa@U!;>nt|A@k>A1t1UJ+Wq5&WD=oGrH^1`( zi*4Yq?f<^TCSR}rYU(VD%^J`rqs26f&F?a%;*!Z0<9C&hPq0`LU&WeNEjC8Hc%;Q9 zm900p*-(ph*z&=^SpzKA`O%fLxAwDG@A-!=J=WXU0h7}kPUvB=x?%wxE!OzLfFJs_ zF}79WjbfWyTI~6Q>pp4E_;a57`0dI~ELOQozBT`;Z?P)h7e2DLp2eR0^3nSW*0$Kw z^k=d96Gc*bd+aAvymMcA z|JrK5Xg^hmUz;@SS>sRb=Vn7k{W<+Pi@kE@`uYo{F!5I^?e%~b?>$9>k!*seI!pgWA%G(xBC|J zh3zuy(Z|~Uzw2**kKTt~hqkxxPoKa2SYM~#_d0i1jgwIfx^lYTsvrug$74~5|6LD8 z5qv8OCW?YB1d>EiiMyleyHQlX^$16fCNRs`v(SGxSTHiE_5@8O}1V6TAx!4Z(3n?^u!EgeDf)2qLzdo=)+w;Mp9O50BZ zY?#`KBS?Ve(Leug{XziRNBsbl7exXH*iW4GK{&m~he*dsbA7xtuERqmhWY5fD&P_K zu`MAI)<<~h&t*URC9L)t1}~kV)(1vWr*%IFp!Tn=#fu1>**+bb~jgIQA`kJx+#ix z3JK*Dh14BUG;jIZwfoMRdZlpmiGUY5N9qU&uOk66A|mjC%G&@^Zv;sBR}h{L!D7>v zoBESF$_L4-`+W%Sym>K<1(QUjhtZw_8ph7LR}a0sA&fNf`ek7ZA;*M4>A5b9Gjku~ zsT)QEO0Y02T39%zzD3oi4EJ$Tgk^`X{1SlN(k_7HzW@!Z2?U5hM76`b2uN?g)12l^p$JWNj${8 zQ~lOSJVg%41eHa>Qz{QmLF$4m^dbPw!YuB+DWkJ6Lu~hxEamSkOciC8KcU{60Ks+F z1h{I7m6@1;$b!pD5|Fw(fe%L(op#`OLe(L^ZiBcnf$_wEH4omnCsq<>dNQI0N%=HQu-v4>_va7Ltk^$3P#Ce3;8bfEA(09R|hZ#E6U!HbtD z$3@anzQQ+ytK|M#hL8H(Z2He{=kBZrxsEWXZ0Cl-Nskkb!R;4j)Z-KuMj;YRVO7t@ zI^~4`593ooL$8G~SES_>!*~Bo8n$oDLQ`M4+lf2F%8y~JZoccU8{u%y`Y8SXF+4%K ztA-D$>rK7iO@-?#9P(9O`U>BD<#%7qS(aMzBE0$wtOf z8aIfb8Q=YZXCjdOp}C%e9$!KP)OVCOw!3)zeG$}taqn~Y;N-qQo?*qI_jZ&a=K|ocT{B!`4a}@)HQzrkke|r+)aQr>KibvGP z3PGbk`-pQEx^mh_ekvPttmC7tzE3~rITC_CMD~=*8{lI>ecu z%HXh$udKrW4T!)4P`fbJ=soR+t^}9?h^|)zs6qn99BZ8@@gcy?fusa-?4gpPK+W&^ zerbF=dTmkwzko>a03<`knd5H+JSGK7uLRif=EbcYhX)Xt2KS_T(*nRk)m8_ zxuE&MAy?^6SH;a$yy_}HyUI_l!VwoOsWEktmxP0>^4LL15;%_1OAZL$4mjva10+Y` zwF53%IDj~c&mG7)lLK;Fdk56kKI>p2H3kl74pGJd&8=?PXiEbp8?$)#yzfUFqI!Y+T7NI%gY$iSW~u4o$;EvHhiK z`01h;Jfz_P={SXdX(&qrf>ebwsSv1Cb*8aY#@$j?{8Ck37%TH56`Mrc(;*ePBvGmQ zG|!bfB~9^d8swZW4b;aMut9xkV_Wf$t@y=Oyz40baKN3#8x;qe#p@n&Rb9DC54j4T zJ=Je`N~eTDGT?dB&)iD_!&mY1LGpUHkHw-rz3%7av99yRqCNOBfShjxC_@Be_=03+ z^9XKj7YP?p`aJ>;<~*Uy?}JECMwDKQsQ$$Cx3WJBFhaD!I|5MKI?v#PuB&1+xCtQV zL#Dr#de4VomHnNo;^D%g zL6PZqCEmJD%R;%T@g6g2sO}UB%a~;!9WYr;8e*-OlGi&Iz3yypsCcQ9SLaKFa~g_bdkz zjvUmaLNoQFkM(u>evNw?2j%?4QM~S`zRAVW4~2tVHD_@_GN^*9_}x``>MHy&_^bV- z=N&qZlBZqJe0Zn9b)6S)5b&EhhO8S;=_wE6=>iXO&g>~2?kSz-DShUF=4o|2HDB@+ z9-8ZPzDoV)p&cc-2bw1r^T3;6uDCtW-1(Zp9jQ0X`lVho`6YSZ%yD$Qwg04kGV##7 zGl~G#QT?U^Q5C7jf{tNRJz~>gx@^$UfOty52^X2-j*` z^%1th0UI>#8e^+I!SK!s5sq%Q!c!ZaX)tNyJwZnuZPi~HemXAN)%NB(`CY|U^EC6j z?2~QKxM75?`Vt#zkfQn<8`P%Gw9$_U*;exaTaJ5{+o10)pWCWWGW(UdZ}?Ro%lYpA zwzUX&@4R2ztK#+F*bvj43Msal8=Cj(^=M0Fw$sKZVvMxL#vCClx0`s-_+zWNU++V& zL$7NGjjc?4_lmL66&ugcpvJ)rAuEbGVANPTpgDeB2h_&5aIi~^L0)hmzneMWB)zNu zWCtYA7CX2i(rB%t<|U4r8yK9BeS!n}4)eQ{gG+L*;ix{tRegpF+3>iK-#?7izw7+3 zM-8_Fl3NK!%@rKT?~0C^J35N*9ME^C^$zH}Rk8#6&M@17{2t_>wkWsPj6Qo;v@!i0 zQ2W@;*s?V6aX{nRE{>W%n)7n5V$NT5#RnvFKb8B74oF=MGI~Y7SFc0sMaj>OnoF7M zl117NG0)L;uH&WkrJS2N5MM9X)pkJR#0pN1UMfieaMV2QPMxXqS?8mix0t%pdQ;=T zx$HWn5A?E;kA$zOb2$gG)tui}a~fOC$!*a0m`XM>S_r{ye%I%#(1%zXG_I&+ewXuM z8#MlGYU3?RAcI31C-gcsF6ez|{jBRw*ORUv`F+FGkDSlhxGcspmrOnV7f>4qKZnnF ArvLx| literal 0 HcmV?d00001 diff --git a/Figures/Data/Parameter_N2.mat b/Figures/Data/Parameter_N2.mat new file mode 100644 index 0000000000000000000000000000000000000000..bce2266d6aaf3e16aac9cb9796e8d2cc2aca858c GIT binary patch literal 352 zcmeZu4DoSvQZUssQ1EpO(M`+DN!3vZ$Vn_o%P-2cQgHY2i*PhE(NSyj5-`93qo*%FknI4(6>}aZCnR_*C_FiagUI(>zr@WYq<<89v@LRygp8VX@+}4u?!pb_Tou z0&_s-IK$29Kr+Xe8Ej6WLsDQ;-k##RK+h~zkU8F9a}=ID__)CJ^53408+k zc7cp?g&W0zY!oBNsKzU5Vye%MwWp{&Q+YN`%-EoS=_pG>MSINa*9T-jzu8+cW6SNM NVYiqWlrHe(0RW)bZh!y) literal 0 HcmV?d00001 diff --git a/Figures/Data/Parameter_N3.mat b/Figures/Data/Parameter_N3.mat new file mode 100644 index 0000000000000000000000000000000000000000..0371740d18a355d01b110b53535117d26735eb5d GIT binary patch literal 337 zcmeZu4DoSvQZUssQ1EpO(M`+DN!3vZ$Vn_o%P-2cQgHY2i*PhE(NSyj5-`93qo*%FknI4(6>}aZCnR_*C_FiagUI(>zr@WYq<<89v@LRygp8VX@+}4u?!pb_Tou z0&_s-n83{`Kr+Xe8Ej6WLsDQ;-k##RK+h~z4_V;=R%u~jW`@!Ytf3&iu5i5^$a)z; zdK<5(iK#w2)}EsBOy$`$F=K-QrlTwkH;PVoZmRv9ezW?cLFV@EwOiO3R9$!$0{{y( BW1;{6 literal 0 HcmV?d00001 diff --git a/Figures/Data_ERP_N3.m b/Figures/Data_ERP_N3.m new file mode 100644 index 0000000..f13e088 --- /dev/null +++ b/Figures/Data_ERP_N3.m @@ -0,0 +1,71 @@ +function Data_ERP_N3() +% This function creates the model data depicted in Figure 7 of +% +% A thalamocortical neural mass model of the EEG during NREM sleep and its +% response to auditory stimulation. +% M Schellenberger Costa, A Weigenand, H-VV Ngo, L Marshall, J Born, T Martinetz, +% JC Claussen. +% PLoS Computational Biology (in review). + +% To ensure availability of the simulation routine and the utility functions +% from fieldtrip it should be called from Create_Data.m + +% Simulation parameters +load('Data/Parameter_N3'); + +% Simulation duration +T = 3600; % duration of the simulation + +% stimulation parameters +% first number is the mode of stimulation +% 0 == none +% 1 == semi-periodic +% 2 == phase dependend + +var_stim = [ 2; % mode of stimulation + 70; % strength of the stimulus in Hz (spikes per second) + 80; % duration of the stimulus in ms + 5; % time between stimulation events in s (ISI) + 0; % range of ISI in s [ISI-range,ISI+range] + 2; % Number of stimuli per event + 1050; % time between stimuli within a event in ms + 450]; % time until stimuli after minimum in ms + + +% Model Output is given as: +% 1. Pyramidal membrane voltage in mV +% 2. Thalamic relay membrane voltage in mV +% 3. Thalamic calcium concentration in mu mol +% 4. Thalamic I_h activation unitless +% 5. Stimulation markers in sampling rate +[Vp, Vt, Ca, ah, Marker_Stim] = TC_mex(T, Param_Cortex, Param_Thalamus, Connectivity, var_stim); + +% Power of spindle bands (BP filtered) +Fs = length(Vp)/T; +Vp_FSP = ft_preproc_hilbert(ft_preproc_bandpassfilter(Vp, Fs, [12, 15], 513, 'fir'), 'abs').^2; +Vp_SSP = ft_preproc_hilbert(ft_preproc_bandpassfilter(Vp, Fs, [09, 12], 513, 'fir'), 'abs').^2; + +% Check whether stimuli are too early +Range_ERP = [-1, 3]; + +Marker_Stim = Marker_Stim(Marker_Stim>-Range_ERP(1)*Fs); +Marker_Stim = Marker_Stim(Marker_Stim< (T-Range_ERP(2))*Fs); + +% Define the matrices +N_Stim = length(Marker_Stim); +time_event = linspace(Range_ERP(1), Range_ERP(2), (Range_ERP(2)-Range_ERP(1))*Fs+1); +Events = zeros(length(time_event), N_Stim); +Events_T = zeros(length(time_event), N_Stim); +Events_FSP = zeros(length(time_event), N_Stim); +Events_SSP = zeros(length(time_event), N_Stim); + +for i=1:N_Stim + Events(:,i) = Vp((Marker_Stim(i)+Range_ERP(1)*Fs)+1:(Marker_Stim(i)+Range_ERP(2)*Fs+1)); + Events_T(:,i) = Vt((Marker_Stim(i)+Range_ERP(1)*Fs)+1:(Marker_Stim(i)+Range_ERP(2)*Fs+1)); + Events_FSP(:,i) = Vp_FSP((Marker_Stim(i)+Range_ERP(1)*Fs)+1:(Marker_Stim(i)+Range_ERP(2)*Fs+1)); + Events_SSP(:,i) = Vp_SSP((Marker_Stim(i)+Range_ERP(1)*Fs)+1:(Marker_Stim(i)+Range_ERP(2)*Fs+1)); +end + +save('Data/ERP_Stim_Model', 'Events', 'Events_T', 'Events_FSP', 'Events_SSP', 'Marker_Stim'); +save('Data/Ca_Depletion', 'Vp', 'Vt', 'Ca', 'ah', 'Marker_Stim'); +end \ No newline at end of file diff --git a/Figures/Data_SO_Average.m b/Figures/Data_SO_Average.m new file mode 100644 index 0000000..404dbeb --- /dev/null +++ b/Figures/Data_SO_Average.m @@ -0,0 +1,73 @@ +function Data_SO_Average(type) +% This function creates the model data depicted in Figure 6 of +% +% A thalamocortical neural mass model of the EEG during NREM sleep and its +% response to auditory stimulation. +% M Schellenberger Costa, A Weigenand, H-VV Ngo, L Marshall, J Born, T Martinetz, +% JC Claussen. +% PLoS Computational Biology (in review). + +% To ensure availability of the simulation routine and the utility functions +% from fieldtrip it should be called from Create_Data.m +if nargin==0 + type = 2; +end + +switch type + case 1; + load_model = 'Data/Time_Series_N2.mat'; + load_data = 'Data/KC_Average_data.mat'; + savename = 'Data/KC_Average.mat'; + case 2; + load_model = 'Data/Time_Series_N3.mat'; + load_data = 'Data/SO_Average_data.mat'; + savename = 'Data/SO_Average.mat'; +end + +% Load data +load(load_model, 'Vp'); +load(load_data); + +% Rename data +mean_ERP_data = mean_ERP_sham; +mean_FSP_data = mean_FSP_sham; +sem_ERP_data = sem_ERP_sham; +sem_FSP_data = sem_FSP_sham; + +% Process model Data +T = 3600; +Fs = length(Vp)/T; + +% Power of fast spindle band (BP filtered) and slow wave activity for peak +% detection +Vp_low = ft_preproc_bandpassfilter(Vp, Fs, [0.25,4], 513, 'fir') + mean(Vp); +Vp_FSP = ft_preproc_hilbert(ft_preproc_bandpassfilter(Vp, Fs, [12, 15], 513, 'fir'), 'abs').^2; + +% Search for peaks +[~, x_SO] = findpeaks(-Vp_low, 'MINPEAKHEIGHT', 68, 'MINPEAKDISTANCE', 0.2*Fs); + +% Remove those events, that are too close to begin/end +x_SO = x_SO(x_SO<(T-2)*Fs); +x_SO = x_SO(x_SO> 2*Fs); + +% Set the variables +N_Stim = length(x_SO); +Range_ERP = [-1.25, 1.25]; +Events = zeros(length(time_events), N_Stim); +Events_FSP = zeros(length(time_events), N_Stim); + +% Segmentation +for i=1:N_Stim + Events(:,i) = Vp ((x_SO(i)+Range_ERP(1)*Fs)+1:(x_SO(i)+Range_ERP(2)*Fs+1)); + Events_FSP(:,i) = Vp_FSP((x_SO(i)+Range_ERP(1)*Fs)+1:(x_SO(i)+Range_ERP(2)*Fs+1)); +end + +% Averaging +mean_ERP_model= mean(Events, 2); %#ok<*NASGU> +mean_FSP_model= mean(Events_FSP,2); +sd_ERP_model = std (Events, 0, 2); +sd_FSP_model = std (Events_FSP,0, 2); + +% Save the data +save(savename, 'N_Stim', 'time_events', 'mean_ERP_data', 'mean_FSP_data', 'sem_ERP_data', 'sem_FSP_data', 'mean_ERP_model', 'mean_FSP_model', 'sd_ERP_model', 'sd_FSP_model'); +end \ No newline at end of file diff --git a/Figures/Data_Time_Series.m b/Figures/Data_Time_Series.m new file mode 100644 index 0000000..2a3d9cd --- /dev/null +++ b/Figures/Data_Time_Series.m @@ -0,0 +1,40 @@ +function Data_Time_Series(type) +% This function creates the model data depicted in Figure 4-8 of +% +% A thalamocortical neural mass model of the EEG during NREM sleep and its +% response to auditory stimulation. +% M Schellenberger Costa, A Weigenand, H-VV Ngo, L Marshall, J Born, T Martinetz, +% JC Claussen. +% PLoS Computational Biology (in review). + +% To ensure availability of the simulation routine and the utility functions +% from fieldtrip it should be called from Create_Data.m + +% Load the parameter settings +if type == 1 + load('Data/Parameter_N2'); + Protocol_Name = 'N2'; +else + load('Data/Parameter_N3'); + Protocol_Name = 'N3'; +end + +% there is no stimulation so set var_stim to zero +var_stim = zeros(8,1); + +% Duration of the simulation +T = 3600; + +% Run the simulation (Can take some time) +[Vp, Vt, Ca, ah, ~] = TC_mex(T, Param_Cortex, Param_Thalamus, Connectivity, var_stim); + +% Save the full data +save(['Data/Time_Series_',Protocol_Name], 'Vp', 'Vt', 'Ca', 'ah'); + +% Save a smaller snipplet for example time series plot +Vp = Vp(1:3000); +Vt = Vt(1:3000); +Ca = Ca(1:3000); +ah = ah(1:3000); +save(['Data/Time_Series_Short_',Protocol_Name], 'Vp', 'Vt', 'Ca', 'ah'); +end \ No newline at end of file diff --git a/Figures/Import_Data.m b/Figures/Import_Data.m new file mode 100644 index 0000000..1912da3 --- /dev/null +++ b/Figures/Import_Data.m @@ -0,0 +1,44 @@ +function Import_Data(Type) +% This function imports the experimental data from Ngo et al 2013 utilized +% in +% +% A thalamocortical neural mass model of the EEG during NREM sleep and its +% response to auditory stimulation. +% M Schellenberger Costa, A Weigenand, H-VV Ngo, L Marshall, J Born, T Martinetz, +% JC Claussen. +% PLoS Computational Biology (in review). + +if nargin ==0 + Type = 1; +end + +load('Data/Orig/Experimental_Data'); + +switch (Type) + case 1; + Data = data.KC_Average; + sn = 'Data/KC_Average_data'; + case 2; + Data = data.SO_Average; + sn = 'Data/SO_Average_data'; + case 3; + Data = data.ERP_Average; + sn = 'Data/ERP_Average_data'; +end + +time_events = Data(:,1); %#ok<*NASGU> + +mean_ERP = Data(:,2); +mean_ERP_sham = Data(:,3); + +sem_ERP = Data(:,4); +sem_ERP_sham = Data(:,5); + +mean_FSP = Data(:,6); +mean_FSP_sham = Data(:,7); + +sem_FSP = Data(:,8); +sem_FSP_sham = Data(:,9); + +save(sn, 'time_events', 'mean_ERP', 'mean_ERP_sham', 'sem_ERP', 'sem_ERP_sham', 'mean_FSP', 'mean_FSP_sham', 'sem_FSP', 'sem_FSP_sham'); +end diff --git a/Figures/Pictogramm_TC.tex b/Figures/Pictogramm_TC.tex new file mode 100644 index 0000000..804772c --- /dev/null +++ b/Figures/Pictogramm_TC.tex @@ -0,0 +1,64 @@ +% command to run this on the terminal is lualatex -shell-escape Pictogram_Sleep_Regulation.tex +% requires pgfplots + +\pdfpxdimen=1in +\divide\pdfpxdimen by 600 + +\documentclass[10pt]{standalone} +\usepackage{tikz} +%has to be compiled with ``lualatex -shell-escape Pictogram_TC.tex'' +\usetikzlibrary{external, positioning, backgrounds, arrows.meta,decorations.pathmorphing} + +% set up externalization +%\tikzexternalize[prefix=tikz/] +% WARNING the name of the Image CANNOT be the same as the tex file +\tikzsetnextfilename{Pictogram_TC} +\tikzset{external/system call={lualatex \tikzexternalcheckshellescape -halt-on-error +-interaction=batchmode -jobname "\image" "\texsource"; +convert -density 600 -compress LZW "\image".pdf -flatten -alpha off -colorspace RGB -depth 8 "\image".tif; +convert -density 600 "\image".tif "\image".jpg;}} +\tikzexternalize + +\begin{document} +\begin{tikzpicture}[ + s/.style = {shorten >=1pt}, + every loop/.style={shorten >=1pt}, + pop/.style = {rectangle, draw, text width=10em, text centered, minimum height=3em, rounded corners=3mm}, + noise_C/.style = {decorate, decoration={zigzag}, s,-{Circle[fill=blue!30]}}, + noise_T/.style = {decorate, decoration={zigzag}, s,-{Circle[fill=green!30]}}, + Ex/.style = {rounded corners, s, -{Circle[fill=blue!30]}}, + Tc/.style = {rounded corners, s, -{Circle[fill=green!30]}}, + In/.style = {rounded corners, s, -|}, + ] + + % Cortex populations + \node (ex) at (0,0) [pop, fill=blue!30] {Pyramidal (p)}; + \node (in) at (6,0) [pop, fill=yellow!30] {Inhibitory (i)}; + % Thalamus Populations + \node (tc) at (6,-2) [pop, fill=green!30] {Thalamocortical (t)}; + \node (re) at (0,-2) [pop, fill=yellow!30] {Thalamic Reticular (r)}; + + % Connections within cortex + \draw + ( ex) edge [Ex, loop above] node[below left, xshift = -0.4em] {$N_{pp}$} ( ex) + ([yshift= 0.2cm] ex.east) edge [Ex] node[above] {$N_{ip}$} ([yshift= 0.2cm]in.west) + ([yshift=-0.2cm] in.west) edge [In] node[below] {$N_{pi}$} ([yshift=-0.2cm]ex.east) + ( in) edge [In, loop above] node[below left, xshift = -0.4em] {$N_{ii}$} ( in) + (0.5,1.22) edge [noise_C] node[right] {$\phi_n$} ([xshift= 0.5cm]ex.north) + (6.5,1.22) edge [noise_C] node[right] {$\phi^{'}_n$} ([xshift= 0.5cm]in.north); + + % Connections within thalamus + \draw + (6.5,-3.22) edge [noise_T] node[right] {$\phi^{''}_{n}$}([xshift=0.5cm] tc.south) + ([yshift= 0.15cm]tc.west) edge [Tc] node[above] {$N_{rt}$} ([yshift= 0.15cm] re.east) + ([yshift=-0.15cm]re.east) edge [In] node[below] {$N_{tr}$} ([yshift=-0.15cm] tc.west) + ( re) edge [In, loop below] node[above left, xshift = -0.4em] {$N_{rr}$} ( re); + + % TC connections + \path + (ex.south) edge [Ex] node[left] {$N_{rp}$} ([yshift= 0.02cm] re.north) + (ex.south east) edge [Ex, bend right=15, yshift =-0.2em] node[at start,below]{$N_{tp}$} (tc.north west) + (tc.north) edge [Tc] node[right] {$N_{it}$} ([yshift=-0.02cm] in.south) + (tc.north west) edge [Tc, bend right=15, yshift = 0.2em] node[at start,above]{$N_{pt}$} (ex.south east); +\end{tikzpicture} +\end{document} \ No newline at end of file diff --git a/Figures/Plot_Calcium_Depletion.m b/Figures/Plot_Calcium_Depletion.m new file mode 100644 index 0000000..b5ed50c --- /dev/null +++ b/Figures/Plot_Calcium_Depletion.m @@ -0,0 +1,86 @@ +function Plot_Calcium_Depletion() +% This function creates Figure 8 of +% +% A thalamocortical neural mass model of the EEG during NREM sleep and its +% response to auditory stimulation. +% M Schellenberger Costa, A Weigenand, H-VV Ngo, L Marshall, J Born, T Martinetz, +% JC Claussen. +% PLoS Computational Biology (in review). + +% Path to utility functions +if(isempty(strfind(path, [pwd, '/Tools']))) + addpath([pwd, '/Tools']); +end + +% Load the endogenous data +load('Data/Time_Series_N3.mat'); +Vp_endo = Vp(10900:11400); +Vt_endo = Vt(10900:11400); +ah_endo = ah(10900:11400); +timeaxis = linspace(0,5,501); + +Range = [-1.25, 3.75]; +Fs = 100; +i = 5; + +% Load the stimulationd data +load('Data/Ca_Depletion.mat'); +Vp_stim = Vp((Marker_Stim(i)+Range(1)*Fs)+1:(Marker_Stim(i)+Range(2)*Fs+1)); +Vt_stim = Vt((Marker_Stim(i)+Range(1)*Fs)+1:(Marker_Stim(i)+Range(2)*Fs+1)); +ah_stim = ah((Marker_Stim(i)+Range(1)*Fs)+1:(Marker_Stim(i)+Range(2)*Fs+1)); + +Lines = [1.26, 2.31; 1.26, 2.31]; + +% Create figure +figure(2) +clf + +% Create panel +p = panel('no-manage-font'); +p.pack(3,2); + +% Set margins +p.de.marginbottom = 10; +p.de.margintop = 5; +p.de.marginleft = 5; +p.de.marginright = 5; + +% Plot the data +p(1,1).select(); +plot(timeaxis, Vp_endo, 'color', 'black'); +set(gca, 'XTick', 1:4, 'XTickLabel', {'','','',''}, 'YTick', [-80, -60, -40], 'YLim', [-80, -40], 'XLim', [0,5]); +p(1,1).ylabel('V_{p} [mV]'); +p(1,1).title('(A) Endogenous SOs'); + +p(2,1).select(); +plot(timeaxis, Vt_endo, 'color', 'black'); +set(gca, 'XTick', 1:4, 'XTickLabel', {'','','',''}, 'YTick', [-70, -55, -40], 'YLim', [-70, -40], 'XLim', [0,5]); +p(2,1).ylabel('V_{t} [mV]'); + +p(3,1).select(); +plot(timeaxis, ah_endo.*51, 'color', 'black'); +set(gca, 'XTick', 1:4, 'XTickLabel', {'','','',''}, 'YTick', [15, 18, 21], 'YLim', [15, 21], 'XLim', [0,5]); +p(3,1).xlabel('Time [s]'); +p(3,1).ylabel('g_{h} [\muS/cm^{2}]'); + +p(1,2).select(); +plot(timeaxis, Vp_stim, 'color', 'black'); +set(gca, 'XTick', 1:4, 'XTickLabel', {'','','',''}, 'YAxisLocation', 'Right', 'YTick', [-80, -60, -40], 'YLim', [-80, -40], 'XLim', [0,5]); +line(Lines, get(gca, 'YLim'), 'color', 'black', 'LineStyle', ':'); +p(1,2).ylabel('V_{p} [mV]'); +p(1,2).title('(B) Closed loop stimulation'); + +p(2,2).select(); +plot(timeaxis, Vt_stim, 'color', 'black'); +set(gca, 'XTick', 1:4, 'XTickLabel', {'','','',''}, 'YAxisLocation', 'Right', 'YTick', [-70, -55, -40], 'YLim', [-70, -40], 'XLim', [0,5]); +line(Lines, get(gca, 'YLim'), 'color', 'black', 'LineStyle', ':'); +p(2,2).ylabel('V_{t} [mV]'); + +p(3,2).select(); +plot(timeaxis, ah_stim.*51, 'color', 'black'); +set(gca, 'XTick', 1:4, 'XTickLabel', {'','','',''}, 'YAxisLocation', 'Right', 'YTick', [15, 18, 21], 'YLim', [15, 21], 'XLim', [0,5]); +line(Lines, get(gca, 'YLim'), 'color', 'black', 'LineStyle', ':'); +p(3,2).xlabel('Time [s]'); +p(3,2).ylabel('g_{h} [\muS/cm^{2}]'); + +end \ No newline at end of file diff --git a/Figures/Plot_Compare_ERP.m b/Figures/Plot_Compare_ERP.m new file mode 100644 index 0000000..b67a31d --- /dev/null +++ b/Figures/Plot_Compare_ERP.m @@ -0,0 +1,111 @@ +function Plot_Compare_ERP() +% This function creates Figure 7 of +% +% A thalamocortical neural mass model of the EEG during NREM sleep and its +% response to auditory stimulation. +% M Schellenberger Costa, A Weigenand, H-VV Ngo, L Marshall, J Born, T Martinetz, +% JC Claussen. +% PLoS Computational Biology (in review). + +% Path to utility functions +if(isempty(strfind(path, [pwd, '/Tools']))) + addpath([pwd, '/Tools']); +end + +% Load the data +load('Data/ERP_Average_data'); +load('Data/ERP_double_model', 'Events', 'Events_FSP'); + +% Calculate average/sem +mean_ERP_model= mean(Events, 2); +mean_FSP_model= mean(Events_FSP,2); +sd_ERP_model = std (Events, 0, 2); +sd_FSP_model = std (Events_FSP,0, 2); + +mean_ERP_data = mean_ERP; +mean_FSP_data = mean_FSP; +sem_ERP_data = sem_ERP; +sem_FSP_data = sem_FSP; + +% Option array for set +Model_Range_ERP = [-75, -45]; +Data_Range_ERP = [-80, 50]; +Model_Range_FSP = [-0.25, 1.25]; +Data_Range_FSP = [2, 8]; +xRange = -1:0.5:3; + +Option_Name = { 'YLim'; + 'YTick'; + 'YTickLabel'; + 'YColor'; + 'XTick'; + 'XLim'}'; + +Option_Model_ERP = {Model_Range_ERP; + -75:10:-40; + -75:10:-40; + 'black'; + xRange; + [xRange(1),xRange(end)]}'; %#ok<*NBRAK> + +Option_Data_ERP = {Data_Range_ERP; + -80:40:40; + -80:40:40; + 'black'; + xRange; + [xRange(1),xRange(end)]}'; + +Option_Model_FSP = {Model_Range_FSP; + linspace(Model_Range_FSP(1), Model_Range_FSP(2), 4); + linspace(Model_Range_FSP(1), Model_Range_FSP(2), 4); + 'black'; + xRange; + [xRange(1),xRange(end)]}'; + +Option_Data_FSP = {Data_Range_FSP; + linspace(Data_Range_FSP(1), Data_Range_FSP(2), 4); + linspace(Data_Range_FSP(1), Data_Range_FSP(2), 4); + 'black'; + xRange; + [xRange(1),xRange(end)]}'; + + +Lines = [0, 1.05; 0, 1.05]; + +% Define handle for plotting +BL_model =@(y,x) boundedline(y,x(:,1), x(:,2), 'alpha', 'transparency', 0.1, 'r'); +BL_data =@(y,x) boundedline(y,x(:,1), x(:,2), 'alpha', 'transparency', 0.1, 'black'); + +% Create figure +figure(1) +clf + +% Create panel +p = panel('no-manage-font'); +p.pack(2,1); + +% set margins +p.de.margintop = 10; + +% ERP +p(1,1).select(); +[AX1, ~, ~] = plotyy(time_events,[mean_ERP_data, sem_ERP_data],time_events,[mean_ERP_model, sd_ERP_model], BL_data, BL_model); +set(AX1(1),Option_Name, Option_Data_ERP, 'box', 'off'); +set(AX1(2),Option_Name, Option_Model_ERP); +line(Lines, get(AX1(1), 'YLim'), 'Color', 'black', 'LineStyle', ':'); +line(Lines, get(AX1(2), 'YLim'), 'Color', 'black', 'LineStyle', ':'); +ylabel(AX1(1),'EEG [\muV]'); +ylabel(AX1(2),'V_{p} [mV]'); + +% Spindle power +p(2,1).select(); +[AX2, ~, ~] = plotyy(time_events,[mean_FSP_data, sem_FSP_data],time_events,[mean_FSP_model, sd_FSP_model], BL_data, BL_model); +set(AX2(1),Option_Name, Option_Data_FSP); +set(AX2(2),Option_Name, Option_Model_FSP); +line(Lines, get(AX1(1), 'YLim'), 'Color', 'black', 'LineStyle', ':'); +line(Lines, get(AX1(2), 'YLim'), 'Color', 'black', 'LineStyle', ':'); +line(Lines, get(gca, 'YLim'), 'Color', 'black', 'LineStyle', ':'); +ylabel(AX2(1),'Spindle Power [\muV^{2}]'); +ylabel(AX2(2),'Spindle Power [mV^{2}cd ../NM_Th]'); +xlabel('Time [s]'); +end \ No newline at end of file diff --git a/Figures/Plot_Compare_SO_Average.m b/Figures/Plot_Compare_SO_Average.m new file mode 100644 index 0000000..ee20b2d --- /dev/null +++ b/Figures/Plot_Compare_SO_Average.m @@ -0,0 +1,154 @@ +function Plot_Compare_SO_Average() +% This function creates Figure 6 of +% +% A thalamocortical neural mass model of the EEG during NREM sleep and its +% response to auditory stimulation. +% M Schellenberger Costa, A Weigenand, H-VV Ngo, L Marshall, J Born, T Martinetz, +% JC Claussen. +% PLoS Computational Biology (in review). + +% Path to utility functions +if(isempty(strfind(path, [pwd, '/Tools']))) + addpath([pwd, '/Tools']); +end + +% Load the data +load('Data/KC_Average'); +load('Data/SO_Average'); + + +% Image settings +Data_Range_ERP = [-120, 60]; +xRange = [-1.25,1.25]; +xTicks = -1:0.5:1; + +Model_Range_ERP_N2 = [-75, -41]; +Model_Ticks_ERP_N2 = linspace(-75,-45,4); +Model_Range_FSP_N2 = [-0.25, 2]; +Data_Range_FSP_N2 = [-1, 5]; + +Model_Range_ERP_N3 = [-75, -45]; +Model_Ticks_ERP_N3 = linspace(-75,-45,4); +Model_Range_FSP_N3 = [-0.5, 1.]; +Data_Range_FSP_N3 = [2, 5]; + +mean_ERP_data_N2 = mean_ERP_data; +mean_ERP_model_N2 = mean_ERP_model; +mean_FSP_data_N2 = mean_FSP_data; +mean_FSP_model_N2 = mean_FSP_model; +sem_ERP_data_N2 = sem_ERP_data; +sem_FSP_data_N2 = sem_FSP_data; +sd_ERP_model_N2 = sd_ERP_model; +sd_FSP_model_N2 = sd_FSP_model; + +mean_ERP_data_N3 = mean_ERP_data; +mean_ERP_model_N3 = mean_ERP_model; +mean_FSP_data_N3 = mean_FSP_data; +mean_FSP_model_N3 = mean_FSP_model; +sem_ERP_data_N3 = sem_ERP_data; +sem_FSP_data_N3 = sem_FSP_data; +sd_ERP_model_N3 = sd_ERP_model; +sd_FSP_model_N3 = sd_FSP_model; + +% Define handle for plotting +BL_model =@(y,x) boundedline(y,x(:,1), x(:,2), 'alpha', 'transparency', 0.1, 'r'); +BL_data =@(y,x) boundedline(y,x(:,1), x(:,2), 'alpha', 'transparency', 0.1, 'black'); + +% Option array for set +Option_Name = { 'YLim'; + 'YTick'; + 'YColor'; + 'XTick'; + 'XLim'}'; +Option_Data_ERP = {Data_Range_ERP; + linspace(Data_Range_ERP(1), Data_Range_ERP(2), 4); + 'black'; + xTicks; + xRange}'; + +Option_Model_ERP_N2 = { Model_Range_ERP_N2; + Model_Ticks_ERP_N2; + 'black'; + xTicks; + xRange}'; %#ok<*NBRAK> + +Option_Model_FSP_N2 = { Model_Range_FSP_N2; + linspace(Model_Range_FSP_N2(1), Model_Range_FSP_N2(2), 4); + 'black'; + xTicks; + xRange}'; + +Option_Data_FSP_N2 = { Data_Range_FSP_N2; + linspace(Data_Range_FSP_N2(1), Data_Range_FSP_N2(2), 4); + 'black'; + xTicks; + xRange}'; + +Option_Model_ERP_N3 = { Model_Range_ERP_N3; + Model_Ticks_ERP_N3; + 'black'; + xTicks; + xRange}'; %#ok<*NBRAK> + +Option_Model_FSP_N3 = { Model_Range_FSP_N3; + linspace(Model_Range_FSP_N3(1), Model_Range_FSP_N3(2), 4); + 'black'; + xTicks; + xRange}'; + +Option_Data_FSP_N3 = { Data_Range_FSP_N3; + linspace(Data_Range_FSP_N3(1), Data_Range_FSP_N3(2), 4); + 'black'; + xTicks; + xRange}'; + +% Create figure +figure(1); +clf, shg + +% Create panel +p = panel('no-manage-font'); +p.pack(2,2); + +% set margins +p.de.marginbottom = 10; +p.de.margintop = 10; +p.de.marginleft = 25; +p.de.marginright = 25; + +% ERP +p(1,1).select(); +[AX1_N2, ~, ~] = plotyy(time_events,[mean_ERP_data_N2, sem_ERP_data_N2],time_events,[mean_ERP_model_N2, sd_ERP_model_N2], BL_data, BL_model); +set(AX1_N2(1), Option_Name, Option_Data_ERP, 'box', 'off', 'XTickLabel', []); +set(AX1_N2(2), Option_Name, Option_Model_ERP_N2, 'XTickLabel', []); +ylabel(AX1_N2(1),'EEG [\muV]'); +ylabel(AX1_N2(2),'V_{p} [mV]'); +p(1,1).title('(A) N2'); + +% Spindle power +p(2,1).select(); +[AX2_N2, ~, ~] = plotyy(time_events,[mean_FSP_data_N2, sem_FSP_data_N2],time_events,[mean_FSP_model_N2, sd_FSP_model_N2], BL_data, BL_model); +set(AX2_N2(1), Option_Name, Option_Data_FSP_N2); +set(AX2_N2(2), Option_Name, Option_Model_FSP_N2); +xlabel('Time [s]'); +ylabel(AX2_N2(1),'Spindle Power [\muV^{2}]'); +ylabel(AX2_N2(2),'Spindle Power [mV^{2}]'); + +% ERP +p(1,2).select(); +[AX1_N3, ~, ~] = plotyy(time_events,[mean_ERP_data_N3, sem_ERP_data_N3],time_events,[mean_ERP_model_N3, sd_ERP_model_N3], BL_data, BL_model); +set(AX1_N3(1), Option_Name, Option_Data_ERP, 'XTickLabel', []); +set(AX1_N3(2), Option_Name, Option_Model_ERP_N3, 'XTickLabel', []); +ylabel(AX1_N3(1),'EEG [\muV]'); +ylabel(AX1_N3(2),'V_{p} [mV]'); +p(1,2).title('(B) N3'); + +% Spindle power +p(2,2).select(); +[AX2_N3, ~, ~] = plotyy(time_events,[mean_FSP_data_N3, sem_FSP_data_N3],time_events,[mean_FSP_model_N3, sd_FSP_model_N3], BL_data, BL_model); +set(AX2_N3(1), Option_Name, Option_Data_FSP_N3); +set(AX2_N3(2), Option_Name, Option_Model_FSP_N3); +xlabel('Time [s]'); +ylabel(AX2_N3(1),'Spindle Power [\muV^{2}]'); +ylabel(AX2_N3(2),'Spindle Power [mV^{2}]'); +end \ No newline at end of file diff --git a/Figures/Plot_Time_Series.m b/Figures/Plot_Time_Series.m new file mode 100644 index 0000000..04ec6ba --- /dev/null +++ b/Figures/Plot_Time_Series.m @@ -0,0 +1,62 @@ +function Plot_Time_Series(type) +% This function creates Figure 4 and 5 of +% +% A thalamocortical neural mass model of the EEG during NREM sleep and its +% response to auditory stimulation. +% M Schellenberger Costa, A Weigenand, H-VV Ngo, L Marshall, J Born, T Martinetz, +% JC Claussen. +% PLoS Computational Biology (in review). + +if nargin==0 + type = 1; +end + + +% Path to utility functions +if(isempty(strfind(path, [pwd, '/Tools']))) + addpath([pwd, '/Tools']); +end + +if type == 1 + Name = 'N2'; +else + Name = 'N3'; +end + +% load the data +load(['Data/Time_Series_Short_',Name]); + +T = 30; +L = max(size(Vt)); +timeaxis = linspace(0,T,L); + +% Create figure +figure(1) +clf + +% Create panel +p = panel('no-manage-font'); +p.pack(2,1); + +% set margins +p.de.margintop = 10; + +p(1,1).select(); +plot(timeaxis, Vp, 'Color', 'black'); +ylabel(gca, 'V_{p} [mV]'); +set(gca, 'XTick', [], 'YTick', [-80, -60, -40], 'YLim', [-80, -40]); + +% thalamic relay membrane voltage +p(2,1).select(); +plot(timeaxis, Vt, 'Color', 'black'); +p(2,1).xlabel('Time [s]'); +ylabel(gca, 'V_{t} [mV]'); +switch (type) + case 1 + ylim([-70, -30]); + set(gca, 'ytick', [-70, -50, -30]); + case 2 + ylim([-70, -40]); + set(gca, 'ytick', [-70, -55, -40]); +end +end \ No newline at end of file diff --git a/Figures/Tools/boundedline.m b/Figures/Tools/boundedline.m new file mode 100644 index 0000000..a1e75e9 --- /dev/null +++ b/Figures/Tools/boundedline.m @@ -0,0 +1,392 @@ +function varargout = boundedline(varargin) +%BOUNDEDLINE Plot a line with shaded error/confidence bounds +% +% [hl, hp] = boundedline(x, y, b) +% [hl, hp] = boundedline(x, y, b, linespec) +% [hl, hp] = boundedline(x1, y1, b1, linespec1, x2, y2, b2, linespec2) +% [hl, hp] = boundedline(..., 'alpha') +% [hl, hp] = boundedline(..., ax) +% [hl, hp] = boundedline(..., 'transparency', trans) +% [hl, hp] = boundedline(..., 'orientation', orient) +% [hl, hp] = boundedline(..., 'cmap', cmap) +% +% Input variables: +% +% x, y: x and y values, either vectors of the same length, matrices +% of the same size, or vector/matrix pair where the row or +% column size of the array matches the length of the vector +% (same requirements as for plot function). +% +% b: npoint x nside x nline array. Distance from line to +% boundary, for each point along the line (dimension 1), for +% each side of the line (lower/upper or left/right, depending +% on orientation) (dimension 2), and for each plotted line +% described by the preceding x-y values (dimension 3). If +% size(b,1) == 1, the bounds will be the same for all points +% along the line. If size(b,2) == 1, the bounds will be +% symmetrical on both sides of the lines. If size(b,3) == 1, +% the same bounds will be applied to all lines described by +% the preceding x-y arrays (only applicable when either x or +% y is an array). Bounds cannot include Inf, -Inf, or NaN, +% +% linespec: line specification that determines line type, marker +% symbol, and color of the plotted lines for the preceding +% x-y values. +% +% 'alpha': if included, the bounded area will be rendered with a +% partially-transparent patch the same color as the +% corresponding line(s). If not included, the bounded area +% will be an opaque patch with a lighter shade of the +% corresponding line color. +% +% ax: handle of axis where lines will be plotted. If not +% included, the current axis will be used. +% +% transp: Scalar between 0 and 1 indicating with the transparency or +% intensity of color of the bounded area patch. Default is +% 0.2. +% +% orient: 'vert': add bounds in vertical (y) direction (default) +% 'horiz': add bounds in horizontal (x) direction +% +% cmap: n x 3 colormap array. If included, lines will be colored +% (in order of plotting) according to this colormap, +% overriding any linespec or default colors. +% +% Output variables: +% +% hl: handles to line objects +% +% hp: handles to patch objects +% +% Example: +% +% x = linspace(0, 2*pi, 50); +% y1 = sin(x); +% y2 = cos(x); +% e1 = rand(size(y1))*.5+.5; +% e2 = [.25 .5]; +% +% ax(1) = subplot(2,2,1); +% [l,p] = boundedline(x, y1, e1, '-b*', x, y2, e2, '--ro'); +% outlinebounds(l,p); +% title('Opaque bounds, with outline'); +% +% ax(2) = subplot(2,2,2); +% boundedline(x, [y1;y2], rand(length(y1),2,2)*.5+.5, 'alpha'); +% title('Transparent bounds'); +% +% ax(3) = subplot(2,2,3); +% boundedline([y1;y2], x, e1(1), 'orientation', 'horiz') +% title('Horizontal bounds'); +% +% ax(4) = subplot(2,2,4); +% boundedline(x, repmat(y1, 4,1), permute(0.5:-0.1:0.2, [3 1 2]), ... +% 'cmap', cool(4), 'transparency', 0.5); +% title('Multiple bounds using colormap'); + + +% Copyright 2010 Kelly Kearney + +%-------------------- +% Parse input +%-------------------- + +% Alpha flag + +isalpha = cellfun(@(x) ischar(x) && strcmp(x, 'alpha'), varargin); +if any(isalpha) + usealpha = true; + varargin = varargin(~isalpha); +else + usealpha = false; +end + +% Axis + +isax = cellfun(@(x) isscalar(x) && ishandle(x) && strcmp('axes', get(x,'type')), varargin); +if any(isax) + hax = varargin{isax}; + varargin = varargin(~isax); +else + hax = gca; +end + +% Transparency + +[found, trans, varargin] = parseparam(varargin, 'transparency'); + +if ~found + trans = 0.2; +end + +if ~isscalar(trans) || trans < 0 || trans > 1 + error('Transparency must be scalar between 0 and 1'); +end + +% Orientation + +[found, orient, varargin] = parseparam(varargin, 'orientation'); + +if ~found + orient = 'vert'; +end + +if strcmp(orient, 'vert') + isvert = true; +elseif strcmp(orient, 'horiz') + isvert = false; +else + error('Orientation must be ''vert'' or ''horiz'''); +end + + +% Colormap + +[hascmap, cmap, varargin] = parseparam(varargin, 'cmap'); + + +% X, Y, E triplets, and linespec + +[x,y,err,linespec] = deal(cell(0)); +while ~isempty(varargin) + if length(varargin) < 3 + error('Unexpected input: should be x, y, bounds triplets'); + end + if all(cellfun(@isnumeric, varargin(1:3))) + x = [x varargin(1)]; + y = [y varargin(2)]; + err = [err varargin(3)]; + varargin(1:3) = []; + else + error('Unexpected input: should be x, y, bounds triplets'); + end + if ~isempty(varargin) && ischar(varargin{1}) + linespec = [linespec varargin(1)]; + varargin(1) = []; + else + linespec = [linespec {[]}]; + end +end + +%-------------------- +% Reformat x and y +% for line and patch +% plotting +%-------------------- + +% Calculate y values for bounding lines + +plotdata = cell(0,7); + +htemp = figure('visible', 'off'); +for ix = 1:length(x) + + % Get full x, y, and linespec data for each line (easier to let plot + % check for properly-sized x and y and expand values than to try to do + % it myself) + + try + if isempty(linespec{ix}) + hltemp = plot(x{ix}, y{ix}); + else + hltemp = plot(x{ix}, y{ix}, linespec{ix}); + end + catch + close(htemp); + error('X and Y matrices and/or linespec not appropriate for line plot'); + end + + linedata = get(hltemp, {'xdata', 'ydata', 'marker', 'linestyle', 'color'}); + + nline = size(linedata,1); + + % Expand bounds matrix if necessary + + if nline > 1 + if ndims(err{ix}) == 3 + err2 = squeeze(num2cell(err{ix},[1 2])); + else + err2 = repmat(err(ix),nline,1); + end + else + err2 = err(ix); + end + + % Figure out upper and lower bounds + + [lo, hi] = deal(cell(nline,1)); + for iln = 1:nline + + x2 = linedata{iln,1}; + y2 = linedata{iln,2}; + nx = length(x2); + + if isvert + lineval = y2; + else + lineval = x2; + end + + sz = size(err2{iln}); + + if isequal(sz, [nx 2]) + lo{iln} = lineval - err2{iln}(:,1)'; + hi{iln} = lineval + err2{iln}(:,2)'; + elseif isequal(sz, [nx 1]) + lo{iln} = lineval - err2{iln}'; + hi{iln} = lineval + err2{iln}'; + elseif isequal(sz, [1 2]) + lo{iln} = lineval - err2{iln}(1); + hi{iln} = lineval + err2{iln}(2); + elseif isequal(sz, [1 1]) + lo{iln} = lineval - err2{iln}; + hi{iln} = lineval + err2{iln}; + elseif isequal(sz, [2 nx]) % not documented, but accepted anyways + lo{iln} = lineval - err2{iln}(:,1); + hi{iln} = lineval + err2{iln}(:,2); + elseif isequal(sz, [1 nx]) % not documented, but accepted anyways + lo{iln} = lineval - err2{iln}; + hi{iln} = lineval + err2{iln}; + elseif isequal(sz, [2 1]) % not documented, but accepted anyways + lo{iln} = lineval - err2{iln}(1); + hi{iln} = lineval + err2{iln}(2); + else + error('Error bounds must be npt x nside x nline array'); + end + + end + + % Combine all data (xline, yline, marker, linestyle, color, lower bound + % (x or y), upper bound (x or y) + + plotdata = [plotdata; linedata lo hi]; + +end +close(htemp); + +% Override colormap + +if hascmap + nd = size(plotdata,1); + cmap = repmat(cmap, ceil(nd/size(cmap,1)), 1); + cmap = cmap(1:nd,:); + plotdata(:,5) = num2cell(cmap,2); +end + + +%-------------------- +% Plot +%-------------------- + +% Setup of x and y, plus line and patch properties + +nline = size(plotdata,1); +[xl, yl, xp, yp, marker, lnsty, lncol, ptchcol, alpha] = deal(cell(nline,1)); + +for iln = 1:nline + xl{iln} = plotdata{iln,1}; + yl{iln} = plotdata{iln,2}; +% if isvert +% xp{iln} = [plotdata{iln,1} fliplr(plotdata{iln,1})]; +% yp{iln} = [plotdata{iln,6} fliplr(plotdata{iln,7})]; +% else +% xp{iln} = [plotdata{iln,6} fliplr(plotdata{iln,7})]; +% yp{iln} = [plotdata{iln,2} fliplr(plotdata{iln,2})]; +% end + + [xp{iln}, yp{iln}] = calcpatch(plotdata{iln,1}, plotdata{iln,2}, isvert, plotdata{iln,6}, plotdata{iln,7}); + + marker{iln} = plotdata{iln,3}; + lnsty{iln} = plotdata{iln,4}; + + if usealpha + lncol{iln} = plotdata{iln,5}; + ptchcol{iln} = plotdata{iln,5}; + alpha{iln} = trans; + else + lncol{iln} = plotdata{iln,5}; + ptchcol{iln} = interp1([0 1], [1 1 1; lncol{iln}], trans); + alpha{iln} = 1; + end +end + +% Plot patches and lines + +if verLessThan('matlab', '8.4.0') + [hp,hl] = deal(zeros(nline,1)); +else + [hp,hl] = deal(gobjects(nline,1)); +end + +axes(hax); +hold all; + +for iln = 1:nline + hp(iln) = patch(xp{iln}, yp{iln}, ptchcol{iln}, 'facealpha', alpha{iln}, 'edgecolor', 'none'); +end + +for iln = 1:nline + hl(iln) = line(xl{iln}, yl{iln}, 'marker', marker{iln}, 'linestyle', lnsty{iln}, 'color', lncol{iln}); +end + +%-------------------- +% Assign output +%-------------------- + +nargchk(0, 2, nargout); + +if nargout >= 1 + varargout{1} = hl; +end + +if nargout == 2 + varargout{2} = hp; +end + +%-------------------- +% Parse optional +% parameters +%-------------------- + +function [found, val, vars] = parseparam(vars, param) + +isvar = cellfun(@(x) ischar(x) && strcmpi(x, param), vars); + +if sum(isvar) > 1 + error('Parameters can only be passed once'); +end + +if any(isvar) + found = true; + idx = find(isvar); + val = vars{idx+1}; + vars([idx idx+1]) = []; +else + found = false; + val = []; +end + +%---------------------------- +% Calculate patch coordinates +%---------------------------- + +function [xp, yp] = calcpatch(xl, yl, isvert, lo, hi) + +ismissing = any(isnan([xl;yl;lo;hi]),2); + +if isvert + xp = [xl fliplr(xl)]; + yp = [lo fliplr(hi)]; +else + xp = [lo fliplr(hi)]; + yp = [yl fliplr(yl)]; +end + +if any(ismissing) + warning('NaNs in bounds; inpainting'); + xp = inpaint_nans(xp'); + yp = inpaint_nans(yp'); +end + + diff --git a/Figures/Tools/inpaint_nans.m b/Figures/Tools/inpaint_nans.m new file mode 100755 index 0000000..86cc8c4 --- /dev/null +++ b/Figures/Tools/inpaint_nans.m @@ -0,0 +1,411 @@ +function B=inpaint_nans(A,method) +% inpaint_nans: in-paints over nans in an array +% usage: B=inpaint_nans(A) +% +% solves approximation to one of several pdes to +% interpolate and extrapolate holes +% +% arguments (input): +% A - nxm array with some NaNs to be filled in +% +% method - (OPTIONAL) scalar numeric flag - specifies +% which approach (or physical metaphor to use +% for the interpolation.) All methods are capable +% of extrapolation, some are better than others. +% There are also speed differences, as well as +% accuracy differences for smooth surfaces. +% +% methods {0,1,2} use a simple plate metaphor +% methods {3} use a better plate equation, +% but will be slower +% methods 4 use a spring metaphor +% +% method == 0 --> (DEFAULT) see method 1, but +% this method does not build as large of a +% linear system in the case of only a few +% NaNs in a large array. +% Extrapolation behavior is linear. +% +% method == 1 --> simple approach, applies del^2 +% over the entire array, then drops those parts +% of the array which do not have any contact with +% NaNs. Uses a least squares approach, but it +% does not touch existing points. +% In the case of small arrays, this method is +% quite fast as it does very little extra work. +% Extrapolation behavior is linear. +% +% method == 2 --> uses del^2, but solving a direct +% linear system of equations for nan elements. +% This method will be the fastest possible for +% large systems since it uses the sparsest +% possible system of equations. Not a least +% squares approach, so it may be least robust +% to noise on the boundaries of any holes. +% This method will also be least able to +% interpolate accurately for smooth surfaces. +% Extrapolation behavior is linear. +% +% method == 3 --+ See method 0, but uses del^4 for +% the interpolating operator. This may result +% in more accurate interpolations, at some cost +% in speed. +% +% method == 4 --+ Uses a spring metaphor. Assumes +% springs (with a nominal length of zero) +% connect each node with every neighbor +% (horizontally, vertically and diagonally) +% Since each node tries to be like its neighbors, +% extrapolation is as a constant function where +% this is consistent with the neighboring nodes. +% +% +% arguments (output): +% B - nxm array with NaNs replaced + +% I always need to know which elements are NaN, +% and what size the array is for any method +[n,m]=size(A); +nm=n*m; +k=isnan(A(:)); + +% list the nodes which are known, and which will +% be interpolated +nan_list=find(k); +known_list=find(~k); + +% how many nans overall +nan_count=length(nan_list); + +% convert NaN indices to (r,c) form +% nan_list==find(k) are the unrolled (linear) indices +% (row,column) form +[nr,nc]=ind2sub([n,m],nan_list); + +% both forms of index in one array: +% column 1 == unrolled index +% column 2 == row index +% column 3 == column index +nan_list=[nan_list,nr,nc]; + +% supply default method +if (nargin<2)|isempty(method) + method = 0; +elseif ~ismember(method,0:4) + error 'If supplied, method must be one of: {0,1,2,3,4}.' +end + +% for different methods +switch method + case 0 + % The same as method == 1, except only work on those + % elements which are NaN, or at least touch a NaN. + + % horizontal and vertical neighbors only + talks_to = [-1 0;0 -1;1 0;0 1]; + neighbors_list=identify_neighbors(n,m,nan_list,talks_to); + + % list of all nodes we have identified + all_list=[nan_list;neighbors_list]; + + % generate sparse array with second partials on row + % variable for each element in either list, but only + % for those nodes which have a row index > 1 or < n + L = find((all_list(:,2) > 1) & (all_list(:,2) < n)); + nl=length(L); + if nl>0 + fda=sparse(repmat(all_list(L,1),1,3), ... + repmat(all_list(L,1),1,3)+repmat([-1 0 1],nl,1), ... + repmat([1 -2 1],nl,1),nm,nm); + else + fda=spalloc(n*m,n*m,size(all_list,1)*5); + end + + % 2nd partials on column index + L = find((all_list(:,3) > 1) & (all_list(:,3) < m)); + nl=length(L); + if nl>0 + fda=fda+sparse(repmat(all_list(L,1),1,3), ... + repmat(all_list(L,1),1,3)+repmat([-n 0 n],nl,1), ... + repmat([1 -2 1],nl,1),nm,nm); + end + + % eliminate knowns + rhs=-fda(:,known_list)*A(known_list); + k=find(any(fda(:,nan_list(:,1)),2)); + + % and solve... + B=A; + B(nan_list(:,1))=fda(k,nan_list(:,1))\rhs(k); + + case 1 + % least squares approach with del^2. Build system + % for every array element as an unknown, and then + % eliminate those which are knowns. + + % Build sparse matrix approximating del^2 for + % every element in A. + % Compute finite difference for second partials + % on row variable first + [i,j]=ndgrid(2:(n-1),1:m); + ind=i(:)+(j(:)-1)*n; + np=(n-2)*m; + fda=sparse(repmat(ind,1,3),[ind-1,ind,ind+1], ... + repmat([1 -2 1],np,1),n*m,n*m); + + % now second partials on column variable + [i,j]=ndgrid(1:n,2:(m-1)); + ind=i(:)+(j(:)-1)*n; + np=n*(m-2); + fda=fda+sparse(repmat(ind,1,3),[ind-n,ind,ind+n], ... + repmat([1 -2 1],np,1),nm,nm); + + % eliminate knowns + rhs=-fda(:,known_list)*A(known_list); + k=find(any(fda(:,nan_list),2)); + + % and solve... + B=A; + B(nan_list(:,1))=fda(k,nan_list(:,1))\rhs(k); + + case 2 + % Direct solve for del^2 BVP across holes + + % generate sparse array with second partials on row + % variable for each nan element, only for those nodes + % which have a row index > 1 or < n + L = find((nan_list(:,2) > 1) & (nan_list(:,2) < n)); + nl=length(L); + if nl>0 + fda=sparse(repmat(nan_list(L,1),1,3), ... + repmat(nan_list(L,1),1,3)+repmat([-1 0 1],nl,1), ... + repmat([1 -2 1],nl,1),n*m,n*m); + else + fda=spalloc(n*m,n*m,size(nan_list,1)*5); + end + + % 2nd partials on column index + L = find((nan_list(:,3) > 1) & (nan_list(:,3) < m)); + nl=length(L); + if nl>0 + fda=fda+sparse(repmat(nan_list(L,1),1,3), ... + repmat(nan_list(L,1),1,3)+repmat([-n 0 n],nl,1), ... + repmat([1 -2 1],nl,1),n*m,n*m); + end + + % fix boundary conditions at extreme corners + % of the array in case there were nans there + if ismember(1,nan_list(:,1)) + fda(1,[1 2 n+1])=[-2 1 1]; + end + if ismember(n,nan_list(:,1)) + fda(n,[n, n-1,n+n])=[-2 1 1]; + end + if ismember(nm-n+1,nan_list(:,1)) + fda(nm-n+1,[nm-n+1,nm-n+2,nm-n])=[-2 1 1]; + end + if ismember(nm,nan_list(:,1)) + fda(nm,[nm,nm-1,nm-n])=[-2 1 1]; + end + + % eliminate knowns + rhs=-fda(:,known_list)*A(known_list); + + % and solve... + B=A; + k=nan_list(:,1); + B(k)=fda(k,k)\rhs(k); + + case 3 + % The same as method == 0, except uses del^4 as the + % interpolating operator. + + % del^4 template of neighbors + talks_to = [-2 0;-1 -1;-1 0;-1 1;0 -2;0 -1; ... + 0 1;0 2;1 -1;1 0;1 1;2 0]; + neighbors_list=identify_neighbors(n,m,nan_list,talks_to); + + % list of all nodes we have identified + all_list=[nan_list;neighbors_list]; + + % generate sparse array with del^4, but only + % for those nodes which have a row & column index + % >= 3 or <= n-2 + L = find( (all_list(:,2) >= 3) & ... + (all_list(:,2) <= (n-2)) & ... + (all_list(:,3) >= 3) & ... + (all_list(:,3) <= (m-2))); + nl=length(L); + if nl>0 + % do the entire template at once + fda=sparse(repmat(all_list(L,1),1,13), ... + repmat(all_list(L,1),1,13) + ... + repmat([-2*n,-n-1,-n,-n+1,-2,-1,0,1,2,n-1,n,n+1,2*n],nl,1), ... + repmat([1 2 -8 2 1 -8 20 -8 1 2 -8 2 1],nl,1),nm,nm); + else + fda=spalloc(n*m,n*m,size(all_list,1)*5); + end + + % on the boundaries, reduce the order around the edges + L = find((((all_list(:,2) == 2) | ... + (all_list(:,2) == (n-1))) & ... + (all_list(:,3) >= 2) & ... + (all_list(:,3) <= (m-1))) | ... + (((all_list(:,3) == 2) | ... + (all_list(:,3) == (m-1))) & ... + (all_list(:,2) >= 2) & ... + (all_list(:,2) <= (n-1)))); + nl=length(L); + if nl>0 + fda=fda+sparse(repmat(all_list(L,1),1,5), ... + repmat(all_list(L,1),1,5) + ... + repmat([-n,-1,0,+1,n],nl,1), ... + repmat([1 1 -4 1 1],nl,1),nm,nm); + end + + L = find( ((all_list(:,2) == 1) | ... + (all_list(:,2) == n)) & ... + (all_list(:,3) >= 2) & ... + (all_list(:,3) <= (m-1))); + nl=length(L); + if nl>0 + fda=fda+sparse(repmat(all_list(L,1),1,3), ... + repmat(all_list(L,1),1,3) + ... + repmat([-n,0,n],nl,1), ... + repmat([1 -2 1],nl,1),nm,nm); + end + + L = find( ((all_list(:,3) == 1) | ... + (all_list(:,3) == m)) & ... + (all_list(:,2) >= 2) & ... + (all_list(:,2) <= (n-1))); + nl=length(L); + if nl>0 + fda=fda+sparse(repmat(all_list(L,1),1,3), ... + repmat(all_list(L,1),1,3) + ... + repmat([-1,0,1],nl,1), ... + repmat([1 -2 1],nl,1),nm,nm); + end + + % eliminate knowns + rhs=-fda(:,known_list)*A(known_list); + k=find(any(fda(:,nan_list(:,1)),2)); + + % and solve... + B=A; + B(nan_list(:,1))=fda(k,nan_list(:,1))\rhs(k); + + case 4 + % Spring analogy + % interpolating operator. + + % list of all springs between a node and a horizontal + % or vertical neighbor + hv_list=[-1 -1 0;1 1 0;-n 0 -1;n 0 1]; + hv_springs=[]; + for i=1:4 + hvs=nan_list+repmat(hv_list(i,:),nan_count,1); + k=(hvs(:,2)>=1) & (hvs(:,2)<=n) & (hvs(:,3)>=1) & (hvs(:,3)<=m); + hv_springs=[hv_springs;[nan_list(k,1),hvs(k,1)]]; + end + + % delete replicate springs + hv_springs=unique(sort(hv_springs,2),'rows'); + + % build sparse matrix of connections, springs + % connecting diagonal neighbors are weaker than + % the horizontal and vertical springs + nhv=size(hv_springs,1); + springs=sparse(repmat((1:nhv)',1,2),hv_springs, ... + repmat([1 -1],nhv,1),nhv,nm); + + % eliminate knowns + rhs=-springs(:,known_list)*A(known_list); + + % and solve... + B=A; + B(nan_list(:,1))=springs(:,nan_list(:,1))\rhs; + +end + +% ==================================================== +% end of main function +% ==================================================== +% ==================================================== +% begin subfunctions +% ==================================================== +function neighbors_list=identify_neighbors(n,m,nan_list,talks_to) +% identify_neighbors: identifies all the neighbors of +% those nodes in nan_list, not including the nans +% themselves +% +% arguments (input): +% n,m - scalar - [n,m]=size(A), where A is the +% array to be interpolated +% nan_list - array - list of every nan element in A +% nan_list(i,1) == linear index of i'th nan element +% nan_list(i,2) == row index of i'th nan element +% nan_list(i,3) == column index of i'th nan element +% talks_to - px2 array - defines which nodes communicate +% with each other, i.e., which nodes are neighbors. +% +% talks_to(i,1) - defines the offset in the row +% dimension of a neighbor +% talks_to(i,2) - defines the offset in the column +% dimension of a neighbor +% +% For example, talks_to = [-1 0;0 -1;1 0;0 1] +% means that each node talks only to its immediate +% neighbors horizontally and vertically. +% +% arguments(output): +% neighbors_list - array - list of all neighbors of +% all the nodes in nan_list + +if ~isempty(nan_list) + % use the definition of a neighbor in talks_to + nan_count=size(nan_list,1); + talk_count=size(talks_to,1); + + nn=zeros(nan_count*talk_count,2); + j=[1,nan_count]; + for i=1:talk_count + nn(j(1):j(2),:)=nan_list(:,2:3) + ... + repmat(talks_to(i,:),nan_count,1); + j=j+nan_count; + end + + % drop those nodes which fall outside the bounds of the + % original array + L = (nn(:,1)<1)|(nn(:,1)>n)|(nn(:,2)<1)|(nn(:,2)>m); + nn(L,:)=[]; + + % form the same format 3 column array as nan_list + neighbors_list=[sub2ind([n,m],nn(:,1),nn(:,2)),nn]; + + % delete replicates in the neighbors list + neighbors_list=unique(neighbors_list,'rows'); + + % and delete those which are also in the list of NaNs. + neighbors_list=setdiff(neighbors_list,nan_list,'rows'); + +else + neighbors_list=[]; +end + + + + + + + + + + + + + + + + + diff --git a/Figures/Tools/license_panel.txt b/Figures/Tools/license_panel.txt new file mode 100644 index 0000000..d494c4b --- /dev/null +++ b/Figures/Tools/license_panel.txt @@ -0,0 +1,24 @@ +Copyright (c) 2015, Ben Mitch +All rights reserved. + +Redistribution and use in source and binary forms, with or without +modification, are permitted provided that the following conditions are +met: + + * Redistributions of source code must retain the above copyright + notice, this list of conditions and the following disclaimer. + * Redistributions in binary form must reproduce the above copyright + notice, this list of conditions and the following disclaimer in + the documentation and/or other materials provided with the distribution + +THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" +AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE +IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE +ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE +LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR +CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF +SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS +INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN +CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) +ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE +POSSIBILITY OF SUCH DAMAGE. diff --git a/Figures/Tools/outlinebounds.m b/Figures/Tools/outlinebounds.m new file mode 100644 index 0000000..49e3d7b --- /dev/null +++ b/Figures/Tools/outlinebounds.m @@ -0,0 +1,31 @@ +function hnew = outlinebounds(hl, hp) +%OUTLINEBOUNDS Outline the patch of a boundedline +% +% hnew = outlinebounds(hl, hp) +% +% This function adds an outline to the patch objects created by +% boundedline, matching the color of the central line associated with each +% patch. +% +% Input variables: +% +% hl: handles to line objects from boundedline +% +% hp: handles to patch objects from boundedline +% +% Output variables: +% +% hnew: handle to new line objects + +% Copyright 2012 Kelly Kearney + + +hnew = zeros(size(hl)); +for il = 1:numel(hp) + col = get(hl(il), 'color'); + xy = get(hp(il), {'xdata','ydata'}); + ax = ancestor(hl(il), 'axes'); + + hnew(il) = line(xy{1}, xy{2}, 'parent', ax, 'linestyle', '-', 'color', col); +end + diff --git a/Figures/Tools/panel.m b/Figures/Tools/panel.m new file mode 100644 index 0000000..c206fd3 --- /dev/null +++ b/Figures/Tools/panel.m @@ -0,0 +1,5086 @@ + +% Panel is an alternative to Matlab's "subplot" function. +% +% INSTALLATION. To install panel, place the file "panel.m" +% on your Matlab path. +% +% DOCUMENTATION. Scan the introductory information in the +% folder "docs". Learn to use panel by working through the +% demonstration scripts in the folder "demo" (list the demos +% by typing "help panel/demo"). Reference information is +% available through "doc panel" or "help panel". For the +% change log, use "edit panel" to view the file "panel.m". + + + +% CHANGE LOG +% +% ############################################################ +% 22/05/2011 +% First Public Release Version 2.0 +% ############################################################ +% +% 23/05/2011 +% Incorporated an LP solver, since the one we were using +% "linprog()" is not available to users who do not have the +% Optimisation Toolbox installed. +% +% 21/06/2011 +% Added -opdf option, and changed PageSize to be equal to +% PaperPosition. +% +% 12/07/2011 +% Made some linprog optimisations, inspired by "Ian" on +% Matlab central. Tested against subplot using +% demopanel2(N=9). Subplot is faster, by about 20%, but +% panel is better :). For my money, 20% isn't much of a hit +% for the extra functionality. NB: Using Jeff Stuart's +% linprog (unoptimised), panel is much slower (especially +% for large N problems); we will probably have to offer a +% faster solver at some point (optimise Jeff's?). +% +% NOTES: You will see a noticeable delay, also, on resize. +% That's the price of using physical units for the layout, +% because we have to recalculate everything when the +% physical canvas size changes. I suppose in the future, we +% could offer an option so that physical units are only used +% during export; that would make resizes fast, and the user +% may not care so much about layout on screen, if they are +% aiming for print figures. Or, you could have the ability +% to turn off auto-refresh on resize(). +% +% ############################################################ +% 20/07/2011 +% Release Version 2.1 +% ############################################################ +% +% 05/10/2011 +% Tidied in-file documentation (panel.m). +% +% 11/12/2011 +% Added flag "no-manage-font" to constructor, as requested +% by Matlab Central user Mukhtar Ullah. +% +% ############################################################ +% 13/12/2011 +% Release Version 2.2 +% ############################################################ +% +% 21/01/2012 +% Fixed bug in explicit height export option "-hX" which +% wasn't working right at all. +% +% 25/01/12 +% Fixed bug in tick label display during print. _Think_ I've +% got it right, this time! Some notes below, search for +% "25/01/12". +% +% 25/01/12 +% Fixed DPI bug in smoothed export figures. Bug was flagged +% up by Jesper at Matlab Central. +% +% ############################################################ +% 26/01/2012 +% Release Version 2.3 +% ############################################################ +% +% 09/03/12 +% Fixed bug whereby re-positioning never got done if only +% one panel was created in an existing figure window. +% +% ############################################################ +% 13/03/2012 +% Release Version 2.4 +% ############################################################ +% +% 15/03/12 +% NB: On 2008b, and possibly later versions, the fact that +% the resizeCallback() and closeCallback() are private makes +% things not work. You can fix this by removing the "Access +% = Private" modifier on that section of "methods". It works +% fine in later versions, they must have changed the access +% rules I guess. +% +% 19/07/12 +% Modified so that more than one object can be managed by +% one axis. Just use p.select([h1 h2 ...]). Added function +% "getAllManagedAxes()" which returns only objects from the +% "object list" (h_object), as it now is, which represent +% axes. Suggested by Brendan Sullivan @ Matlab Central. +% +% 19/07/12 +% Added support for zlabel() call (not applicable to parent +% panels, since they are implicitly 2D for axis labelling). +% +% 19/07/12 +% Fixed another export bug - how did this one not get +% noticed? XLimMode (etc.) was not getting locked during +% export, so that axes without manual limits might get +% re-dimensioned during export, which is bad news. Added +% locking of limits as well as ticks, in storeAxisState(). +% Hope this has no side effects! +% +% ############################################################ +% 19/07/12 +% Release Version 2.5 +% +% NB: Owing to the introduction of management of multiple +% objects by each panel, this release should be considered +% possibly flaky. Revert to 2.4 if you have problems with +% 2.5. +% ############################################################ +% +% 23/07/12 +% Improved documentation for figure export in demopanelA. +% +% 24/07/12 +% Added support for export to SVG, using "plot2svg" (Matlab +% Central File Exchange) as the renderer. Along the way, +% tidied the behaviour of export() a little, and improved +% reporting to the user. Changed default DPI for EPS to 600, +% since otherwise the output files are pretty shoddy, and +% the filesize is relatively unaffected. +% +% 24/07/12 +% Updated documentation, particularly HTML pages and +% associated figures. Bit nicer, now. +% +% ############################################################ +% 24/07/12 +% Release Version 2.6 +% ############################################################ +% +% 22/09/12 +% Added demopanelH, which illustrates how to do insets. Kudos +% to Ann Hickox for the idea. +% +% 20/03/13 +% Added panel.plot() to work around poor rendering of dashed +% lines, etc. Added demopanelI to illustrate its use. +% +% 20/03/13 +% Renamed setCallback to addCallback, so we can have more +% than one. Added "userdata" argument to addCallback(), and +% "event" field (and "userdata" field) to "data" passed when +% callback is fired. +% +% ############################################################ +% 21/03/13 +% Release Version 2.7 +% ############################################################ +% +% 21/03/13 +% Fixed bug in panel.plot() which did not handle solid lines +% correctly. +% +% 12/04/13 +% Added back setCallback() with appropriate semantics, for +% the use of legacy code (or, really, future code, these +% semantics might be useful to someone). Also added the +% function clearCallbacks(). +% +% 12/04/13 +% Removed panel.plot() because it just seemed to be too hard +% to manage. Instead, we'll let the user plot things in the +% usual way, but during export (when things are well under +% our control), we'll fix up any dashed lines that the user +% has requested using the call fixdash(). Thus, we apply the +% fix only where it's needed, when printing to an image +% file, and save all the faffing with resize callbacks. +% +% ############################################################ +% 12/04/13 +% Release Version 2.8 +% ############################################################ +% +% 13/04/13 +% Changed panel.export() to infer image format from file +% extension, in the case that it is not explicitly specified +% and the passed filename has an extension. +% +% 03/05/13 +% Changed term "render", where misused, to "layout", so as +% not to confuse users of the help/docs. Changed name of +% callback event from "render-complete" to "layout-updated", +% is the only functional effect. +% +% 03/05/13 +% Added argument to panel constructor so that units can be +% set there, rather than through a separate call to the +% "units" property. +% +% 03/05/13 +% Added set descriptor "family" to go with "children" and +% "descendants". This one should be of particular use for +% the construct p.fa.margin = 0. +% +% 03/05/13 +% Updated demopanel9 to be a walkthrough of how to set +% margins. Will be useful to point users at this if they ask +% "how do I do margins?". +% +% 03/05/13 +% Added panel.version(). +% +% 03/05/13 +% Added page size "LNCS" to export. +% +% ############################################################ +% 03/05/13 +% Release Version 2.9 +% ############################################################ +% +% 10/05/13 +% Removed linprog solution in favour of recursive +% computation. This should speed things up for people who +% don't have the optimisation toolbox. +% +% 10/05/13 +% Added support for panels of fixed physical size. See new +% documentation for panel/pack(). +% +% 10/05/13 +% Added support for packing into panels packed in absolute +% mode, which wasn't previously possible. +% +% 10/05/13 +% Removed advertisement for 'defer' flag, since I suspect +% it's no longer needed now we've moved away from LP. There +% may be some optimisation required before this is true - +% defer still functions as before, it's just not advertised. +% +% ############################################################ +% 10/05/13 +% Release Version 2.10 +% ############################################################ +% +% 14/05/13 +% Some minor optimisations, so now panel is not slower than +% subplot (see demopanelK). +% +% 14/01/15 +% Various fixes to work correctly under R2014b. Essentially, +% checked the demos, added retirement notes to fixdash(), and +% added function "fignum()". +% +% ############################################################ +% 14/01/15 +% Release Version 2.11 +% ############################################################ +% +% 28/03/15 +% Changed export() logic slightly so that if either -h or -w option is +% specified, direct sizing model is selected (and, therefore, /all/ +% options from the paper sizing model are ignored). Thus, either -w or +% -h can be specified, with -a, and intuitively-correct behaviour +% results. +% +% 02/04/15 +% Changed functions x/y/zlabel and title to return a handle to the +% referenced object so that caller can access its properties. +% +% ############################################################ +% 02/04/15 +% Release Version 2.12 +% ############################################################ + + + +classdef (Sealed = true) panel < handle + + + + %% ---- PROPERTIES ---- + + properties (Constant = true, Hidden = true) + + PANEL_TYPE_UNCOMMITTED = 0; + PANEL_TYPE_PARENT = 1; + PANEL_TYPE_OBJECT = 2; + + end + + properties (Constant = true) + + LAYOUT_MODE_NORMAL = 0; + LAYOUT_MODE_PREPRINT = 1; + LAYOUT_MODE_POSTPRINT = 2; + + end + + properties + + % these properties are only here for documentation. they + % are actually stored in "prop". it's just some subsref + % madness. + + % font name to use for axis text (inherited) + % + % access: read/write + % default: normal + fontname + + % font size to use for axis text (inherited) + % + % access: read/write + % default: normal + fontsize + + % font weight to use for axis text (inherited) + % + % access: read/write + % default: normal + fontweight + + % the units that are used when reading/writing margins + % + % units can be set to any of 'mm', 'cm', 'in' or 'pt'. + % it only affects the read/write interface; values + % stored already are not re-interpreted. + % + % access: read/write + % default: mm + units + + % the panel's margin vector in the form [left bottom right top] + % + % the margin is key to the layout process. the layout + % algorithm makes all panels as large as possible whilst + % not violating margin constraints. margins are + % respected between panels within their parent and + % between the root panel and the edges of the canvas + % (figure or image file). + % + % access: read/write + % default: [12 10 2 2] (mm) + % + % see also: marginleft, marginbottom, marginright, margintop + margin + + % one element of the margin vector + % + % access: read/write + % default: see margin + % + % see also: margin + marginleft + + % one element of the margin vector + % + % access: read/write + % default: see margin + % + % see also: margin + marginbottom + + % one element of the margin vector + % + % access: read/write + % default: see margin + % + % see also: margin + marginright + + % one element of the margin vector + % + % access: read/write + % default: see margin + % + % see also: margin + margintop + + % return position of panel + % + % return the panel's position in normalized + % coordinates (normalized to the figure window that + % is associated with the panel). note that parent + % panels have positions too, even though nothing is + % usually rendered. uncommitted panels, too. + % + % access: read only + position + + % return handle of associated figure + % + % access: read only + figure + + % return handle of associated axis + % + % if the panel is not an axis panel, empty is returned. + % object includes axis, but axis does not include + % object. + % + % access: read only + % + % see also: object + axis + + % return handle of associated object + % + % if the panel is not an object panel, empty is + % returned. object includes axis, but axis does not + % include object. + % + % access: read only + % + % see also: axis + object + + % access properties of panel's children + % + % if the panel is a parent panel, "children" gives + % access to some properties of its children (direct + % descendants). "children" can be abbreviated "ch". + % properties that can be accessed are as follows. + % + % axis: read-only, returns an array + % object: read-only, returns an array + % + % margin: write-only + % fontname: write-only + % fontsize: write-only + % fontweight: write-only + % + % EXAMPLE: + % h = p.ch.axis; + % p.ch.margin = 3; + % + % see also: descendants, family + children + + % access properties of panel's descendants + % + % if the panel is a parent panel, "descendants" gives + % access to some properties of its descendants + % (children, grandchildren, etc.). "descendants" can be + % abbreviated "de". properties that can be accessed are + % as follows. + % + % axis: read-only, returns an array + % object: read-only, returns an array + % + % margin: write-only + % fontname: write-only + % fontsize: write-only + % fontweight: write-only + % + % EXAMPLE: + % h = p.de.axis; + % p.de.margin = 3; + % + % see also: children, family + descendants + + % access properties of panel's family + % + % if the panel is a parent panel, "family" gives access + % to some properties of its family (self, children, + % grandchildren, etc.). "family" can be abbreviated + % "fa". properties that can be accessed are as follows. + % + % axis: read-only, returns an array + % object: read-only, returns an array + % + % margin: write-only + % fontname: write-only + % fontsize: write-only + % fontweight: write-only + % + % EXAMPLE: + % h = p.fa.axis; + % p.fa.margin = 3; + % + % see also: children, descendants + family + + end + + properties (Access = private) + + % associated figure window + h_figure + + % parent graphics object + h_parent + + % this is empty for the root PANEL, populated for all others + parent + + % this is always the root panel associated with this + m_root + + % packing specifier + % + % empty: relative positioning mode (stretch) + % scalar fraction: relative positioning mode + % scalar percentage: relative positioning mode + % 1x4 dimension: absolute positioning mode + packspec + + % packing dimension of children + % + % 1 : horizontal + % 2 : vertical + packdim + + % panel type + m_panelType + + % fixdash lines + m_fixdash + m_fixdash_restore + + % associated managed graphics object (usually, an axis) + h_object + + % show axis (only the root has this extra axis, if show() is active) + h_showAxis + + % children (only a parent panel has non-empty, here) + m_children + + % callback (any functions listed in this cell array are called when events occur) + m_callback + + % local properties (actual properties is this overlaid on inherited/default properties) + % + % see getPropertyValue() + prop + + % state + % + % private state information used during various processing + state + + % layout context for this panel + % + % this is the layout context for the panel. this is + % computed in the function recomputeLayout(), and used + % to reposition the panel in applyLayout(). storage of + % this data means that we can call applyLayout() to + % layout only a branch of the panel tree without having + % to recompute the whole thing. however, I don't know + % how efficient this system is, might need some work. + m_context + + end + + + + + + + + %% ---- PUBLIC CTOR/DTOR ---- + + methods + + function p = panel(varargin) + + % create a new panel + % + % p = panel(...) + % create a new root panel. optional arguments listed + % below can be supplied in any order. if "h_parent" + % is not supplied, it is set to gcf - that is, the + % panel fills the current figure. + % + % initially, the root panel is an "uncommitted + % panel". calling pack() or select() on it will + % commit it as a "parent panel" or an "object + % panel", respectively. the following arguments may + % be passed, in any order. + % + % h_parent + % a handle to a graphics object that will act as the + % parent of the new panel. this is usually a figure + % handle, but may be a handle to any graphics + % object, in principle. currently, an error is + % raised unless it's a figure or a uipanel - if you + % want to try other object types, edit the code + % where the error is raised, and let me know if you + % have positive results so I can update panel to + % allow other object types. + % + % 'add' + % usually, when you attach a new root panel to a + % figure, any existing attached root panels are + % first deleted to make way for it. if you pass this + % argument, this is not done, so that you can attach + % more than one root panel to the same figure. see + % demopanelE for an example of this use. + % + % 'no-manage-font' + % by default, a panel will manage fonts of titles + % and axis labels. this prevents the user from + % setting individual fonts on those items. pass this + % flag to disable font management for this panel. + % + % 'mm', 'cm', 'in', 'pt' + % by default, panel uses mm as the unit of + % communication with the user over margin sizes. + % pass any of these to change this (you can achieve + % the same effect after creating a panel by setting + % the property "units"). + % + % see also: panel (overview), pack(), select() + + % PRIVATE DOCUMENTATION + % + % 'defer' + % THIS IS NO LONGER ADVERTISED since we replaced the + % LP solution with a procedural solution, but still + % functions as before, to provide legacy support. + % the panel will be created with layout disabled. + % the layout computations take a little while when + % large numbers of panels are involved, and are + % re-run every time you add a panel or change a + % margin, by default. this is tedious if you are + % preparing a complex layout; pass 'defer', and + % layout will not be computed at all until you call + % refresh() or export() on the root panel. + % + % 'pack' + % this constructor is called internally from pack() + % to create new panels when packing them into + % parents. the first argument is passed as 'pack' in + % this case, which allows us to do slightly quicker + % parsing of the arguments, since we know the + % calling convention (see immediately below). + + % default state + p.state = []; + p.state.name = ''; + p.state.defer = 0; + p.state.manage_font = 1; + p.m_callback = {}; + p.m_fixdash = {}; + p.packspec = []; + p.packdim = 2; + p.m_panelType = p.PANEL_TYPE_UNCOMMITTED; + p.prop = panel.getPropertyInitialState(); + + % handle call from pack() aqap + if nargin && isequal(varargin{1}, 'pack') + + % since we know the calling convention, in this + % case, we can handle this as quickly as possible, + % so that large recursive layouts do not get held up + % by spurious code, here. + + % parent is a panel + passed_h_parent = varargin{2}; + + % become its child + indexInParent = int2str(length(passed_h_parent.m_children)+1); + if passed_h_parent.isRoot() + p.state.name = ['(' indexInParent ')']; + else + p.state.name = [passed_h_parent.state.name(1:end-1) ',' indexInParent ')']; + end + p.h_parent = passed_h_parent.h_parent; + p.h_figure = passed_h_parent.h_figure; + p.parent = passed_h_parent; + p.m_root = passed_h_parent.m_root; + + % done! + return + + end + + % default condition + passed_h_parent = []; + add = false; + + % peel off args + while ~isempty(varargin) + + % get arg + arg = varargin{1}; + varargin = varargin(2:end); + + % handle text + if ischar(arg) + + switch arg + + case 'add' + add = true; + continue; + + case 'defer' + p.state.defer = 1; + continue; + + case 'no-manage-font' + p.state.manage_font = 0; + continue; + + case {'mm' 'cm' 'in' 'pt'} + p.setPropertyValue('units', arg); + continue; + + otherwise + error('panel:InvalidArgument', ['unrecognised text argument "' arg '"']); + + end + + end + + % handle parent + if isscalar(arg) && ishandle(arg) + passed_h_parent = arg; + continue; + end + + % error + error('panel:InvalidArgument', 'unrecognised argument to panel constructor'); + + end + + % attach to current figure if no parent supplied + if isempty(passed_h_parent) + passed_h_parent = gcf; + + % this might cause a figure to be created - if so, + % give it time to display now so we don't get a (or + % two, in fact!) resize event(s) later + drawnow + end + + % we are a root panel + p.state.name = 'root'; + p.parent = []; + p.m_root = p; + + % get parent type + parentType = get(passed_h_parent, 'type'); + + % set handles + switch parentType + + case 'uipanel' + p.h_parent = passed_h_parent; + p.h_figure = getParentFigure(passed_h_parent); + + case 'figure' + p.h_parent = passed_h_parent; + p.h_figure = passed_h_parent; + + otherwise + error('panel:InvalidArgument', ... + ['panel() cannot be attached to an object of type "' parentType '"']); + + end + + % lay in callbacks + addHandleCallback(p.h_figure, 'CloseRequestFcn', @panel.closeCallback); + addHandleCallback(p.h_parent, 'ResizeFcn', @panel.resizeCallback); + + % register for callbacks + if add + panel.callbackDispatcher('registerNoClear', p); + else + panel.callbackDispatcher('register', p); + end + + % lock class in memory (prevent persistent from being cleared by clear all) + panel.lockClass(); + + end + + function delete(p) + + % destroy a panel + % + % delete(p) + % destroy the passed panel, deleting all associated + % graphics objects. + % + % NB: you won't usually have to call this explicitly. + % it is called automatically for all attached panels + % when you close the associated figure. + + % debug output +% panel.debugmsg(['deleting "' p.state.name '"...']); + + % delete managed graphics objects + for n = 1:length(p.h_object) + h = p.h_object(n); + if ishandle(h) + delete(h); + end + end + + % delete associated show axis + if ~isempty(p.h_showAxis) && ishandle(p.h_showAxis) + delete(p.h_showAxis); + end + + % delete all children (child will remove itself from + % "m_children" on delete()) + while ~isempty(p.m_children) + delete(p.m_children(end)); + end + + % unregister... + if p.isRoot() + + % ...for callbacks + panel.callbackDispatcher('unregister', p); + + else + + % ...from parent + p.parent.removeChild(p); + + end + + % debug output +% panel.debugmsg(['deleted "' p.state.name '"!']); + + end + + end + + + + + + + + + + + + + + + %% ---- PUBLIC DISPLAY ---- + + methods (Hidden = true) + + function disp(p) + + display(p); + + end + + function display(p, indent) + + % default + if nargin < 2 + indent = ''; + end + + % handle non-scalar (should not exist!) + nels = numel(p); + if nels > 1 + sz = size(p); + sz = sprintf('%dx', sz); + disp([sz(1:end-1) ' array of panel objects']); + return + end + + % header + header = indent; + if p.isObject() + header = [header 'Object ' p.state.name ': ']; + elseif p.isParent() + header = [header 'Parent ' p.state.name ': ']; + else + header = [header 'Uncommitted ' p.state.name ': ']; + end + if p.isRoot() + pp = ['attached to Figure ' panel.fignum(p.h_figure)]; + else + if isempty(p.packspec) + pp = 'stretch'; + elseif iscell(p.packspec) + units = p.getPropertyValue('units'); + val = panel.resolveUnits({p.packspec{1} 'mm'}, units); + pp = sprintf('%.1f%s', val, units); + elseif isscalar(p.packspec) + if p.packspec > 1 + pp = sprintf('%.0f%%', p.packspec); + else + pp = sprintf('%.3f', p.packspec); + end + else + pp = sprintf('%.3f ', p.packspec); + pp = pp(1:end-1); + end + end + header = [header '[' pp]; + if p.isParent() + edges = {'hori' 'vert'}; + header = [header ', ' edges{p.packdim}]; + end + header = [header ']']; + + % margin + header = rpad(header, 60); + header = [header '[ margin ' sprintf('%.3g ', p.getPropertyValue('margin')) p.getPropertyValue('units') ']']; + +% % index +% if isfield(p.state, 'index') +% header = [header ' (' int2str(p.state.index) ')']; +% end + + % display + disp(header); + + % children + for c = 1:length(p.m_children) + p.m_children(c).display([indent ' ']); + end + + end + + end + + + + + + + + + + + %% ---- PUBLIC METHODS ---- + + methods + + function h = xlabel(p, text) + + % apply an xlabel to the panel (or group) + % + % p.xlabel(...) + % behaves just like xlabel() at the prompt (you can + % use that as an alternative) when called on an axis + % panel. when called on a parent panel, however, the + % group of objects within that parent have a label + % applied. when called on a non-axis object panel, + % an error is raised. + + h = get(p.getOrCreateAxis(), 'xlabel'); + set(h, 'string', text); + if p.isParent() + set(h, 'visible', 'on'); + end + + end + + function h = ylabel(p, text) + + % apply a ylabel to the panel (or group) + % + % p.ylabel(...) + % behaves just like ylabel() at the prompt (you can + % use that as an alternative) when called on an axis + % panel. when called on a parent panel, however, the + % group of objects within that parent have a label + % applied. when called on a non-axis object panel, + % an error is raised. + + h = get(p.getOrCreateAxis(), 'ylabel'); + set(h, 'string', text); + if p.isParent() + set(h, 'visible', 'on'); + end + + end + + function h = zlabel(p, text) + + % apply a zlabel to the panel (or group) + % + % p.zlabel(...) + % behaves just like zlabel() at the prompt (you can + % use that as an alternative) when called on an axis + % panel. when called on a parent panel, however, + % this method raises an error, since parent panels + % are assumed to be 2D, with respect to axes. + + if p.isParent() + error('panel:ZLabelOnParentAxis', 'can only call zlabel() on an object panel'); + end + + h = get(p.getOrCreateAxis(), 'zlabel'); + set(h, 'string', text); + + end + + function h = title(p, text) + + % apply a title to the panel (or group) + % + % p.title(...) + % behaves just like title() at the prompt (you can + % use that as an alternative) when called on an axis + % panel. when called on a parent panel, however, the + % group of objects within that parent have a title + % applied. when called on a non-axis object panel, + % an error is raised. + + h = title(p.getOrCreateAxis(), text); + if p.isParent() + set(h, 'visible', 'on'); + end + + end + + function hold(p, state) + + % set the hold on/off state of the associated axis + % + % p.hold('on' / 'off') + % you can use matlab's "hold" function with plots in + % panel, just like any other plot. there is, + % however, a very minor gotcha that is somewhat + % unlikely to ever come up, but for completeness + % this is the problem and the solutions: + % + % if you create a panel "p", change its font using + % panel, e.g. "p.fontname = 'Courier New'", then call + % "hold on", then "hold off", then plot into it, the + % font is not respected. this situation is unlikely to + % arise because there's usually no reason to call + % "hold off" on a plot. however, there are three + % solutions to get round it, if it does: + % + % a) call p.refresh() when you're finished, to + % update all fonts to managed values. + % + % b) if you're going to call p.export() anyway, + % fonts will get updated when you do. + % + % c) if for some reason you can't do (a) OR (b) (I + % can't think why), you can use the hold() function + % provided by panel instead of that provided by + % Matlab. this will not affect your fonts. for + % example, call "p(2).hold('on')". + + % because the matlab "hold off" command sets an axis's + % nextplot state to "replace", we lose control over + % the axis properties (such as fontname). we set + % nextplot to "replacechildren" when we create an + % axis, but if the user does a "hold on, hold off" + % cycle, we lose that. therefore, we offer this + % alternative. + + % check + if ~p.isObject() + error('panel:HoldWhenNotObjectPanel', 'can only call hold() on an object panel'); + end + + % check + h_axes = p.getAllManagedAxes(); + if isempty(h_axes) + error('panel:HoldWhenNoAxes', 'can only call hold() on a panel that manages one or more axes'); + end + + % switch + switch state + case {'on' true 1} + set(h_axes, 'nextplot', 'add'); + case {'off' false 0} + set(h_axes, 'nextplot', 'replacechildren'); + otherwise + error('panel:InvalidArgument', 'argument to hold() must be ''on'', ''off'', or boolean'); + end + + end + + function fixdash(p, hs, linestyle) + + % pass dashed lines to be fixed up during export + % + % NB: Matlab's difficulty with dotted/dashed lines on export + % seems to be fixed in R2014b, so if using this version or a + % later one, this functionality of panel will be of no + % interest. Text below was from pre R2014b. + % + % p.fixdash(h, linestyle) + % add the lines specified as handles in "h" to the + % list of lines to be "fixed up" during export. + % panel will attempt to get the lines to look right + % during export to all formats where they would + % usually get mussed up. see demopanelI for an + % example of how it works. + % + % the above is the usual usage of fixdash(), but + % you can get more control over linestyle by + % specifying the additional argument, "linestyle". + % if "linestyle" is supplied, it is used as the + % linestyle; if not, the current linestyle of the + % line (-, --, -., :) is used. "linestyle" can + % either be a text string or a series of numbers, as + % described below. + % + % '-' solid + % '--' dashed, equal to [2 0.75] + % '-.' dash-dot, equal to [2 0.75 0.5 0.75] + % ':', '.' dotted, equal to [0.5 0.5] + % + % a number series should be 1xN, where N is a + % multiple of 2, as in the examples above, and + % specifies the lengths of any number of dash + % components that are used before being repeated. + % for instance, '-.' generates a 2 unit segment + % (dash), a 0.75 unit gap, then a 0.5 unit segment + % (dot) and a final 0.75 unit gap. at present, the + % units are always millimetres. this system is + % extensible, so that the following examples are + % also valid: + % + % '--..' dash-dash-dot-dot + % '-..-.' dash-dot-dot-dash-dot + % [2 1 4 1 6 1] 2 dash, 4 dash, 6 dash + + % default + if nargin < 3 + linestyle = []; + end + + % bubble up to root + if ~p.isRoot() + p.m_root.fixdash(hs, linestyle); + return + end + + % for each passed handle + for h = (hs(:)') + + % check it's still a handle + if ~ishandle(h) + continue + end + + % check it's a line + if ~isequal(get(h, 'type'), 'line') + continue + end + + % update if in list + found = false; + for i = 1:length(p.m_fixdash) + if h == p.m_fixdash{i}.h + p.m_fixdash{i}.linestyle = linestyle; + found = true; + break + end + end + + % else add to list + if ~found + p.m_fixdash{end+1} = struct('h', h, 'linestyle', linestyle); + end + + end + + end + + function show(p) + + % highlight the outline of the panel + % + % p.show() + % make the outline of the panel "p" show up in red + % in the figure window. this is useful for + % understanding a complex layout. + % + % see also: identify() + + r = p.getObjectPosition(); + + if ~isempty(r) + h = p.getShowAxis(); + delete(get(h, 'children')); + xdata = [r(1) r(1)+r(3) r(1)+r(3) r(1) r(1)]; + ydata = [r(2) r(2) r(2)+r(4) r(2)+r(4) r(2)]; + plot(h, xdata, ydata, 'r-', 'linewidth', 5); + axis([0 1 0 1]) + end + + end + + function export(p, varargin) + + % to export the root panel to an image file + % + % p.export(filename, ...) + % + % export the figure containing panel "p" to an image file. + % you must supply the filename of this output file, with or + % without a file extension. any further arguments must be + % option strings starting with the dash ("-") character. "p" + % should be the root panel. + % + % if the filename does not include an extension, the + % appropriate extension will be added. if it does, the + % output format will be inferred from it, unless overridden + % by the "-o" option, described below. + % + % if you are targeting a print publication, you may find it + % easiest to size your output using the "paper sizing model" + % (below). if you prefer, you can use the "direct sizing + % model", instead. these two sizing models are described + % below. underneath these are listed the options unrelated + % to sizing (which apply regardless of which sizing model + % you use). + % + % + % + % PAPER SIZING MODEL: + % + % using the paper sizing model, you specify your target as a + % region of a piece of paper, and the actual size in + % millimeters is calculated for you. this is usually very + % convenient, but if you find it unsuitable, the direct + % sizing model (next section) is provided as an alternative. + % + % to specify the region, you specify the type (size) of + % paper, the orientation, the number of columns, and the + % aspect ratio of the figure (or the fraction of a column to + % fill). usually, the remaining options can be left as + % defaults. + % + % -pX + % X is the paper type, A2-A6 or letter (default is A4). + % NB: you can also specify paper type LNCS (Lecture Notes + % in Computer Science), using "-pLNCS". If you do this, + % the margins are also adjusted to match LNCS format. + % + % -l + % specify landscape mode (default is portrait). + % + % -mX + % X is the paper margins in mm. you can provide a scalar + % (same margins all round) or a comma-separated list of + % four values, specifying the left, bottom, right, top + % margins separately (default is 20mm all round). + % + % -iX + % X is the inter-column space in mm (default is + % 5mm). + % + % -cX + % X is the number of columns (default is 1). + % + % NB: the following two options represent two ways to + % specify the height of the figure relative to the space + % defined by the above options. if you supply both, + % whichever comes second will be used. + % + % -aX + % X is the aspect ratio of the resulting image file (width + % is set by the paper model). X can be one of the strings: + % s (square), g (landscape golden ratio), gp (portrait + % golden ratio), h (half-height), d (double-height); or, a + % number greater than zero, to specify the aspect ratio + % explicitly. note that, if using the numeric form, the + % ratio is expressed as the quotient of width over height, + % in the usual way. ratios greater than 10 or less than + % 0.1 are disallowed, since these can cause a very large + % figure file to be created accidentally. default is to + % use the landscape golden ratio. + % + % -fX + % X is the fraction of the column (or page, if there are + % not columns) to fill. X can be one of the following + % strings - a (all), tt (two thirds), h (half), t (third), + % q (quarter) - or a fraction between 0 and 1, to specify + % the fraction of the space to fill as a number. default + % is to use aspect ratio, not fill fraction. + % + % + % + % DIRECT SIZING MODEL: + % + % if one of these two options is set, the output image is + % sized according to that option and the aspect ratio (see + % above) and the paper model is not used. if both are set, + % the aspect ratio is not used either. + % + % -wX + % X is the explicit width in mm. + % + % -hX + % X is the explicit height in mm. + % + % + % + % NON-SIZING OPTIONS: + % + % finally, a few options are provided to control how + % the prepared figure is exported. note that DPI below + % 150 is only recommended for sizing drafts, since + % font and line sizes are not rendered even vaguely + % accurately in some cases. at the other end, DPI + % above 600 is unlikely to be useful except when + % submitting camera-ready copy. + % + % -rX + % X is the resolution (DPI) at which to produce the + % output file. X can be one of the following strings + % - d (draft, 75DPI), n (normal, 150DPI), h (high, + % 300DPI), p (publication quality, 600DPI), x + % (extremely high quality, 1200DPI) - or just + % the DPI as a number (must be in 75-2400). the + % default depends on the output format (see below). + % + % -rX/S + % X is the DPI, S is the smoothing factor, which can + % be 2 or 4. the output file is produced at S times + % the specified DPI, and then reduced in size to the + % specified DPI by averaging. thus, the hard edges + % produced by the renderer are smoothed - the effect + % is somewhat like "anti-aliasing". + % + % NB: the DPI setting might be expected to have no + % effect on vector formats. this is true for SVG, but + % not for EPS, where the DPI affects the numerical + % precision used as well as the size of some image + % elements, but has little effect on file size. for + % this reason, the default DPI is 150 for bitmap + % formats but 600 for vector formats. + % + % -s + % print sideways (default is to print upright) + % + % -oX + % X is the output format - choose from most of the + % built-in image device drivers supported by "print" + % (try "help print"). this includes "png", "jpg", + % "tif", "eps" and "pdf". note that "eps"/"ps" + % resolve to "epsc2"/"psc2", for convenience. to use + % the "eps"/"ps" devices, use "-oeps!"/"-ops!". you + % may also specify "svg", if you have the tool + % "plot2svg" on your path (available at Matlab + % Central File Exchange). the default output format + % is inferred from the file extension, or "png" if + % the passed filename has no extension. + % + % + % + % EXAMPLES: + % + % default export of 'myfig', creates 'myfig.png' at a + % size of 170x105mm (1004x620px). this size comes + % from: A4 (210mm wide), minus two 20mm margins + % (170mm), and using the golden aspect ratio to give a + % height of 105mm, and finally 150DPI to give the + % pixel size. + % + % p.export('myfig') + % + % when producing the final camera-ready image for a + % square figure that will sit in one of the two + % columns of a letter-size paper journal with default + % margins and inter-column space, we might use this: + % + % p.export('myfig', '-pletter', '-c2', '-as', '-rp'); + + % LEGACY + % + % (this is legacy since the 'defer' flag is no longer + % needed - though it is still supported) + % + % NB: if you pass 'defer' to the constructor, calling + % export() both exports the panel and releases the + % defer mode. future changes to properties (e.g. + % margins) will cause immediate recomputation of the + % layout. + + % check + if ~p.isRoot() + error('panel:ExportWhenNotRoot', 'cannot export() this panel - it is not the root panel'); + end + + % used below + default_margin = 20; + + % parameters + pars = []; + pars.filename = ''; + pars.fmt = ''; + pars.ext = ''; + pars.dpi = []; + pars.smooth = 1; + pars.paper = 'A4'; + pars.landscape = false; + pars.fill = -1.618; + pars.cols = 1; + pars.intercolumnspacing = 5; + pars.margin = default_margin; + pars.sideways = false; + pars.width = 0; + pars.height = 0; + invalid = false; + + % interpret args + for a = 1:length(varargin) + + % extract + arg = varargin{a}; + + % all arguments must be non-empty strings + if ~isstring(arg) + error('panel:InvalidArgument', ... + 'all arguments to export() must be non-empty strings'); + end + + % is filename? + if arg(1) ~= '-' + + % error if already set + if ~isempty(pars.filename) + error('panel:InvalidArgument', ... + ['at argument "' arg '", the filename is already set ("' pars.filename '")']); + end + + % ok, continue + pars.filename = arg; + continue + + end + + % split off option key and option value + if length(arg) < 2 + error('panel:InvalidArgument', ... + ['at argument "' arg '", no option specified']); + end + key = arg(2); + val = arg(3:end); + + % switch on option key + switch key + + case 'p' + pars.paper = validate_par(val, arg, {'A2' 'A3' 'A4' 'A5' 'A6' 'letter' 'LNCS'}); + + case 'l' + pars.landscape = true; + validate_par(val, arg, 'empty'); + + case 'm' + pars.margin = validate_par(str2num(val), arg, 'dimension', 'nonneg'); + + case 'i' + pars.intercolumnspacing = validate_par(str2num(val), arg, 'scalar', 'nonneg'); + + case 'c' + pars.cols = validate_par(str2num(val), arg, 'scalar', 'integer'); + + case 'f' + switch val + case 'a', pars.fill = 1; % all + case 'w', pars.fill = 1; % whole (legacy, not documented) + case 'tt', pars.fill = 2/3; % two thirds + case 'h', pars.fill = 1/2; % half + case 't', pars.fill = 1/3; % third + case 'q', pars.fill = 1/4; % quarter + otherwise + pars.fill = validate_par(str2num(val), arg, 'scalar', [0 1]); + end + + case 'a' + switch val + case 's', pars.fill = -1; % square + case 'g', pars.fill = -1.618; % golden ratio (landscape) + case 'gp', pars.fill = -1/1.618; % golden ratio (portrait) + case 'h', pars.fill = -2; % half height + case 'd', pars.fill = -0.5; % double height + otherwise + pars.fill = -validate_par(str2num(val), arg, 'scalar', [0.1 10]); + end + + case 'w' + pars.width = validate_par(str2num(val), arg, 'scalar', 'nonneg', [10 Inf]); + + case 'h' + pars.height = validate_par(str2num(val), arg, 'scalar', 'nonneg', [10 Inf]); + + case 'r' + % peel off smoothing ("/...") + if any(val == '/') + f = find(val == '/', 1); + switch val(f+1:end) + case '2', pars.smooth = 2; + case '4', pars.smooth = 4; + otherwise, error('panel:InvalidArgument', ... + ['invalid argument "' arg '", part after / must be "2" or "4"']); + end + val = val(1:end-2); + end + + switch val + case 'd', pars.dpi = 75; % draft + case 'n', pars.dpi = 150; % normal + case 'h', pars.dpi = 300; % high + case 'p', pars.dpi = 600; % publication quality + case 'x', pars.dpi = 1200; % extremely high quality + otherwise + pars.dpi = validate_par(str2num(val), arg, 'scalar', [75 2400]); + end + + case 's' + pars.sideways = true; + validate_par(val, arg, 'empty'); + + case 'o' + fmts = { + 'png' 'png' 'png' + 'tif' 'tiff' 'tif' + 'tiff' 'tiff' 'tif' + 'jpg' 'jpeg' 'jpg' + 'jpeg' 'jpeg' 'jpg' + 'ps' 'psc2' 'ps' + 'ps!' 'psc' 'ps' + 'psc' 'psc' 'ps' + 'ps2' 'ps2' 'ps' + 'psc2' 'psc2' 'ps' + 'eps' 'epsc2' 'eps' + 'eps!' 'eps' 'eps' + 'epsc' 'epsc' 'eps' + 'eps2' 'eps2' 'eps' + 'epsc2' 'epsc2' 'eps' + 'pdf' 'pdf' 'pdf' + 'svg' 'svg' 'svg' + }; + validate_par(val, arg, fmts(:, 1)'); + index = isin(fmts(:, 1), val); + pars.fmt = fmts(index, 2:3); + + otherwise + error('panel:InvalidArgument', ... + ['invalid argument "' argtext '", option is not recognised']); + + end + + end + + % if not specified, infer format from filename + if isempty(pars.fmt) + [path, base, ext] = fileparts(pars.filename); + if ~isempty(ext) + ext = ext(2:end); + end + switch ext + case {'tif' 'tiff'} + pars.fmt = {'tiff' 'tif'}; + case {'jpg' 'jpeg'} + pars.fmt = {'jpeg' 'jpg'}; + case 'eps' + pars.fmt = {'epsc2' 'eps'}; + case {'png' 'pdf' 'svg'} + pars.fmt = {ext ext}; + case '' + pars.fmt = {'png' 'png'}; + otherwise + warning('panel:CannotInferImageFormat', ... + ['unable to infer image format from file extension "' ext '" (PNG assumed)']); + pars.fmt = {'png' 'png'}; + end + end + + % extract + pars.ext = pars.fmt{2}; + pars.fmt = pars.fmt{1}; + + % extract + is_bitmap = ismember(pars.fmt, {'png' 'jpeg' 'tiff'}); + + % default DPI + if isempty(pars.dpi) + if is_bitmap + pars.dpi = 150; + else + pars.dpi = 600; + end + end + + % validate + if isequal(pars.fmt, 'svg') && isempty(which('plot2svg')) + error('panel:Plot2SVGMissing', 'export to SVG requires plot2svg (Matlab Central File Exchange)'); + end + + % validate + if ~is_bitmap && pars.smooth ~= 1 + pars.smooth = 1; + warning('panel:NoSmoothVectorFormat', 'requested smoothing will not be performed (chosen export format is not a bitmap format)'); + end + + % validate + if isempty(pars.filename) + error('panel:InvalidArgument', 'filename not supplied'); + end + + % make sure filename has extension + if ~any(pars.filename == '.') + pars.filename = [pars.filename '.' pars.ext]; + end + + + +%%%% GET TARGET DIMENSIONS (BEGIN) + + % get space for figure + switch pars.paper + case 'A0', sz = [841 1189]; + case 'A1', sz = [594 841]; + case 'A2', sz = [420 594]; + case 'A3', sz = [297 420]; + case 'A4', sz = [210 297]; + case 'A5', sz = [148 210]; + case 'A6', sz = [105 148]; + case 'letter', sz = [216 279]; + case 'LNCS', sz = [152 235]; + % if margin is still at default, set it to LNCS + % margin size + if isequal(pars.margin, default_margin) + pars.margin = [15 22 15 20]; + end + otherwise + error(['unrecognised paper size "' pars.paper '"']) + end + + % orientation of paper + if pars.landscape + sz = sz([2 1]); + end + + % paper margins (scalar or quad) + if isscalar(pars.margin) + pars.margin = pars.margin * [1 1 1 1]; + end + sz = sz - pars.margin(1:2) - pars.margin(3:4); + + % divide by columns + w = (sz(1) + pars.intercolumnspacing) / pars.cols - pars.intercolumnspacing; + sz(1) = w; + + % apply fill / aspect ratio + if pars.fill > 0 + % fill fraction + sz(2) = sz(2) * pars.fill; + elseif pars.fill < 0 + % aspect ratio + sz(2) = sz(1) * (-1 / pars.fill); + end + + % direct sizing model is used if either of width or height + % is set + if pars.width || pars.height + + % use aspect ratio to fill in either one that is not + % specified + if ~pars.width || ~pars.height + + % aspect ratio must not be a fill + if pars.fill >= 0 + error('cannot use fill fraction with direct sizing model'); + end + + % compute width + if ~pars.width + pars.width = pars.height * -pars.fill; + end + + % compute height + if ~pars.height + pars.height = pars.width / -pars.fill; + end + + end + + % store + sz = [pars.width pars.height]; + + end + +%%%% GET TARGET DIMENSIONS (END) + + + + % orientation of figure is upright, unless printing + % sideways, in which case the printing space is rotated too + if pars.sideways + set(p.h_figure, 'PaperOrientation', 'landscape') + sz = fliplr(sz); + else + set(p.h_figure, 'PaperOrientation', 'portrait') + end + + % report export size + msg = ['exporting to ' int2str(sz(1)) 'x' int2str(sz(2)) 'mm']; + if is_bitmap + psz = sz / 25.4 * pars.dpi; + msg = [msg ' (' int2str(psz(1)) 'x' int2str(psz(2)) 'px @ ' int2str(pars.dpi) 'DPI)']; + else + msg = [msg ' (vector format @ ' int2str(pars.dpi) 'DPI)']; + end + disp(msg); + + % if we are in defer state, we need to do a clean + % recompute first so that axes get positioned so that + % axis ticks get set correctly (if they are in + % automatic mode), since the LAYOUT_MODE_PREPRINT + % recompute will store the tick states. + if p.state.defer + p.state.defer = 0; + p.recomputeLayout([]); + end + + % turn off defer, if it is on + p.state.defer = 0; + + % do a pre-print layout + context.mode = panel.LAYOUT_MODE_PREPRINT; + context.size_in_mm = sz; + context.rect = [0 0 1 1]; + p.recomputeLayout(context); + + % need also to disable the warning that we should set + % PaperPositionMode to auto during this operation - + % we're setting it explicitly. + w = warning('off', 'MATLAB:Print:CustomResizeFcnInPrint'); + + % handle smoothing + pars.write_dpi = pars.dpi; + if pars.smooth > 1 + pars.write_dpi = pars.write_dpi * pars.smooth; + print_filename = [pars.filename '-temp']; + else + print_filename = pars.filename; + end + + % disable layout so it doesn't get computed during any + % figure resize operations that occur during printing. + p.state.defer = 1; + + % set size of figure now. it's important we do this + % after the pre-print layout, because in SVG export + % mode the on-screen figure size is changed and that + % would otherwise affect ticks and limits. + switch pars.fmt + + case 'svg' + % plot2svg (our current SVG export mechanism) uses + % 'Units' and 'Position' (i.e. on-screen position) + % rather than the Paper- prefixed ones used by the + % Matlab export functions. + + % store old on-screen position + svg_units = get(p.h_figure, 'Units'); + svg_pos = get(p.h_figure, 'Position'); + + % update on-screen position + set(p.h_figure, 'Units', 'centimeters'); + pos = get(p.h_figure, 'Position'); + pos(3:4) = sz / 10; + set(p.h_figure, 'Position', pos); + + otherwise + set(p.h_figure, ... + 'PaperUnits', 'centimeters', ... + 'PaperPosition', [0 0 sz] / 10, ... + 'PaperSize', sz / 10 ... % * 1.5 / 10 ... % CHANGED 21/06/2011 so that -opdf works correctly - why was this * 1.5, anyway? presumably was spurious... + ); + + end + + % do fixdash (not for SVG, since plot2svg does a nice + % job of dashed lines without our meddling...) + if ~isequal(pars.fmt, 'svg') + p.do_fixdash(context); + end + + % do the export + switch pars.fmt + case 'svg' + plot2svg(print_filename, p.h_figure); + otherwise + print(p.h_figure, '-loose', ['-d' pars.fmt], ['-r' int2str(pars.write_dpi)], print_filename) + end + + % undo fixdash + if ~isequal(pars.fmt, 'svg') + p.do_fixdash([]); + end + + % set on-screen figure size back to what it was, if it + % was changed. + switch pars.fmt + case 'svg' + set(p.h_figure, 'Units', svg_units); + set(p.h_figure, 'Position', svg_pos); + end + + % enable layout again (it was disabled, above, during + % printing). + p.state.defer = 0; + + % enable warnings + warning(w); + + % do a post-print layout + context.mode = panel.LAYOUT_MODE_POSTPRINT; + context.size_in_mm = []; + context.rect = [0 0 1 1]; + p.recomputeLayout(context); + + % handle smoothing + if pars.smooth > 1 + psz = sz * pars.smooth / 25.4 * pars.dpi; + msg = [' (reducing from ' int2str(psz(1)) 'x' int2str(psz(2)) 'px)']; + disp(['smoothing by factor ' int2str(pars.smooth) msg]); + im1 = imread(print_filename); + delete(print_filename); + sz = size(im1); + sz = [sz(1)-mod(sz(1),pars.smooth) sz(2)-mod(sz(2),pars.smooth)] / pars.smooth; + im = zeros(sz(1), sz(2), 3); + mm = 1:pars.smooth:(sz(1) * pars.smooth); + nn = 1:pars.smooth:(sz(2) * pars.smooth); + for m = 0:pars.smooth-1 + for n = 0:pars.smooth-1 + im = im + double(im1(mm+m, nn+n, :)); + end + end + im = uint8(im / (pars.smooth^2)); + + % set the DPI correctly in the new file + switch pars.fmt + case 'png' + dpm = pars.dpi / 25.4 * 1000; + imwrite(im, pars.filename, ... + 'XResolution', dpm, ... + 'YResolution', dpm, ... + 'ResolutionUnit', 'meter'); + case 'tiff' + imwrite(im, pars.filename, ... + 'Resolution', pars.dpi * [1 1]); + otherwise + imwrite(im, pars.filename); + end + end + + end + + function clearCallbacks(p) + + % clear all callback functions for the panel + % + % p.clearCallbacks() + p.m_callback = {}; + + end + + function setCallback(p, func, userdata) + + % set the callback function for the panel + % + % p.setCallback(myCallbackFunction, userdata) + % + % NB: this function clears all current callbacks, then + % calls addCallback(myCallbackFunction, userdata). + p.clearCallbacks(); + p.addCallback(func, userdata); + + end + + function addCallback(p, func, userdata) + + % attach a callback function to receive panel events + % + % p.addCallback(myCallbackFunction, userdata) + % register myCallbackFunction() to be called when + % events occur on the panel. at present, the only + % event is "layout-updated", which usually occurs + % after the figure is resized. myCallbackFunction() + % should accept one argument, "data", which will + % have the following fields. + % + % "userdata": the userdata passed to this function, if + % any was supplied, else empty. + % + % "panel": a reference to the panel on which the + % callback was set. this object can be queried in + % the usual way. + % + % "event": name of event (currently only + % "layout-updated"). + % + % "context": the layout context for the panel. this + % includes a field "size_in_mm" which is the + % physical size of the rendering surface (screen + % real estate, or image file) and "rect" which is + % the relative size of the rectangle within that + % occupied by the panel which is the context of + % the callback (in [left, bottom, width, height] + % format). + + invalid = ~isscalar(func) || ~isa(func, 'function_handle'); + if invalid + error('panel:InvalidArgument', 'argument to callback() must be a function handle'); + end + if nargin == 2 + p.m_callback{end+1} = {func []}; + else + p.m_callback{end+1} = {func userdata}; + end + + end + + function identify(p) + + % add annotations to help identify individual panels + % + % p.identify() + % when creating a complex layout, it can become + % confusing as to which panel is which. this + % function adds a text label to each axis panel + % indicating how to reference the axis panel through + % the root panel. for instance, if "(2, 3)" is + % indicated, you can find that panel at p(2, 3). + % + % see also: show() + + if p.isObject() + + % get managed axes + h_axes = p.getAllManagedAxes(); + + % if no axes, ignore + if isempty(h_axes) + return + end + + % mark first axis + h_axes = h_axes(1); + cla(h_axes); + text(0.5, 0.5, p.state.name, 'fontsize', 12, 'hori', 'center', 'parent', h_axes); + axis(h_axes, [0 1 0 1]); + grid(h_axes, 'off') + + else + + % recurse + for c = 1:length(p.m_children) + p.m_children(c).identify(); + end + + end + + end + + function repack(p, packspec) + + % change the packing specifier for an existing panel + % + % p.repack(packspec) + % repack() is a convenience function provided to + % allow easy development of a layout from the + % command prompt. packspec can be any packing + % specifier accepted by pack(). + % + % see also: pack() + + % deny repack() on root + if p.isRoot() + + % let's deny this. I'm not sure it makes anyway. you + % could always pack into root with a panel with + % absolute positioning, so let's deny first, and + % allow later if we're sure it's a good idea. + error('panel:InvalidArgument', 'root panel cannot be repack()ed'); + + end + + % validate + validate_packspec(packspec); + + % handle units + if iscell(packspec) + units = p.getPropertyValue('units'); + packspec{1} = panel.resolveUnits({packspec{1} units}, 'mm'); + end + + % update the packspec + p.packspec = packspec; + + % and recomputeLayout + p.recomputeLayout([]); + + end + + function pack(p, varargin) + + % add (pack) child panel(s) into an existing panel + % + % p.pack(...) + % add children to the panel "p", committing it as a + % "parent" panel (if it is not already). new (child) + % panels are created using this call - they start as + % "uncommitted" panels. if the parent already has + % children, the new children are appended. The + % following arguments are understood: + % + % 'h'/'v' - switch "p" to pack in the horizontal or + % vertical packing dimension for relative packing + % mode (default for new panels is vertical). + % + % {a, b, c, ...} (a cell row vector) - pack panels + % into "p" with "packing specifiers" a, b, c, etc. + % packing specifiers are detailed below. + % + % PACKING MODES + % panels can be packed into their parent in two + % modes, dependent on their packing specifier. you + % can see a visual representation of these modes on + % the HTML documentation page "Layout". + % + % (i) Relative Mode - panels are packed into the space + % occupied by their parent. size along the parent's + % "packing dimension" is dictated by the packing + % specifier; along the other dimension size matches + % the parent. the following packing specifiers + % indicate Relative Mode. + % + % a) Fixed Size: the specifier is a scalar double in + % a cell {d}. The panel will be of size d in the + % current units of "p" (see the property "p.units" + % for details, but default units are mm). + % + % b) Fractional Size: the specifier is a scalar + % double between 0 and 1 (or between 1 and 100, as a + % percentage). The panel is sized as a fraction of + % the space remaining in its parent after Fixed Size + % panels and inter-panel margins have been subtracted. + % + % c) Stretchable: the specifier is the empty matrix + % []. remaining space in the parent after Fixed and + % Fractional Size panels have been subtracted is + % shared out amongst Stretchable Size panels. + % + % (ii) Absolute Mode - panels hover above their + % parent and do not take up space, as if using the + % position:absolute property in CSS. The packing + % specifier is a 1x4 double vector indicating the + % [left bottom width height] of the panel in + % normalised coordinates of its parent. for example, + % the specifier [0 0 1 1] generates a child panel + % that fills its parent. + % + % SHORTCUTS + % + % ** a small scalar integer, N, (1 to 32) is expanded + % to {[], [], ... []}, with N entries. that is, it + % packs N panels in Relative Mode (Stretchable) and + % shares the available space between them. + % + % ** the call to pack() is recursive, so following a + % packing specifier list, an additional list will + % be used to generate a separate call to pack() on + % each of the children created by the first. hence: + % + % p.pack({[] []}, {[] []}) + % + % will create a 2x2 grid of panels that share the + % space of their parent, "p". since the argument + % "2" expands to {[] []} (see above), the same grid + % can be created using: + % + % p.pack(2, 2) + % + % which is a common idiom in the demos. NB: on + % recursion, the packing dimension is flipped + % automatically, so that a grid is formed. + % + % ** if no arguments are passed at all, a single + % argument {[]} is assumed, so that a single + % additional panel is packed into the parent in + % relative packing mode and with stretchable size. + % + % see also: panel (overview), panel/panel(), select() + % + % LEGACY + % + % the interface to pack() was changed at release + % 2.10 to add support for panels of fixed physical + % size. the interface offered at 2.9 and earlier is + % still available (look inside panel.m - search for + % text "LEGACY" - for details). + + % LEGACY + % + % releases of panel prior to 2.10 did not support + % panels of fixed physical size, and therefore had + % developed a different argument form to that used in + % 2.10 and beyond. specifically, the following + % additional arguments are accepted, for legacy + % support: + % + % 'abs' + % the next argument will be an absolute position, as + % described below. you should avoid using absolute + % positioning mode, in general, since this does not + % take advantage of panel's automatic layout. + % however, on occasion, you may need to take manual + % control of the position of one or more panels. see + % demopanelH for an example. + % + % 1xN row vector (without 'abs') + % pack N new panels along the packing dimension in + % relative mode, with the relative size of each + % given by the elements of the vector. -1 can be + % passed for any elements to mark those panel as + % 'stretchable', so that they fill available space + % left over by other panels packed alongside. the + % sum of the vector (apart from any -1 entries) + % should not come to more than 1, or a warning will + % be generated during laying out. an example would + % be [1/4 1/4 -1], to pack 3 panels, at 25, 25 and + % 50% relative sizes. though, NB, you can use + % percentages instead of fractions if you prefer, in + % which case they should not sum to over 100. so + % that same pack() would be [25 25 -1]. + % + % 1x4 row vector (after 'abs') + % pack 1 new panel using absolute positioning. the + % argument indicates the [left bottom width height] + % of the new panel, in normalised coordinates, as a + % fraction of its parent's position. panels using + % absolute positioning mode are ignored for the sake + % of layout, much like items using + % 'position:absolute' in CSS. + + % handle legacy, parse arguments from varargin into args + args = {}; + while ~isempty(varargin) + + % peel + arg = varargin{1}; + varargin = varargin(2:end); + + % handle shortcut (small integer) on current interface + if isa(arg, 'double') && isscalar(arg) && round(arg) == arg && arg >= 1 && arg <= 32 + arg = cell(1, arg); + end + + % handle current interface - note that the argument + % "recursive" is private and not advertised to the + % user. + if isequal(arg, 'h') || isequal(arg, 'v') || (iscell(arg) && isrow(arg)) || isequal(arg, 'recursive') + args{end+1} = arg; + continue + end + + % report (DEBUG) +% panel.debugmsg('use of LEGACY interface to pack()', 1); + + % handle legacy case + if isequal(arg, 'abs') + if length(varargin) ~= 1 || ~isnumeric(varargin{1}) || ~isofsize(varargin{1}, [1 4]) + error('panel:LegacyAbsNotFollowedBy1x4', 'the argument "abs" on the legacy interface should be followed by a [1x4] row vector'); + end + abs = varargin{1}; + varargin = varargin(2:end); + args{end+1} = {abs}; + continue + end + + % handle legacy case + if isa(arg, 'double') && isrow(arg) + arg_ = {}; + for a = 1:length(arg) + aa = arg(a); + if isequal(aa, -1) + arg_{end+1} = []; + else + arg_{end+1} = aa; + end + end + args{end+1} = arg_; + continue + end + + % unrecognised argument + error('panel:InvalidArgument', 'argument to pack() not recognised'); + + end + + % check m_panelType + if p.isObject() + error('panel:PackWhenObjectPanel', ... + 'cannot pack() into this panel - it is already committed as an object panel'); + end + + % if no arguments, simulate an argument of [], to pack + % a single panel of stretchable size + if isempty(args) + args = {{[]}}; + end + + % state + recursive = false; + + % handle arguments one by one + while ~isempty(args) && ischar(args{1}) + + % extract + arg = args{1}; + args = args(2:end); + + % handle string arguments + switch arg + case 'h' + p.packdim = 1; + case 'v' + p.packdim = 2; + case 'recursive' + recursive = true; + otherwise + error('panel:InvalidArgument', ['pack() did not recognise the argument "' arg '"']); + end + + end + + % if no more arguments that's weird but not bad + if isempty(args) + return + end + + % next argument now must be a cell + arg = args{1}; + args = args(2:end); + if ~iscell(arg) + panel.error('InternalError'); + end + + % commit as parent + p.commitAsParent(); + + % for each element + for i = 1:length(arg) + + % get packspec + packspec = arg{i}; + + % validate + validate_packspec(packspec); + + % handle units + if iscell(packspec) + units = p.getPropertyValue('units'); + packspec{1} = panel.resolveUnits({packspec{1} units}, 'mm'); + end + + % create a child + child = panel('pack', p); + child.packspec = packspec; + + % store it in the parent + if isempty(p.m_children) + p.m_children = child; + else + p.m_children(end+1) = child; + end + + % recurse (further argumens are passed on) + if ~isempty(args) + child_packdim = flippackdim(p.packdim); + edges = 'hv'; + child.pack('recursive', edges(child_packdim), args{:}); + end + + end + + % this must generate a recomputeLayout(), since the + % addition of new panels may affect the layout. any + % recursive call passes 'recursive', so that only the + % root call actually bothers doing a layout. + if ~recursive + p.recomputeLayout([]); + end + + end + + function h_out = select(p, h_object) + + % create or select an axis or object panel + % + % h = p.select(h) + % this call will return the handle of the object + % associated with the panel. if the panel is not yet + % committed, this will involve first committing it + % as an "object panel". if a list of objects ("h") + % is passed, these are the objects associated with + % the panel; if not, a new axis is created by the + % panel when this function is called. + % + % if the object list includes axes, then the "object + % panel" is also known as an "axis panel". in this + % case, the call to select() will make the (first) + % axis current, unless an output argument is + % requested, in which case the handle of the axes + % are returned but no axis is made current. + % + % the passed objects can be user-created axes (e.g. + % a colorbar) or any graphics object that is to have + % its position managed (e.g. a uipanel). your + % mileage may vary with different types of graphics + % object, please let me know. + % + % see also: panel (overview), panel/panel(), pack() + + % handle "all" and "data" + if nargin == 2 && isstring(h_object) && (strcmp(h_object, 'all') || strcmp(h_object, 'data')) + + % collect + h_out = []; + + % commit all uncommitted panels as axis panels by + % selecting them once + if p.isParent() + + % recurse + for c = 1:length(p.m_children) + h_out = [h_out p.m_children(c).select(h_object)]; + end + + elseif p.isUncommitted() + + % select in an axis + h_out = p.select(); + + % plot some data + if strcmp(h_object, 'data') + plot(h_out, randn(100, 1), 'k-'); + end + + end + + % ok + return + + end + + % check m_panelType + if p.isParent() + error('panel:SelectWhenParent', 'cannot select() this panel - it is already committed as a parent panel'); + end + + % commit as object + p.commitAsObject(); + + % assume not a new object + newObject = false; + + % use passed graphics object + if nargin >= 2 + + % validate + if ~all(ishandle(h_object)) + error('panel:InvalidArgument', 'argument to select() must be a list of handles to graphics objects'); + end + + % validate + if ~isempty(p.h_object) + error('panel:SelectWithObjectWhenObject', 'cannot select() new objects into this panel - it is already managing objects'); + end + + % store + p.h_object = h_object; + newObject = true; + + % make sure it has the correct parent - this doesn't + % seem to affect axes, so we do it for all + set(p.h_object, 'parent', p.h_parent); + + end + + % create new axis if necessary + if isempty(p.h_object) + % 'NextPlot', 'replacechildren' + % make sure fonts etc. don't get changed when user + % plots into it + p.h_object = axes( ... + 'Parent', p.h_parent, ... + 'NextPlot', 'replacechildren' ... + ); + newObject = true; + end + + % if wrapped objects include an axis, and no output args, make it current + h_axes = p.getAllManagedAxes(); + if ~isempty(h_axes) && ~nargout + set(p.h_figure, 'CurrentAxes', h_axes(1)); + + % 12/07/11: this call is slow, because it implies "drawnow" +% figure(p.h_figure); + + % 12/07/11: this call is fast, because it doesn't + set(0, 'CurrentFigure', p.h_figure); + + end + + % and return object list + if nargout + h_out = p.h_object; + end + + % this must generate a applyLayout(), since the axis + % will need positioning appropriately + if newObject + % 09/03/12 mitch + % if there isn't a context yet, we'll have to + % recomputeLayout(), in fact, to generate a context first. + % this will happen, for instance, if a single panel + % is generated in a window that was already open + % (no resize event will fire, and since pack() is + % not called, it will not call recomputeLayout() either). + % nonetheless, we have to reposition this object, so + % this forces us to recomputeLayout() now and generate + % that context we need. + if isempty(p.m_context) + p.recomputeLayout([]); + else + p.applyLayout(); + end + end + + end + + end + + + + + + + + + + + + + + %% ---- HIDDEN OVERLOADS ---- + + methods (Hidden = true) + + function out = vertcat(p, q) + error('panel2:MethodNotImplemented', 'concatenation is not supported by panel (use a cell array instead)'); + end + + function out = horzcat(p, q) + error('panel2:MethodNotImplemented', 'concatenation is not supported by panel (use a cell array instead)'); + end + + function out = cat(dim, p, q) + error('panel2:MethodNotImplemented', 'concatenation is not supported by panel (use a cell array instead)'); + end + + function out = ge(p, q) + error('panel2:MethodNotImplemented', 'inequality operators are not supported by panel'); + end + + function out = le(p, q) + error('panel2:MethodNotImplemented', 'inequality operators are not supported by panel'); + end + + function out = lt(p, q) + error('panel2:MethodNotImplemented', 'inequality operators are not supported by panel'); + end + + function out = gt(p, q) + error('panel2:MethodNotImplemented', 'inequality operators are not supported by panel'); + end + + function out = eq(p, q) + out = eq@handle(p, q); + end + + function out = ne(p, q) + out = ne@handle(p, q); + end + + end + + + + + + + + + + + + + + + + + + + + + + + + + + %% ---- PUBLIC HIDDEN GET/SET ---- + + methods (Hidden = true) + + function p = descend(p, indices) + + while ~isempty(indices) + + % validate + if numel(p) > 1 + error('panel:InvalidIndexing', 'you can only use () on a single (scalar) panel'); + end + + % validate + if ~p.isParent() + error('panel:InvalidIndexing', 'you can only use () on a parent panel'); + end + + % extract + index = indices{1}; + indices = indices(2:end); + + % only accept numeric + if ~isnumeric(index) || ~isscalar(index) + error('panel:InvalidIndexing', 'you can only use () with scalar indices'); + end + + % do the reference + p = p.m_children(index); + + end + + end + + function p_out = subsasgn(p, refs, value) + + % output is always subject + p_out = p; + + % handle () indexing + if strcmp(refs(1).type, '()') + p = p.descend(refs(1).subs); + refs = refs(2:end); + end + + % is that it? + if isempty(refs) + error('panel:InvalidIndexing', 'you cannot assign to a child panel'); + end + + % next ref must be . + if ~strcmp(refs(1).type, '.') + panel.error('InvalidIndexing'); + end + + % either one (.X) or two (.ch.X) + switch numel(refs) + + case 2 + + % validate + if ~strcmp(refs(2).type, '.') + panel.error('InvalidIndexing'); + end + + % validate + switch refs(2).subs + case {'fontname' 'fontsize' 'fontweight'} + case {'margin' 'marginleft' 'marginbottom' 'marginright' 'margintop'} + otherwise + panel.error('InvalidIndexing'); + end + + % avoid computing layout whilst setting descendant + % properties + p.defer(); + + % recurse + switch refs(1).subs + case {'children' 'ch'} + cs = p.m_children; + for c = 1:length(cs) + subsasgn(cs(c), refs(2:end), value); + end + case {'descendants' 'de'} + cs = p.getPanels('*'); + for c = 1:length(cs) + if cs{c} ~= p + subsasgn(cs{c}, refs(2:end), value); + end + end + case {'family' 'fa'} + cs = p.getPanels('*'); + for c = 1:length(cs) + subsasgn(cs{c}, refs(2:end), value); + end + end + + % release for laying out + p.undefer(); + + % mark for appropriate updates + refs(1).subs = refs(2).subs; + + case 1 + + % delegate + p.setPropertyValue(refs(1).subs, value); + + otherwise + panel.error('InvalidIndexing'); + + end + + % update layout as necessary + switch refs(1).subs + case {'fontname' 'fontsize' 'fontweight'} + p.applyLayout('recurse'); + case {'margin' 'marginleft' 'marginbottom' 'marginright' 'margintop'} + p.recomputeLayout([]); + end + + end + + function out = subsref(p, refs) + + % handle () indexing + if strcmp(refs(1).type, '()') + p = p.descend(refs(1).subs); + refs = refs(2:end); + end + + % is that it? + if isempty(refs) + out = p; + return + end + + % next ref must be . + if ~strcmp(refs(1).type, '.') + panel.error('InvalidIndexing'); + end + + % switch on "fieldname" + switch refs(1).subs + + case { ... + 'fontname' 'fontsize' 'fontweight' ... + 'margin' 'marginleft' ... + 'marginbottom' 'marginright' 'margintop' ... + 'units' ... + } + + % delegate this property get + out = p.getPropertyValue(refs(1).subs); + + case 'position' + out = p.getObjectPosition(); + + case 'figure' + out = p.h_figure; + + case 'packspec' + out = p.packspec; + + case 'axis' + if p.isObject() + out = p.getAllManagedAxes(); + else + out = []; + end + + case 'object' + if p.isObject() + h = p.h_object; + ih = ishandle(h); + out = h(ih); + else + out = []; + end + + case {'ch' 'children' 'de' 'descendants' 'fa' 'family'} + + % get the set + switch refs(1).subs + case {'children' 'ch'} + out = p.m_children; + case {'descendants' 'de'} + out = p.getPanels('*'); + for c = 1:length(out) + if out{c} == p + out = out([1:c-1 c+1:end]); + break + end + end + case {'family' 'fa'} + out = p.getPanels('*'); + end + + % we handle a special case of deeper reference + % here, because we are abusing matlab's syntax to + % do it. other cases (non-abusing) will be handled + % recursively, as usual. this is when we go: + % + % p.ch.axis + % + % which isn't syntactically sound since p.ch is a + % cell array (and potentially a non-singular one + % at that). we re-interpret this to mean + % [p.ch{1}.axis p.ch{2}.axis ...], as follows. + if length(refs) > 1 && isequal(refs(2).type, '.') + switch refs(2).subs + case {'axis' 'object'} + pp = out; + out = []; + for i = 1:length(pp) + out = cat(2, out, subsref(pp{i}, refs(2))); + end + refs = refs(2:end); % used up! + otherwise + % give an informative error message + panel.error('InvalidIndexing'); + end + end + + case { ... + 'addCallback' 'setCallback' 'clearCallbacks' ... + 'hold' ... + 'refresh' 'export' ... + 'pack' 'repack' ... + 'identify' 'show' ... + } + + % validate + if length(refs) ~= 2 || ~strcmp(refs(2).type, '()') + error('panel:InvalidIndexing', ['"' refs(1).subs '" is a function (try "help panel/' refs(1).subs '")']); + end + + % delegate this function call with no output + builtin('subsref', p, refs); + return + + case { ... + 'select' 'fixdash' ... + 'xlabel' 'ylabel' 'zlabel' 'title' ... + } + + % validate + if length(refs) ~= 2 || ~strcmp(refs(2).type, '()') + error('panel:InvalidIndexing', ['"' refs(1).subs '" is a function (try "help panel/' refs(1).subs '")']); + end + + % delegate this function call with output + if nargout + out = builtin('subsref', p, refs); + else + builtin('subsref', p, refs); + end + return + + otherwise + panel.error('InvalidIndexing'); + + end + + % continue + if length(refs) > 1 + out = subsref(out, refs(2:end)); + end + + end + + end + + + + + + + + + + + + + + %% ---- UTILITY METHODS ---- + + methods (Access = private) + + function b = ismanagefont(p) + + % ask root + b = p.m_root.state.manage_font; + + end + + function b = isdefer(p) + + % ask root + b = p.m_root.state.defer ~= 0; + + end + + function defer(p) + + % increment + p.m_root.state.defer = p.m_root.state.defer + 1; + + end + + function undefer(p) + + % decrement + p.m_root.state.defer = p.m_root.state.defer - 1; + + end + + function cs = getPanels(p, panelTypes, edgespec, all) + + % return all the panels that match the specification. + % + % panelTypes "*": return all panels + % panelTypes "s": return all sizeable panels (parent, + % object and uncommitted) + % panelTypes "p": return only physical panels (object + % or uncommitted) + % panelTypes "o": return only object panels + % + % if edgespec/all is specified, only panels matching + % the edgespec are returned (all of them if "all" is + % true, or any of them - the first one, in fact - if + % "all" is false). + + cs = {}; + + % do not include any that use absolute positioning - + % they stand outside of the sizing model + skip = (numel(p.packspec) == 4) && ~any(panelTypes == '*'); + + if p.isParent() + + % return if appropriate type + if any(panelTypes == '*s') && ~skip + cs = {p}; + end + + % if edgespec was supplied + if nargin == 4 + + % if we are perpendicular to the specified edge + if p.packdim ~= edgespec(1) + + if all + + % return all matching + for c = 1:length(p.m_children) + ppp = p.m_children(c).getPanels(panelTypes, edgespec, all); + cs = cat(2, cs, ppp); + end + + else + + % return only the first one + cs = cat(2, cs, p.m_children(1).getPanels(panelTypes, edgespec, all)); + + end + + else + + % if we are parallel to the specified edge + if edgespec(2) == 2 + + % use last + ppp = p.m_children(end).getPanels(panelTypes, edgespec, all); + cs = cat(2, cs, ppp); + + else + + % use first + cs = cat(2, cs, p.m_children(1).getPanels(panelTypes, edgespec, all)); + + end + + end + + else + + % else, return all + for c = 1:length(p.m_children) + ppp = p.m_children(c).getPanels(panelTypes); + cs = cat(2, cs, ppp); + end + + end + + elseif p.isObject() + + % return if appropriate type + if any(panelTypes == '*spo') && ~skip + cs = {p}; + end + + else + + % return if appropriate type + if any(panelTypes == '*sp') && ~skip + cs = {p}; + end + + end + + end + + function commitAsParent(p) + + if p.isUncommitted() + p.m_panelType = p.PANEL_TYPE_PARENT; + elseif p.isObject() + error('panel:AlreadyCommitted', 'cannot make this panel a parent panel, it is already an object panel'); + end + + end + + function commitAsObject(p) + + if p.isUncommitted() + p.m_panelType = p.PANEL_TYPE_OBJECT; + elseif p.isParent() + error('panel:AlreadyCommitted', 'cannot make this panel an object panel, it is already a parent panel'); + end + + end + + function b = isRoot(p) + + b = isempty(p.parent); + + end + + function b = isParent(p) + + b = p.m_panelType == p.PANEL_TYPE_PARENT; + + end + + function b = isObject(p) + + b = p.m_panelType == p.PANEL_TYPE_OBJECT; + + end + + function b = isUncommitted(p) + + b = p.m_panelType == p.PANEL_TYPE_UNCOMMITTED; + + end + + function h_axes = getAllManagedAxes(p) + + h_axes = []; + for n = 1:length(p.h_object) + h = p.h_object(n); + if isaxis(h) + h_axes = [h_axes h]; + end + end + + end + + function h_object = getOrCreateAxis(p) + + switch p.m_panelType + + case p.PANEL_TYPE_PARENT + + % create if not present + if isempty(p.h_object) + + % 'Visible', 'off' + % this is the hidden axis of a parent panel, + % used for displaying a parent panel's xlabel, + % ylabel and title, but not as a plotting axis + % + % 'NextPlot', 'replacechildren' + % make sure fonts etc. don't get changed when user + % plots into it + p.h_object = axes( ... + 'Parent', p.h_parent, ... + 'Visible', 'off', ... + 'NextPlot', 'replacechildren' ... + ); + + % make sure it's unitary, to help us in + % positioning labels and title + axis(p.h_object, [0 1 0 1]); + + % refresh this axis position + p.applyLayout(); + + end + + % ok + h_object = p.h_object; + + case p.PANEL_TYPE_OBJECT + + % ok + h_object = p.getAllManagedAxes(); + if isempty(h_object) + error('panel:ManagedObjectNotAnAxis', 'this object panel does not manage an axis'); + end + + case p.PANEL_TYPE_UNCOMMITTED + + panel.error('PanelUncommitted'); + + end + + end + + function removeChild(p, child) + + % if not a parent, fail but warn (shouldn't happen) + if ~p.isParent() + warning('panel:NotParentOnRemoveChild', 'i am not a parent (in removeChild())'); + return + end + + % remove from children + for c = 1:length(p.m_children) + if p.m_children(c) == child + p.m_children = p.m_children([1:c-1 c+1:end]); + return + end + end + + % warn + warning('panel:ChildAbsentOnRemoveChild', 'child not found (in removeChild())'); + + end + + function h = getShowAxis(p) + + if p.isRoot() + if isempty(p.h_showAxis) + + % create + p.h_showAxis = axes( ... + 'Parent', p.h_parent, ... + 'units', 'normalized', ... + 'position', [0 0 1 1] ... + ); + + % move to bottom + c = get(p.h_parent, 'children'); + c = [c(2:end); c(1)]; + set(p.h_parent, 'children', c); + + % finalise axis + set(p.h_showAxis, ... + 'xtick', [], 'ytick', [], ... + 'color', 'none', 'box', 'off' ... + ); + axis(p.h_showAxis, [0 1 0 1]); + + % hold + hold(p.h_showAxis, 'on'); + + end + + % return it + h = p.h_showAxis; + + else + h = p.parent.getShowAxis(); + end + + end + + function fireCallbacks(p, event) + + % for each attached callback + for c = 1:length(p.m_callback) + + % extract + callback = p.m_callback{c}; + func = callback{1}; + userdata = callback{2}; + + % fire + data = []; + data.panel = p; + data.event = event; + data.context = p.m_context; + data.userdata = userdata; + func(data); + + end + + end + + end + + + + + + + + + + + + + %% ---- LAYOUT METHODS ---- + + methods + + function refresh(p) + + % recompute layout of all panels + % + % p.refresh() + % recompute the layout of all panels from scratch. + % this should not usually be required, and is + % provided primarily for legacy support. + + % LEGACY + % + % NB: if you pass 'defer' to the constructor, calling + % refresh() both recomputes the layout and releases + % the defer mode. future changes to properties (e.g. + % margins) will cause immediate recomputation of the + % layout, so only call refresh() when you're done. + + % bubble up to root + if ~p.isRoot() + p.m_root.refresh(); + return + end + + % release defer + p.state.defer = 0; + + % debug output +% panel.debugmsg(['refresh "' p.state.name '"...']); + + % call recomputeLayout + p.recomputeLayout([]); + + end + + end + + methods (Access = private) + + function do_fixdash(p, context) + + % if context is [], this is _after_ the layout for + % export, so we need to restore + if isempty(context) + + % restore lines we changed to their original state + for r = 1:length(p.m_fixdash_restore) + + % get + restore = p.m_fixdash_restore{r}; + + % if empty, no change was made + if ~isempty(restore) + set(restore.h_line, ... + 'xdata', restore.xdata, 'ydata', restore.ydata); + delete([restore.h_supp restore.h_mark]); + end + + end + + else + +% % get handles to objects that still exist +% h_lines = p.m_fixdash(ishandle(p.m_fixdash)); + + % no restores + p.m_fixdash_restore = {}; + + % for each line + for i = 1:length(p.m_fixdash) + + % get + fix = p.m_fixdash{i}; + + % final check + if ~ishandle(fix.h) || ~isequal(get(fix.h, 'type'), 'line') + continue + end + + % apply dashstyle + p.m_fixdash_restore{end+1} = dashstyle_line(fix, context); + + end + + end + + end + + function p = recomputeLayout(p, context) + + % this function recomputes the layout from scratch. + % this means calculating the sizes of the root panel + % and all descendant panels. after this is completed, + % the function calls applyLayout to effect the new + % layout. + + % if not root, bubble up to root + if ~p.isRoot() + p.m_root.recomputeLayout(context); + return + end + + % if in defer mode, do not compute layout + if p.isdefer() + return + end + + % if no context supplied (e.g. on resize events), use + % the figure window (a context is supplied if + % exporting to an image file). + if isempty(context) + context.mode = panel.LAYOUT_MODE_NORMAL; + context.size_in_mm = []; + context.rect = [0 0 1 1]; + end + + % debug output +% panel.debugmsg(['recomputeLayout "' p.state.name '"...']); + +% % root may have a packspec of its own +% if ~isempty(p.packspec) +% if isscalar(p.packspec) +% % this should never happen, because it should be +% % caught when the packspec is set in repack() +% warning('panel:RootPanelCannotUseRelativeMode', 'the root panel uses relative positioning mode - this is ignored'); +% else +% context.rect = p.packspec; +% end +% end + + % if not given a context size, use the size on screen + % of the parent figure + if isempty(context.size_in_mm) + + % get context (whole parent) size in its units + pp = get(p.h_figure, 'position'); + context_size = pp(3:4); + + % defaults, in case this fails for any reason + screen_size = [1280 1024]; + if ismac + screen_dpi = 72; + else + screen_dpi = 96; + end + + % get screen DPI + try + local_screen_dpi = get(0, 'ScreenPixelsPerInch'); + if ~isempty(local_screen_dpi) + screen_dpi = local_screen_dpi; + end + end + + % get screen size + try + local_screen_size = get(0, 'ScreenSize'); + if ~isempty(local_screen_size) + screen_size = local_screen_size; + end + end + + % get figure width and height on screen + switch get(p.h_figure, 'Units') + + case 'points' + points_per_inch = 72; + context.size_in_mm = context_size / points_per_inch * 25.4; + + case 'inches' + context.size_in_mm = context_size * 25.4; + + case 'centimeters' + context.size_in_mm = context_size * 10.0; + + case 'pixels' + context.size_in_mm = context_size / screen_dpi * 25.4; + + case 'characters' + context_size = context_size .* [5 13]; % convert to pixels (based on empirical measurement) + context.size_in_mm = context_size / screen_dpi * 25.4; + + case 'normalized' + context_size = context_size .* screen_size(3:4); % convert to pixels (based on screen size) + context.size_in_mm = context_size / screen_dpi * 25.4; + + otherwise + error('panel:CaseNotCoded', ['case not coded, (Parent Units are ' get(p.h_figure, 'Units') ')']); + + end + + end + + % that's the figure size, now we need the size of our + % parent, if it's not the figure too + if p.h_parent ~= p.h_figure + units = get(p.h_parent, 'units'); + set(p.h_parent, 'units', 'normalized'); + pos = get(p.h_parent, 'position'); + set(p.h_parent, 'units', units); + context.size_in_mm = context.size_in_mm .* pos(3:4); + end + + % for the root, we apply the margins here, since it's + % a special case because there's always exactly one of + % it + margin = p.getPropertyValue('margin', 'mm'); + m = margin([1 3]) / context.size_in_mm(1); + context.rect = context.rect + [m(1) 0 -sum(m) 0]; + m = margin([2 4]) / context.size_in_mm(2); + context.rect = context.rect + [0 m(1) 0 -sum(m)]; + + % now, recurse + p.recurseComputeLayout(context); + + % clear h_showAxis when we recompute the layout + if ~isempty(p.h_showAxis) + delete(p.h_showAxis); + p.h_showAxis = []; + end + + % having computed the layout, we now apply it, + % starting at the root panel. + p.applyLayout('recurse'); + + end + + function recurseComputeLayout(p, context) + + % store context + p.m_context = context; + + % if no children, do nothing further + if isempty(p.m_children) + return + end + + % else, we're going to recompute the layout for our + % children + margins = []; + + % get size to pack into + mm_canvas = context.size_in_mm(p.packdim); + mm_context = mm_canvas * context.rect(2+p.packdim); + + % get list of children that are packed relative - we + % do this because the computation only handles these + % relative children; absolute packed children are + % ignored through the computation, and are just packed + % as specified when the time comes. + rel_list = []; + + % for each child + for i = 1:length(p.m_children) + + % get child + c = p.m_children(i); + + % is it packed abs? + if isofsize(c.packspec, [1 4]) + continue + end + + % if not, it's packed relative, so add to list + rel_list(end+1) = i; + + end + + % array of actual sizes as fraction of parent (note we + % only represent the rel_list). + zz = zeros(1, length(rel_list)); + sz_phys = zz; + sz_frac = zz; + i_stretch = zz; + + % for each child that is packed relative + for i = 1:length(rel_list) + + % get child + c = p.m_children(rel_list(i)); + + % get internal margin + margin = c.getPropertyValue('margin', 'mm'); + if p.packdim == 2 + margin = margin([2 4]); + margin = fliplr(margin); % doclink FLIP_PACKDIM_2 - same reason, here! + else + margin = margin([1 3]); + end + margins(i:i+1, i) = margin'; + + % subtract fixed size packspec from packing size + if iscell(c.packspec) + % NB: fixed size is always _stored_ in mm! + sz_phys(i) = c.packspec{1}; + end + + % get relative packing sizes + if isnumeric(c.packspec) && isscalar(c.packspec) + % NB: relative size is a scalar numeric + sz_frac(i) = c.packspec; + % convert perc to frac + if sz_frac(i) > 1 + sz_frac(i) = sz_frac(i) / 100; + end + end + + % get stretch packing size + if isempty(c.packspec) + % NB: these will be filled later + i_stretch(i) = 1; + end + + % else, it's an abs packing size, and we can ignore + % it for this phase of layout + + end + + % finalise internal margins (that is, the margin at + % each boundary between two adjacent relative packed + % panels is the maximum of the margins specified by + % each of the pair). + margins = max(margins, [], 2); + margins = margins(2:end-1)'; + + % subtract internal margins to give available space + % for objects (in mm) + mm_objects = mm_context - sum(margins); + + % now, subtract physically sized objects to give + % available space to share out amongst panels that + % specify their size as a fraction. + mm_share = mm_objects - sum(sz_phys); + + % and now stretch items can be given their actual + % fractional size, since we now know who they are + % sharing space with. + sz_frac(find(i_stretch)) = (1 - sum(sz_frac)) / sum(i_stretch); + + % and we can now get the real physical size of all the + % fractionally-sized panels in mm. + sz_frac = sz_frac * mm_share; + + % finally, we've got the physical boundaries of + % everything; let's just tidy that up. + sz = [[sz_phys + sz_frac]; margins 0]; + sz = sz(1:end-1); + + % and let's normalise the physical boundaries, because + % we're actually going to specify them to matlab in + % normalised form, even though we computed them in mm. + if ~isempty(sz) + + % do it + sz_norm = reshape([0 cumsum(sz / mm_context)]', 2, [])'; + + % for packdim 2, we pack from the top, whereas + % matlab's position property packs from the bottom, so + % we have to flip these. doclink FLIP_PACKDIM_2. + if p.packdim == 2 + sz_norm = fliplr(1 - sz_norm); + end + + end + + % recurse + for i = 1:length(p.m_children) + + % get child + c = p.m_children(i); + + % handle abs packed panels + if isofsize(c.packspec, [1 4]) + + % child context + child_context = context; + rect = child_context.rect; + rect([1 3]) = c.packspec([1 3]) * rect(3) + [rect(1) 0]; + rect([2 4]) = c.packspec([2 4]) * rect(4) + [rect(2) 0]; + child_context.rect = rect; + + else + + % child context + child_context = context; + rr = sz_norm(1, :); + sz_norm = sz_norm(2:end, :); % sz_norm has only as many entries as there are rel-packed panels + ri = p.packdim + [0 2]; + a = child_context.rect(ri(1)); + b = child_context.rect(ri(2)); + child_context.rect(ri) = [a+rr(1)*b diff(rr)*b]; + + end + + % recurse + c.recurseComputeLayout(child_context); + + end + + end + + function applyLayout(p, varargin) + + % this function applies the layout that is stored in + % each panel objects "m_context" member, and fixes up + % the position of any associated objects (such as axis + % group labels). + + % skip if disabled + if p.isdefer() + return + end + + % debug output +% panel.debugmsg(['applyLayout "' p.state.name '"...']); + + % defaults + recurse = false; + + % handle arguments + while ~isempty(varargin) + + % get + arg = varargin{1}; + varargin = varargin(2:end); + + % handle + switch arg + case 'recurse' + recurse = true; + otherwise + panel.error('InternalError'); + end + + end + + % recurse + if recurse + pp = p.getPanels('*'); + else + pp = {p}; + end + + % why do we have to split the applyLayout() operation + % into two? + % + % because the "group labels" are positioned with + % respect to the axes in their group depending on + % whether those axes have tick labels, and what those + % tick labels are. if those tick labels are in + % automatic mode (as they usually are), they may + % change when those axes are positioned. since an axis + % group may contain many of these nested deep, we have + % to position all axes (step 1) first, then (step 2) + % position any group labels. + + % step 1 + for pi = 1:length(pp) + pp{pi}.applyLayout1(); + end + + % step 2 + for pi = 1:length(pp) + pp{pi}.applyLayout2(); + end + + % callbacks + for pi = 1:length(pp) + fireCallbacks(pp{pi}, 'layout-updated'); + end + + end + + function r = getObjectPosition(p) + + % get packed position + r = p.m_context.rect; + + % if empty, must be absolute position + if isempty(r) + r = p.packspec; + pp = getObjectPosition(p.parent); + r = panel.getRectangleOfRectangle(pp, r); + end + + end + + function applyLayout1(p) + + % if no context yet, skip this call + if isempty(p.m_context) + return + end + + % if no managed objects, skip this call + if isempty(p.h_object) + return + end + + % debug output +% panel.debugmsg(['applyLayout1 "' p.state.name '"...']); + + % handle LAYOUT_MODE + switch p.m_context.mode + + case panel.LAYOUT_MODE_PREPRINT + + % if in LAYOUT_MODE_PREPRINT, store current axis + % layout (ticks and ticklabels) and lock them into + % manual mode so they don't get changed during the + % print operation + h_axes = p.getAllManagedAxes(); + for n = 1:length(h_axes) + p.state.store{n} = storeAxisState(h_axes(n)); + end + + case panel.LAYOUT_MODE_POSTPRINT + + % if in LAYOUT_MODE_POSTPRINT, restore axis + % layout, leaving it as it was before we ran + % export + h_axes = p.getAllManagedAxes(); + for n = 1:length(h_axes) + restoreAxisState(h_axes(n), p.state.store{n}); + end + + end + + % position it + try + set(p.h_object, 'position', p.getObjectPosition(), 'units', 'normalized'); + catch err + if strcmp(err.identifier, 'MATLAB:hg:set_chck:DimensionsOutsideRange') + w = warning('query', 'backtrace'); + warning off backtrace + warning('panel:PanelZeroSize', 'a panel had zero size, and the managed object was hidden'); + set(p.h_object, 'position', [-0.3 -0.3 0.2 0.2]); + if strcmp(w.state, 'on') + warning on backtrace + end + elseif strcmp(err.identifier, 'MATLAB:class:InvalidHandle') + % this will happen if the user deletes the managed + % objects manually. an obvious way that this + % happens is if the user select()s some panels so + % that axes get created, then calls clf. it would + % be nice if we could clear the panels attached to + % a figure in response to a clf call, but there + % doesn't seem any obvious way to pick up the clf + % call, only the delete(objects) that follows, and + % this is indistinguishable from a call by the + % user to delete(my_axis), for instance. how are + % we to respond if the user deletes the axis the + % panel is managing? it's not clear. so, we'll + % just fail silently, for now, and these panels + % will either never be used again (and will be + % destroyed when the figure is closed) or will be + % destroyed when the user creates a new panel on + % this figure. either way, i think, no real harm + % done. +% w = warning('query', 'backtrace'); +% warning off backtrace +% warning('panel:PanelObjectDestroyed', 'the object managed by a panel has been destroyed'); +% if strcmp(w.state, 'on') +% warning on backtrace +% end +% panel.debugmsg('***WARNING*** the object managed by a panel has been destroyed'); + return + else + rethrow(err) + end + end + + % if managing fonts + if p.ismanagefont() + + % apply properties to objects + h = p.h_object; + + % get those which are axes + h_axes = p.getAllManagedAxes(); + + % and labels/title objects, for any that are axes + for n = 1:length(h_axes) + h = [h ... + get(h_axes(n), 'xlabel') ... + get(h_axes(n), 'ylabel') ... + get(h_axes(n), 'zlabel') ... + get(h_axes(n), 'title') ... + ]; + end + + % apply font properties + set(h, ... + 'fontname', p.getPropertyValue('fontname'), ... + 'fontsize', p.getPropertyValue('fontsize'), ... + 'fontweight', p.getPropertyValue('fontweight') ... + ); + + end + + end + + function applyLayout2(p) + + % if no context yet, skip this call + if isempty(p.m_context) + return + end + + % if no object, skip this call + if isempty(p.h_object) + return + end + + % if not a parent, skip this call + if ~p.isParent() + return + end + + % if not an axis, skip this call - NB: this is not a + % displayed and managed object, rather it is the + % invisible axis used to display parent labels/titles. + % we checked above if this panel is a parent. thus, + % the member h_object must be scalar, if it is + % non-empty. + if ~isaxis(p.h_object) + return + end + + % debug output +% panel.debugmsg(['applyLayout2 "' p.state.name '"...']); + + % matlab moves x/ylabels around depending on + % whether the axis in question has any x/yticks, + % so that the label is always "near" the axis. + % we try to do the same, but it's hack-o-rama. + + % calibration offsets - i measured these + % empirically, what a load of shit + font_fudge = [2 1/3]; + nofont_fudge = [2 0]; + + % do xlabel + cs = p.getPanels('o', [2 2], true); + y = 0; + for c = 1:length(cs) + ch = cs{c}; + h_axes = ch.getAllManagedAxes(); + for h_axis = h_axes + % only if there are some tick labels, and they're + % at the bottom... + if ~isempty(get(h_axis, 'xticklabel')) && ~isempty(get(h_axis, 'xtick')) ... + && strcmp(get(h_axis, 'xaxislocation'), 'bottom') + fontoffset_mm = get(h_axis, 'fontsize') * font_fudge(2) + font_fudge(1); + y = max(y, fontoffset_mm); + end + end + end + y = max(y, get(p.h_object, 'fontsize') * nofont_fudge(2) + nofont_fudge(1)); + + % convert and lay in + axisheight_mm = p.m_context.size_in_mm(2) * p.m_context.rect(4); + y = y / axisheight_mm; + set(get(p.h_object, 'xlabel'), ... + 'VerticalAlignment', 'Cap', ... + 'Units', 'Normalized', ... + 'Position', [0.5 -y 1]); + + % calibration offsets - i measured these + % empirically, what a load of shit + font_fudge = [3 1/6]; + nofont_fudge = [2 0]; + + % do ylabel + cs = p.getPanels('o', [1 1], true); + x = 0; + for c = 1:length(cs) + ch = cs{c}; + h_axes = ch.getAllManagedAxes(); + for h_axis = h_axes + % only if there are some tick labels, and they're + % at the left... + if ~isempty(get(h_axis, 'yticklabel')) && ~isempty(get(h_axis, 'ytick')) ... + && strcmp(get(h_axis, 'yaxislocation'), 'left') + yt = get(h_axis, 'yticklabel'); + if ischar(yt) + ml = size(yt, 2); + else + ml = 0; + for i = 1:length(yt) + ml = max(ml, length(yt{i})); + end + end + fontoffset_mm = get(h_axis, 'fontsize') * ml * font_fudge(2) + font_fudge(1); + x = max(x, fontoffset_mm); + end + end + end + x = max(x, get(p.h_object, 'fontsize') * nofont_fudge(2) + nofont_fudge(1)); + + % convert and lay in + axisheight_mm = p.m_context.size_in_mm(1) * p.m_context.rect(3); + x = x / axisheight_mm; + set(get(p.h_object, 'ylabel'), ... + 'VerticalAlignment', 'Bottom', ... + 'Units', 'Normalized', ... + 'Position', [-x 0.5 1]); + + % calibration offsets - made up based on the + % ones i measured for the labels + nofont_fudge = [2 0]; + + % get y position + y = max(y, get(p.h_object, 'fontsize') * nofont_fudge(2) + nofont_fudge(1)); + + % convert and lay in + axisheight_mm = p.m_context.size_in_mm(2) * p.m_context.rect(4); + y = y / axisheight_mm; + set(get(p.h_object, 'title'), ... + 'VerticalAlignment', 'Bottom', ... + 'Position', [0.5 1+y 1]); + + end + + end + + + + + + + %% ---- PROPERTY METHODS ---- + + methods (Access = private) + + function value = getPropertyValue(p, key, units) + + value = p.prop.(key); + + if isempty(value) + + % inherit + if ~isempty(p.parent) + switch key + case {'fontname' 'fontsize' 'fontweight' 'margin' 'units'} + if nargin == 3 + value = p.parent.getPropertyValue(key, units); + else + value = p.parent.getPropertyValue(key); + end + return + end + end + + % default + if isempty(value) + value = panel.getPropertyDefault(key); + end + + end + + % translate dimensions + switch key + case {'margin'} + if nargin < 3 + units = p.getPropertyValue('units'); + end + value = panel.resolveUnits(value, units); + end + + end + + function setPropertyValue(p, key, value) + + % root properties + switch key + case 'units' + if ~isempty(p.parent) + p.parent.setPropertyValue(key, value); + return + end + end + + % value validation + switch key + case 'units' + invalid = ~( (isstring(value) && isin({'mm' 'in' 'cm' 'pt'}, value)) || isempty(value) ); + case 'fontname' + invalid = ~( isstring(value) || isempty(value) ); + case 'fontsize' + invalid = ~( (isnumeric(value) && isscalar(value) && value >= 4 && value <= 60) || isempty(value) ); + case 'fontweight' + invalid = ~( (isstring(value) && isin({'normal' 'bold'}, value)) || isempty(value) ); + case 'margin' + invalid = ~( (isdimension(value)) || isempty(value) ); + case {'marginleft' 'marginbottom' 'marginright' 'margintop'} + invalid = ~isscalardimension(value); + otherwise + error('panel:UnrecognisedProperty', ['unrecognised property "' key '"']); + end + + % value validation + if invalid + error('panel:InvalidValueForProperty', ['invalid value for property "' key '"']); + end + + % marginX properties + switch key + case {'marginleft' 'marginbottom' 'marginright' 'margintop'} + index = isin({'left' 'bottom' 'right' 'top'}, key(7:end)); + element = value; + value = p.getPropertyValue('margin'); + value(index) = element; + key = 'margin'; + end + + % translate dimensions + switch key + case {'margin'} + if isscalar(value) + value = value * [1 1 1 1]; + end + if ~isempty(value) + units = p.getPropertyValue('units'); + value = {panel.resolveUnits({value units}, 'mm') 'mm'}; + end + end + + % lay in + p.prop.(key) = value; + + end + + end + + methods (Static = true, Access = private) + + function s = fignum(h) + + % handled differently pre/post 2014b + if isa(h, 'matlab.ui.Figure') + % R2014b + s = num2str(h.Number); + else + % pre-R2014b + s = num2str(h); + end + end + + function prop = getPropertyInitialState() + + prop = panel.getPropertyDefaults(); + for key = fieldnames(prop)' + prop.(key{1}) = []; + end + + end + + function value = getPropertyDefault(key) + + persistent defprop + + if isempty(defprop) + defprop = panel.getPropertyDefaults(); + end + + value = defprop.(key); + + end + + function defprop = getPropertyDefaults() + + % root properties + defprop.units = 'mm'; + + % inherited properties + defprop.fontname = get(0, 'defaultAxesFontName'); + defprop.fontsize = get(0, 'defaultAxesFontSize'); + defprop.fontweight = 'normal'; + defprop.margin = {[15 15 5 5] 'mm'}; + + % not inherited properties + % CURRENTLY, NONE! +% defprop.align = false; + + end + + end + + + + + + + + %% ---- STATIC PUBLIC METHODS ---- + + methods (Static = true) + + function p = recover(h_figure) + + % get a handle to the root panel associated with a figure + % + % p = recover(h_fig) + % if you have not got a handle to the root panel of + % the figure h_fig, this call will retrieve it. if + % h_fig is not supplied, gcf is used. + + if nargin < 1 + h_figure = gcf; + end + + p = panel.callbackDispatcher('recover', h_figure); + + end + + function version() + + % report the version of panel that is active + % + % panel.version() + + fid = fopen(which('panel')); + tag = '% Release Version'; + ltag = length(tag); + tagline = 'Unable to determine Release Version'; + while 1 + line = fgetl(fid); + if ~ischar(line) + break + end + if length(line) > ltag && strcmp(line(1:ltag), tag) + tagline = line(3:end); + end + end + fclose(fid); + disp(tagline) + + end + + function panic() + + % call delete on all children of the global workspace, + % to recover from bugs that leave us with uncloseable + % figures. call this as "panel.panic()". + % + % NB: if you have to call panic(), something has gone + % wrong. if you are able to reproduce the problem, + % please contact me to report the bug. + delete(allchild(0)); + + end + + end + + + + + + + %% ---- STATIC PRIVATE METHODS ---- + + methods (Static = true, Access = private) + + function error(id) + + switch id + case 'PanelUncommitted' + throwAsCaller(MException('panel:PanelUncommitted', 'this action cannot be performed on an uncommitted panel')); + case 'InvalidIndexing' + throwAsCaller(MException('panel:InvalidIndexing', 'you cannot index a panel object in this way')); + case 'InternalError' + throwAsCaller(MException('panel:InternalError', 'an internal error occurred')); + otherwise + throwAsCaller(MException('panel:UnknownError', ['an unknown error was generated with id "' id '"'])); + end + + end + + function lockClass() + + persistent hasLocked + if isempty(hasLocked) + + % only lock if not in debug mode + if ~panel.isDebug() + % in production code, must mlock() file at this point, + % to avoid persistent variables being cleared by user + if strcmp(getenv('USERDOMAIN'), 'BERGEN') + % my machine, do nothing + else + mlock + end + end + + % mark that we've handled this + hasLocked = true; + + end + + end + + function debugmsg(msg, focus) + + % focus can be supplied to force only focussed + % messages to be shown + if nargin < 2 + focus = 1; + end + + % display, if in debug mode + if focus + if panel.isDebug() + disp(msg); + end + end + + end + + function state = isDebug() + + % persistent + persistent debug + + % create + if isempty(debug) + try + debug = panel_debug_state(); + catch + debug = false; + end + end + + % ok + state = debug; + + end + + function r = getFractionOfRectangle(r, dim, range) + + switch dim + case 1 + r = [r(1)+range(1)*r(3) r(2) range(2)*r(3) r(4)]; + case 2 + r = [r(1) r(2)+(1-sum(range))*r(4) r(3) range(2)*r(4)]; + otherwise + error('panel:CaseNotCoded', ['case not coded, dim = ' dim ' (internal error)']); + end + + end + + function r = getRectangleOfRectangle(r, s) + + w = r(3); + h = r(4); + r = [r(1)+s(1)*w r(2)+s(2)*h s(3)*w s(4)*h]; + + end + + function a = getUnionRect(a, b) + + if isempty(a) + a = b; + end + if ~isempty(b) + d = a(1) - b(1); + if d > 0 + a(1) = a(1) - d; + a(3) = a(3) + d; + end + d = a(2) - b(2); + if d > 0 + a(2) = a(2) - d; + a(4) = a(4) + d; + end + d = b(1) + b(3) - (a(1) + a(3)); + if d > 0 + a(3) = a(3) + d; + end + d = b(2) + b(4) - (a(2) + a(4)); + if d > 0 + a(4) = a(4) + d; + end + end + + end + + function r = reduceRectangle(r, margin) + + r(1:2) = r(1:2) + margin(1:2); + r(3:4) = r(3:4) - margin(1:2) - margin(3:4); + + end + + function v = normaliseDimension(v, space_size_in_mm) + + v = v ./ [space_size_in_mm space_size_in_mm]; + + end + + function v = resolveUnits(d, units) + + % first, convert into mm + v = d{1}; + switch d{2} + case 'mm' + % ok + case 'cm' + v = v * 10.0; + case 'in' + v = v * 25.4; + case 'pt' + v = v / 72.0 * 25.4; + otherwise + error('panel:CaseNotCoded', ['case not coded, storage units = ' units ' (internal error)']); + end + + % then, convert to specified units + switch units + case 'mm' + % ok + case 'cm' + v = v / 10.0; + case 'in' + v = v / 25.4; + case 'pt' + v = v / 25.4 * 72.0; + otherwise + error('panel:CaseNotCoded', ['case not coded, requested units = ' units ' (internal error)']); + end + + end + + function resizeCallback(obj, evt) + + panel.callbackDispatcher('resize', obj); + + end + + function closeCallback(obj, evt) + + panel.callbackDispatcher('delete', obj); + delete(obj); + + end + + function out = callbackDispatcher(op, data) + + % debug output +% panel.debugmsg(['callbackDispatcher(' op ')...']) + + % persistent store + persistent registeredPanels + + % switch on operation + switch op + + case {'register' 'registerNoClear'} + + % if a root panel is already attached to this + % figure, we could throw an error and refuse to + % create the new object, we could delete the + % existing panel, or we could allow multiple + % panels to be attached to the same figure. + % + % we should allow multiple panels, because they + % may have different parents within the same + % figure (e.g. uipanels). but by default we don't, + % unless the panel.add() static constructor is + % used. + + if strcmp(op, 'register') + + argument_h_figure = data.h_figure; + i = 0; + while i < length(registeredPanels) + i = i + 1; + if registeredPanels(i).h_figure == argument_h_figure + delete(registeredPanels(i)); + i = 0; + end + end + + end + + % register the new panel + if isempty(registeredPanels) + registeredPanels = data; + else + registeredPanels(end+1) = data; + end + + % debug output +% panel.debugmsg(['panel registered (' int2str(length(registeredPanels)) ' now registered)']); + + case 'unregister' + + % debug output +% panel.debugmsg(['on unregister, ' int2str(length(registeredPanels)) ' registered']); + + for r = 1:length(registeredPanels) + if registeredPanels(r) == data + registeredPanels = registeredPanels([1:r-1 r+1:end]); + + % debug output +% panel.debugmsg(['panel unregistered (' int2str(length(registeredPanels)) ' now registered)']); + + return + end + end + + % warn + warning('panel:AbsentOnCallbacksUnregister', 'panel was absent from the callbacks register when it tried to unregister itself'); + + case 'resize' + + argument_h_parent = data; + for r = 1:length(registeredPanels) + if registeredPanels(r).h_parent == argument_h_parent + registeredPanels(r).recomputeLayout([]); + end + end + + case 'recover' + + argument_h_figure = data; + out = []; + for r = 1:length(registeredPanels) + if registeredPanels(r).h_figure == argument_h_figure + if isempty(out) + out = registeredPanels(r); + else + out(end+1) = registeredPanels(r); + end + end + end + + case 'delete' + + argument_h_figure = data; + i = 0; + while i < length(registeredPanels) + i = i + 1; + if registeredPanels(i).h_figure == argument_h_figure + delete(registeredPanels(i)); + i = 0; + end + end + + end + + end + + end + + + + +end + + + + + + + + + + + + + + + + + +% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% +% HELPERS +% +% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +function restore = dashstyle_line(fix, context) + +% get axis size in mm +h_line = fix.h; +h_axis = get(h_line, 'parent'); +u = get(h_axis, 'units'); +set(h_axis, 'units', 'norm'); +pos = get(h_axis, 'position'); +set(h_axis, 'units', u); +axis_in_mm = pos(3:4) .* context.size_in_mm; + +% recover data +xdata = get(h_line, 'xdata'); +ydata = get(h_line, 'ydata'); +zdata = get(h_line, 'zdata'); +linestyle = get(h_line, 'linestyle'); +marker = get(h_line, 'marker'); + +% empty restore +restore = []; + +% do not handle 3D +if ~isempty(zdata) + warning('panel:NoFixdash3D', 'panel cannot fixdash() a 3D line - no action taken'); + return +end + +% get range of axis +ax = axis(h_axis); + +% get scale in each dimension (mm per unit) +sc = axis_in_mm ./ (ax([2 4]) - ax([1 3])); + +% create empty line +data = NaN; + +% override linestyle +if ~isempty(fix.linestyle) + linestyle = fix.linestyle; +end + +% transcribe linestyle +linestyle = dashstyle_parse_linestyle(linestyle); +if isempty(linestyle) + return +end + +% scale +scale = 1; +dashes = linestyle * scale; + +% store for restore +restore.h_line = h_line; +restore.xdata = xdata; +restore.ydata = ydata; + +% create another, separate, line to overlay on the original +% line and render the fixed-up dashes. +restore.h_supp = copyobj(h_line, h_axis); + +% if the original line has markers, we'll have to create yet +% another separate line instance to represent them, because +% they shouldn't be "dashed", as it were. note that we don't +% currently attempt to get the z-order right for these +% new lines. +if ~isequal(marker, 'none') + restore.h_mark = copyobj(h_line, h_axis); + set(restore.h_mark, 'linestyle', 'none'); + set(restore.h_supp, 'marker', 'none'); +else + restore.h_mark = []; +end + +% hide the original line. this line remains in existence so +% that if there is a legend, it doesn't get messed up. +set(h_line, 'xdata', NaN, 'ydata', NaN); + +% extract pattern length +patlen = sum(dashes); + +% position within pattern is initially zero +pos = 0; + +% linedata +line_xy = complex(xdata, ydata); + +% for each line segment +while length(line_xy) > 1 + + % get line segment + xy = line_xy(1:2); + line_xy = line_xy(2:end); + + % any NaNs, and we're outta here + if any(isnan(xy)) + continue + end + + % get start etc. + O = xy(1); + V = xy(2) - xy(1); + + % get mm length of this line segment + d = sqrt(sum(([real(V) imag(V)] .* sc) .^ 2)); + + % and mm unit vector + U = V / d; + + % generate a long-enough pattern for this segment + n = ceil((pos + d) / patlen); + pat = [0 cumsum(repmat(dashes, [1 n]))] - pos; + pos = d - (pat(end) - patlen); + pat = [pat(1:2:end-1); pat(2:2:end)]; + + % trim spurious segments + pat = pat(:, any(pat >= 0) & any(pat <= d)); + + % skip if that's it + if isempty(pat) + continue + end + + % and reduce ones that are oversized + pat(1) = max(pat(1), 0); + pat(end) = min(pat(end), d); + + % finally, add these segments to the line data + seg = [O + pat * U; NaN(1, size(pat, 2))]; + data = [data seg(:).']; + +end + +% update line +set(restore.h_supp, 'xdata', real(data), 'ydata', imag(data), ... + 'linestyle', '-'); + +end + + +function linestyle = dashstyle_parse_linestyle(linestyle) + +if isequal(linestyle, 'none') || isequal(linestyle, '-') + linestyle = []; + return +end + +while 1 + + % if numbers + if isnumeric(linestyle) + if ~isa(linestyle, 'double') || ~isrow(linestyle) || mod(length(linestyle), 2) ~= 0 + break + end + % no need to parse + return + end + + % else, must be char + if ~ischar(linestyle) || ~isrow(linestyle) + break + end + + % translate matlab non-standard codes into codes we can + % easily parse + switch linestyle + case ':' + linestyle = '.'; + case '--' + linestyle = '-'; + end + + % must be only - and . + if any(linestyle ~= '.' & linestyle ~= '-') + break + end + + % transcribe + c = linestyle; + linestyle = []; + for l = c + switch l + case '-' + linestyle = [linestyle 2 0.75]; + case '.' + linestyle = [linestyle 0.5 0.75]; + end + end + return + +end + +warning('panel:BadFixdashLinestyle', 'unusable linestyle in fixdash()'); +linestyle = []; + +end + + + +% MISCELLANEOUS + +function index = isin(list, value) + +for i = 1:length(list) + if strcmp(value, list{i}) + index = i; + return + end +end + +index = 0; + +end + +function dim = flippackdim(dim) + +% this function, used between arguments in a recursive call, +% causes the dim to be switched with each recurse, so that +% we build a grid, rather than a long, long row. +dim = 3 - dim; + +end + + + +% STRING PADDING FUNCTIONS + +function s = rpad(s, l) + +if nargin < 2 + l = 16; +end + +if length(s) < l + s = [s repmat(' ', 1, l - length(s))]; +end + +end + +function s = lpad(s, l) + +if nargin < 2 + l = 16; +end + +if length(s) < l + s = [repmat(' ', 1, l - length(s)) s]; +end + +end + + + +% HANDLE GRAPHICS HELPERS + +function h = getParentFigure(h) + +if strcmp(get(h, 'type'), 'figure') + return +else + h = getParentFigure(get(h, 'parent')); +end + +end + +function addHandleCallback(h, name, func) + +% % get current list of callbacks +% callbacks = get(h, name); +% +% % if empty, turn into a cell +% if isempty(callbacks) +% callbacks = {}; +% elseif iscell(callbacks) +% % only add ourselves once +% for c = 1:length(callbacks) +% if callbacks{c} == func +% return +% end +% end +% else +% callbacks = {callbacks}; +% end +% +% % and add ours (this is friendly, in case someone else has a +% % callback attached) +% callbacks{end+1} = func; +% +% % lay in +% set(h, name, callbacks); + +% the above isn't as simple as i thought - for now, we'll +% just stamp on any existing callbacks +set(h, name, func); + +end + +function store = storeAxisState(h) + +% LOCK TICKS AND LIMITS +% +% (LOCK TICKS) +% +% lock state so that the ticks and labels do not change when +% the figure is resized for printing. this is what the user +% will expect, which is why we go through this palaver. +% +% however, for fuck's sake. the following code illustrates +% an idiosyncrasy of matlab (i would call this an +% inconsistency, myself, but there you go). +% +% figure +% axis([0 1 0 1]) +% set(gca, 'ytick', [-1 0 1 2]) +% get(gca, 'yticklabel') +% set(gca, 'yticklabelmode', 'manual') +% +% now, resize the figure window. at least in R2011b, the +% tick labels change on the first resize event. presumably, +% this is because matlab treats the ticklabel value +% differently depending on if the ticklabelmode is auto or +% manual. if it's manual, the value is used as documented, +% and [0 1] is used to label [-1 0 1 2], cyclically. +% however, if the ticklabelmode is auto, and the ticks +% extend outside the figure, then the ticklabels are set +% sensibly, but the _value_ of ticklabel is not consistent +% with what it would need to be to get this tick labelling +% were the mode manual. and, in a final bizarre twist, this +% doesn't become evident until the resize event. i think +% this is a bug, no other way of looking at it; at best it's +% an inconsistency that is either tedious or impossible to +% work around in the general case. +% +% in any case, we have to lock the ticks to manual as we go +% through the print cycle, so that the ticks do not get +% changed if they were in automatic mode. but we mustn't fix +% the tick labels to manual, since if we do we may encounter +% this inconsistency and end up with the wrong tick labels +% in the print out. i can't, at time of writing, think of a +% case where we'd have to fix the tick labels to manual too. +% the possible cases are: +% +% ticks auto, labels auto: in this case, fixing the ticks to +% manual should be enough. +% +% ticks manual, labels auto: leave as is. +% +% ticks manual, labels manual: leave as is. +% +% the only other case is ticks auto, labels manual, which is +% a risky case to use, but in any case we can also fix the +% ticks to manual in that case. thus, our preferred solution +% is to always switch the ticks to manual, if they're not +% already, and otherwise leave things be. +% +% (LOCK LIMITS) +% +% the other thing that may get modified, if the user hasn't +% fixed it, is the axis limits. so we lock them too, any +% that are set to auto, and mark them for unlocking when the +% print is complete. + +store = ''; + +% manual-ise ticks on any axis where they are currently +% automatic, and indicate that we need to switch them back +% afterwards. +if strcmp(get(h, 'XTickMode'), 'auto') + store = [store 'X']; + set(h, 'XTickMode', 'manual'); +end +if strcmp(get(h, 'YTickMode'), 'auto') + store = [store 'Y']; + set(h, 'YTickMode', 'manual'); +end +if strcmp(get(h, 'ZTickMode'), 'auto') + store = [store 'Z']; + set(h, 'ZTickMode', 'manual'); +end + +% manual-ise limits on any axis where they are currently +% automatic, and indicate that we need to switch them back +% afterwards. +if strcmp(get(h, 'XLimMode'), 'auto') + store = [store 'x']; + set(h, 'XLimMode', 'manual'); +end +if strcmp(get(h, 'YLimMode'), 'auto') + store = [store 'y']; + set(h, 'YLimMode', 'manual'); +end +if strcmp(get(h, 'ZLimMode'), 'auto') + store = [store 'z']; + set(h, 'ZLimMode', 'manual'); +end + +% % OLD CODE OBSOLETED 25/01/12 - see notes above +% +% % store current state +% store.XTick = get(h, 'XTick'); +% store.XTickMode = get(h, 'XTickMode'); +% store.XTickLabel = get(h, 'XTickLabel'); +% store.XTickLabelMode = get(h, 'XTickLabelMode'); +% store.YTickMode = get(h, 'YTickMode'); +% store.YTick = get(h, 'YTick'); +% store.YTickLabel = get(h, 'YTickLabel'); +% store.YTickLabelMode = get(h, 'YTickLabelMode'); +% store.ZTick = get(h, 'ZTick'); +% store.ZTickMode = get(h, 'ZTickMode'); +% store.ZTickLabel = get(h, 'ZTickLabel'); +% store.ZTickLabelMode = get(h, 'ZTickLabelMode'); +% +% % lock state to manual +% set(h, 'XTickLabelMode', 'manual'); +% set(h, 'XTickMode', 'manual'); +% set(h, 'YTickLabelMode', 'manual'); +% set(h, 'YTickMode', 'manual'); +% set(h, 'ZTickLabelMode', 'manual'); +% set(h, 'ZTickMode', 'manual'); + +end + +function restoreAxisState(h, store) + +% unmanualise +for item = store + switch item + case {'X' 'Y' 'Z'} + set(h, [item 'TickMode'], 'auto'); + case {'x' 'y' 'z'} + set(h, [upper(item) 'TickMode'], 'auto'); + end +end + +% % OLD CODE OBSOLETED 25/01/12 - see notes above +% +% % restore passed state +% set(h, 'XTick', store.XTick); +% set(h, 'XTickMode', store.XTickMode); +% set(h, 'XTickLabel', store.XTickLabel); +% set(h, 'XTickLabelMode', store.XTickLabelMode); +% set(h, 'YTick', store.YTick); +% set(h, 'YTickMode', store.YTickMode); +% set(h, 'YTickLabel', store.YTickLabel); +% set(h, 'YTickLabelMode', store.YTickLabelMode); +% set(h, 'ZTick', store.ZTick); +% set(h, 'ZTickMode', store.ZTickMode); +% set(h, 'ZTickLabel', store.ZTickLabel); +% set(h, 'ZTickLabelMode', store.ZTickLabelMode); + +end + + + +% DIM AND EDGE HANDLING + +% we describe each edge of a panel in terms of "dim" (1 or +% 2, horizontal or vertical) and "edge" (1 or 2, former or +% latter). together, [dim edge] is an "edgespec". + +function s = edgestr(edgespec) + +s = 'lbrt'; +s = s(edgeindex(edgespec)); + +end + +function i = edgeindex(edgespec) + +% edge indices. margins are stored as [l b r t] but +% dims are packed left to right and top to bottom, so +% relationship between 'dim' and 'end' and index into +% margin is non-trivial. we call the index into the margin +% the "edgeindex". an "edgespec" is just [dim end], in a +% single array. +i = [1 3; 4 2]; +i = i(edgespec(1), edgespec(2)); + +end + + + +% VARIABLE TYPE HELPERS + +function val = validate_par(val, argtext, varargin) + +% this helper validates arguments to some functions in the +% main body + +for n = 1:length(varargin) + + % get validation constraint + arg = varargin{n}; + + % handle string list + if iscell(arg) + % string list + if ~isin(arg, val) + error('panel:InvalidArgument', ... + ['invalid argument "' argtext '", "' val '" is not a recognised data value for this option']); + end + continue; + end + + % handle strings + if isstring(arg) + switch arg + case 'empty' + if ~isempty(val) + error('panel:InvalidArgument', ... + ['invalid argument "' argtext '", option does not expect any data']); + end + case 'dimension' + if ~isdimension(val) + error('panel:InvalidArgument', ... + ['invalid argument "' argtext '", option expects a dimension']); + end + case 'scalar' + if ~(isnumeric(val) && isscalar(val) && ~isnan(val)) + error('panel:InvalidArgument', ... + ['invalid argument "' argtext '", option expects a scalar value']); + end + case 'nonneg' + if any(val(:) < 0) + error('panel:InvalidArgument', ... + ['invalid argument "' argtext '", option expects non-negative values only']); + end + case 'integer' + if any(val(:) ~= round(val(:))) + error('panel:InvalidArgument', ... + ['invalid argument "' argtext '", option expects integer values only']); + end + end + continue; + end + + % handle numeric range + if isnumeric(arg) && isofsize(arg, [1 2]) + if any(val(:) < arg(1)) || any(val(:) > arg(2)) + error('panel:InvalidArgument', ... + ['invalid argument "' argtext '", option data must be between ' num2str(arg(1)) ' and ' num2str(arg(2))]); + end + continue; + end + + % not recognised + arg + error('panel:InternalError', 'internal error - bad argument to validate_par (above)'); + +end + +end + +function b = checkpar(value, mn, mx) + +b = isscalar(value) && isnumeric(value) && ~isnan(value); +if b + if nargin >= 2 + b = b && value >= mn; + end + if nargin >= 3 + b = b && value <= mx; + end +end + +end + +function b = isintegral(v) + +b = all(all(v == round(v))); + +end + +function b = isstring(value) + +sz = size(value); +b = ischar(value) && length(sz) == 2 && sz(1) == 1 && sz(2) >= 1; + +end + +function b = isdimension(value) + +b = isa(value, 'double') && (isscalar(value) || isofsize(value, [1 4])); + +end + +function b = isscalardimension(value) + +b = isa(value, 'double') && isscalar(value); + +end + +function b = isofsize(value, siz) + +sz = size(value); +b = length(sz) == length(siz) && all(sz == siz); + +end + +function b = isaxis(h) + +b = ishandle(h) && strcmp(get(h, 'type'), 'axes'); + +end + +function validate_packspec(packspec) + + % stretchable + if isempty(packspec) + return + end + + % scalar + if isa(packspec, 'double') && isscalar(packspec) + + % fraction + if packspec > 0 && packspec <= 1 + return + end + + % percentage + if packspec > 1 && packspec <= 100 + return + end + + end + + % fixed + if iscell(packspec) && isscalar(packspec) + + % delve + d = packspec{1}; + if isa(d, 'double') && isscalar(d) && d > 0 + return + end + + end + + % abs + if isa(packspec, 'double') && isofsize(packspec, [1 4]) && all(packspec(3:4)>0) + return + end + + % otherwise, bad form + error('panel:BadPackingSpecifier', 'the packing specifier was not valid - see help panel/pack'); + +end + + + diff --git a/Plots.m b/Plots.m deleted file mode 100644 index dfcf2da..0000000 --- a/Plots.m +++ /dev/null @@ -1,100 +0,0 @@ -% mex command is given by: -% mex CXXFLAGS="\$CXXFLAGS -std=c++11 -O3" TC_mex.cpp Cortical_Column.cpp Thalamic_Column.cpp - -function Plots(type) -if nargin == 0 - type = 2; -end - -%mex CXXFLAGS="\$CXXFLAGS -std=c++11 -O3" TC_mex.cpp Cortical_Column.cpp Thalamic_Column.cpp - -if type == 1 - Param_Cortex = [4.7; % sigma_p - 1.33; % g_KNa - 2.]; % dphi - - Param_Thalamus = [0.034; % g_LK - 0.052]; % g_h -else - Param_Cortex = [5.8; % sigma_p - 1.8; % g_KNa - 2.]; % dphi - - Param_Thalamus = [0.038; % g_LK - 0.052]; % g_h -end - -Connectivity = [ 2.5; % N_tp - 2.5; % N_rp - 5; % N_pt - 10]; % N_it - -% stimulation parameters -% first number is the mode of stimulation -% 0 == none -% 1 == semi-periodic -% 2 == phase dependend - -var_stim = [ 2; % mode of stimulation - 200; % strength of the stimulus in Hz (spikes per second) - 100; % duration of the stimulus in ms - 5; % time between stimulation events in s (ISI) - 0; % range of ISI in s [ISI-range,ISI+range] - 3; % Number of stimuli per event - 1050; % time between stimuli within a event in ms - 400]; % time until stimuli after minimum in ms - -T = 30; % duration of the simulation - -[Vp, Vt, Ca, ah, Marker_Stim] = TC_mex(T, Param_Cortex, Param_Thalamus, Connectivity, var_stim); - -L = length(Vt); -timeaxis = linspace(0,T,L); - -figure(2) -subplot(411), plot(timeaxis,Vp) -title('Pyramidal membrane voltage'), -xlabel('Time in s'), -ylabel('V_{p} in mV') -xlim([0,30]); -ylim([-80, -40]); -% vertical line for markers -for i=1:var_stim(6) - %hx = graph2d.constantline(Marker_Stim/1E2+(i-1)*var_stim(7)/1E3,'ydata', get(gca,'ylim'),'xdata', get(gca,'xlim'), 'color', 'black', 'LineStyle', ':'); - %changedependvar(hx,'x'); -end - -subplot(412), plot(timeaxis,Vt) -title('Thalamic relay membrane voltage'), -xlabel('Time in s'), -ylabel('V_{t} in mV') -xlim([0,30]); -%ylim([-80,-35]); -% vertical line for markers - - -subplot(413), plot(timeaxis,ah) -title('Thalamic relay I_h activation'), -xlabel('Time in s'), -ylabel('m_h in mV') -xlim([0,30]); -ylim(get(gca,'ylim')); -% vertical line for markers - - -subplot(414), plot(timeaxis,Ca) -title('Thalamic relay ca concentration'), -xlabel('Time in s'), -ylabel('[Ca] in \mu mol') -xlim([0,30]); -ylim(get(gca,'ylim')); -% vertical line for markers - -% [Pxx,f] = pwelch(Ve-mean(Ve),hamming(10*L/T), 2*L/T, [], L/T); -% n = find(f<=30, 1, 'last' ); -% -% figure(2) -% plot(f(1:n),log(Pxx(1:n))) -% title('Powerspectrum with pwelch'), xlabel('frequency in Hz'), ylabel('Power (log)') -%save('Timeseries', 'Ve', 'Vt'); -end \ No newline at end of file