From 4e93655e9a266e10efecb6210a2333cbce9b826c Mon Sep 17 00:00:00 2001
From: Josue De Santiago <joselotl@gmail.com>
Date: Mon, 30 Oct 2023 15:30:21 -0600
Subject: [PATCH] 45 parallelize point est hist and naive stack (#64)

* First working parallel version for Naive_Stack

* Simplifying the logic in the iterator

* Added Output_BPZ_lite.hdf5 with the new format for pq files
---
 src/rail/core/data.py                         |  12 ++++
 src/rail/core/stage.py                        |  21 +++---
 src/rail/estimation/algos/naive_stack.py      |  62 ++++++++++++++----
 .../testdata/output_BPZ_lite.hdf5             | Bin 0 -> 37192 bytes
 tests/core/test_core.py                       |   2 +-
 tests/estimation/test_classifier.py           |   2 +-
 tests/estimation/test_summarizers.py          |   2 +-
 7 files changed, 79 insertions(+), 22 deletions(-)
 create mode 100644 src/rail/examples_data/testdata/output_BPZ_lite.hdf5

diff --git a/src/rail/core/data.py b/src/rail/core/data.py
index 74dd10df..1e336b58 100644
--- a/src/rail/core/data.py
+++ b/src/rail/core/data.py
@@ -35,6 +35,7 @@ def __init__(self, tag, data=None, path=None, creator=None):
         self.fileObj = None
         self.groups = None
         self.partial = False
+        self.lenght = None
 
     def open(self, **kwargs):
         """Open and return the associated file
@@ -309,6 +310,17 @@ def _validate_data(cls, data):
         if not isinstance(data, qp.Ensemble):
             raise TypeError(f"Expected `data` to be a `qp.Ensemble`, but {type(data)} was provided. Perhaps you meant to use `TableHandle`?")
 
+    # @classmethod
+    def _size(cls, path, **kwargs):
+        if path == 'None':
+            return cls.data.npdf
+        return tables_io.io.getInputDataLengthHdf5(path, groupname='data')
+
+    @classmethod
+    def _iterator(cls, path, **kwargs):
+        """Iterate over the data"""
+        kwargs.pop('groupname','None')
+        return qp.iterator(path, **kwargs)
 
 def default_model_read(modelfile):
     """Default function to read model files, simply used pickle.load"""
diff --git a/src/rail/core/stage.py b/src/rail/core/stage.py
index c2a1ca1c..c8ca0ea0 100644
--- a/src/rail/core/stage.py
+++ b/src/rail/core/stage.py
@@ -326,10 +326,17 @@ def input_iterator(self, tag, **kwargs):
             These will be passed to the Handle's iterator method
         """
         handle = self.get_handle(tag, allow_missing=True)
-        if self.config.hdf5_groupname and handle.path: 
-            self._input_length = handle.size(groupname=self.config.hdf5_groupname)
+
+        try:
+            self.config.hdf5_groupname
+        except:
+            self.config.hdf5_groupname = None
+        self._input_length = handle.size(groupname=self.config.hdf5_groupname)
+
+        if handle.path and handle.path!='None':
             total_chunks_needed = ceil(self._input_length/self.config.chunk_size)
-            if total_chunks_needed<self.size:  #pragma: no cover
+            # If the number of process is larger than we need, we wemove some of them
+            if total_chunks_needed < self.size:  #pragma: no cover
                 color = self.rank+1 <= total_chunks_needed
                 newcomm = self.comm.Split(color=color,key=self.rank)
                 if color:
@@ -342,17 +349,15 @@ def input_iterator(self, tag, **kwargs):
                           parallel_size=self.size)
             kwcopy.update(**kwargs)
             return handle.iterator(**kwcopy)
-                # If data is in memory and not in a file, it means is small enough to process it
-                # in a single chunk.
+        # If data is in memory and not in a file, it means is small enough to process it
+        # in a single chunk.
         else:  #pragma: no cover
             if self.config.hdf5_groupname:
                 test_data = self.get_data('input')[self.config.hdf5_groupname]
             else:
                 test_data = self.get_data('input')
             s = 0
-            e = len(list(test_data.items())[0][1])
-            self._input_length=e
-            iterator=[[s, e, test_data]]
+            iterator=[[s, self._input_length, test_data]]
             return iterator
 
     def connect_input(self, other, inputTag=None, outputTag=None):
diff --git a/src/rail/estimation/algos/naive_stack.py b/src/rail/estimation/algos/naive_stack.py
index f6b629b4..9d3c29d4 100644
--- a/src/rail/estimation/algos/naive_stack.py
+++ b/src/rail/estimation/algos/naive_stack.py
@@ -44,18 +44,58 @@ def __init__(self, args, comm=None):
         self.zgrid = None
 
     def run(self):
-        rng = np.random.default_rng(seed=self.config.seed)
-        test_data = self.get_data('input')
+        iterator = self.input_iterator('input')
         self.zgrid = np.linspace(self.config.zmin, self.config.zmax, self.config.nzbins + 1)
-        pdf_vals = test_data.pdf(self.zgrid)
-        yvals = np.expand_dims(np.sum(np.where(np.isfinite(pdf_vals), pdf_vals, 0.), axis=0), 0)
-        qp_d = qp.Ensemble(qp.interp, data=dict(xvals=self.zgrid, yvals=yvals))
+        # Initiallizing the stacking pdf's
+        yvals = np.zeros((1, len(self.zgrid)))
+        bvals = np.zeros((self.config.nsamples, len(self.zgrid)))
+        bootstrap_matrix  = self.broadcast_bootstrap_matrix()
+
+        first = True
+        for s, e, test_data in iterator:
+            print(f"Process {self.rank} running estimator on chunk {s} - {e}")
+            self._process_chunk(s, e, test_data, first, bootstrap_matrix, yvals, bvals)
+            first = False
+        if self.comm is not None:  # pragma: no cover
+            bvals, yvals = self.join_histograms(bvals, yvals)
+
+        if self.rank == 0:
+            sample_ens = qp.Ensemble(qp.interp, data=dict(xvals=self.zgrid, yvals=bvals))
+            qp_d = qp.Ensemble(qp.interp, data=dict(xvals=self.zgrid, yvals=yvals))
+            self.add_data('output', sample_ens)
+            self.add_data('single_NZ', qp_d)
 
-        bvals = np.empty((self.config.nsamples, len(self.zgrid)))
+
+    def  _process_chunk(self, start, end, data, first, bootstrap_matrix, yvals, bvals):
+        pdf_vals = data.pdf(self.zgrid)
+        yvals += np.expand_dims(np.sum(np.where(np.isfinite(pdf_vals), pdf_vals, 0.), axis=0), 0)
+        # qp_d is the normalized probability of the stack, we need to know how many galaxies were
         for i in range(self.config.nsamples):
-            bootstrap_draws = rng.integers(low=0, high=test_data.npdf, size=test_data.npdf)
-            bvals[i] = np.sum(pdf_vals[bootstrap_draws], axis=0)
-        sample_ens = qp.Ensemble(qp.interp, data=dict(xvals=self.zgrid, yvals=bvals))
+            bootstrap_draws = bootstrap_matrix[:, i]
+            # Neither all of the bootstrap_draws are in this chunk nor the index starts at "start"
+            mask = (bootstrap_draws>=start) & (bootstrap_draws<end)
+            bootstrap_draws = bootstrap_draws[mask] - start
+            bvals[i] += np.sum(pdf_vals[bootstrap_draws], axis=0)
+
+    def broadcast_bootstrap_matrix(self):
+        rng = np.random.default_rng(seed=self.config.seed)
+        # Only one of the nodes needs to produce the bootstrap indices
+        ngal = self._input_length
+        print('i am the rank with number of galaxies',self.rank,ngal)
+        if self.rank == 0:
+            bootstrap_matrix = rng.integers(low=0, high=ngal, size=(ngal,self.config.nsamples))
+        else:  # pragma: no cover
+            bootstrap_matrix = None
+        if self.comm is not None:  # pragma: no cover
+            self.comm.Barrier()
+            bootstrap_matrix = self.comm.bcast(bootstrap_matrix, root = 0)
+        return bootstrap_matrix
+
+    def join_histograms(self, bvals, yvals):
+        bvals_r = self.comm.reduce(bvals)
+        yvals_r = self.comm.reduce(yvals)
+        return(bvals_r, yvals_r)
+
+
+
 
-        self.add_data('output', sample_ens)
-        self.add_data('single_NZ', qp_d)
diff --git a/src/rail/examples_data/testdata/output_BPZ_lite.hdf5 b/src/rail/examples_data/testdata/output_BPZ_lite.hdf5
new file mode 100644
index 0000000000000000000000000000000000000000..83faec2fc63d62fa80af2be4e9cee8d35aeadf7f
GIT binary patch
literal 37192
zcmeI52UHcynzlC?1&Im<q6rig5mW^AZAB3yAfku~6>~roB^(2aqGCW0l%OC2ii!xE
zEJy|cL2?!(XC&t|=q}E>bN+Sio%^p@GiUDI&017dcXfT=^LBN2?W)?l?rN-4pEPm$
zMB2}do14Z?6ZjP!|2IlD4Lio3@$~p{0}^2&;T;nCureYnw4Z-8+ElV%kd*(&d^@(U
zUr!s$IbQyFNKT`bb1-7u<2#t)zm)<S>(^-g1lA+eQY1_!^^7EGXO5U0{24bq$OsRf
zJ$&NySbDsEDvxLUFNS|se>K<`1so^!@*mY-Wm0iv|BH5++5WW@_+$NDOz`fI`ZG>@
z&hY4gvj@){`P=lRBj+!iFgp8dm)WI*r!S1f#_RJxs!uxGG}p1v*nWTh_T6EWJbu2^
zT}GO^#|X!t{`E8eb3TDtjP&9;jF41kniRPY#>@NTeZk26fB!!Vjf0^%0c9bIvam3A
zjPs9&Y&7m)`Gv=JaE{?~kTfeR>#r=_<a{<3%Lxp1e&+hy&FOCwf4TU0y&A6vvRr@P
z#qsyQrQqMsH(s9~3H|+i|0+jMoHaRe{@lL`4U>bZ08;^`0!#&%3NRJ;uPE@(w<`sH
zwJQb28eWXH7xl}`pSLfqN&9lVr8e%>qw&aJD*S)-V6+HX7|rJ2Kldyo4ejyv_IUc|
zC!L@B1f}u2K1w_^rBU)Z<<E7poW>C1IMHid7PIjPGh`~jRDh`fQvs#||LqFkrbWdD
zzsbMBfrM^h+D*44!r(>HZ=8Q_F7PA!1AZRWZblF}9y)Q=InRQL2q*a>WxL~S_s&BU
z$saRzaMd|4#FF#k6@P6=Am=BJ<(BP^$7fPb@-Oq-kwVH#{pp@X1=2{l=|8QlGR9X@
zekMZ_WD+}Bzp`r}i`dKlS$v!=a){lWpYC`TAeY$B8>4BuC5b%ZNB+^D+c7Y*O+>*@
z>l^$~Nc<}Ni6`5GBH~~1FZ|UgCVrOul&9$yz7c=FrTyH}QA+$S`&HZ&lo9{G|8j3H
z$_Ym~u~SZX%85Nj@fxTgb}I;9C9z*gI4g-CRfM;S_)|r=tBGILguj~jS3`8DA%4~n
zJ!*-+wM3U%;&&anUUkI(I-*lO;ixAZ^@O9Ia5RwX(?B>H2uCB~X#B<T2#tiJiE#YN
zHl}M6;b<lt&4i<waI_GP7Q)d&IDQb0AB5uv;b<ist%ReMaI_JQHe#oZ@U#(oZG@|h
z*li<x?ZkdN;cO><v=iQT;!iu_?jU}35dIG0UkACK9mLNLqDLq3x0C47N&N04*SVAU
z-$``pA{<?Wql<8K5sogx(M34A2uC;J=q4Q9grl2q{K__V-*yv@9>UQ>IC=<258>z`
z96f}imvHnFj$XpiOE`K7M=#;%BOHB%qmOX(5sp5>(MLG?2uDBR=qDWggrlEu^poq+
zPdEk$#{l6NARGgvUJekB0m3mzI0gyFAmJD!9D{^okZ=qUjv>M^L^y^B#}MHdA{;}6
zW0-IZ6OLiRF-$my3CA$u7$F=Zgkywoj1Z0y!ZAWPMhV9#;TR<xql9CWaEvlIC<X?N
zqA~1HGzJevW7wl;3@(cH%kEK%#^9r9jQT*)7@QQ1;Ri)y@KQ8}KNO9@O|dZiqF5OG
z6br*YiiM#A#lrBDVqxe(u`v9lSQxrcEDXOX7KT0)3-O<Yp%X>MaVb{9!Adw-2?s0T
zU}fk>u@Vke!of;7SP2I!;b0{ktb~J&aIg^$Hp0P1IM@gW8{uFh9BhPxjc~9L4mQHU
zMmX3R7ldLb9PEUHop7)d4tB!9PB_>J2Rq?lCmigAgPm}&6Ali-!9h4U2nPq@;2<0v
zgoA@{a1ag-!ofi}I0y#^;UL3N6er=}BpjTCgOhM@5)Mwn!AUqc2?r<P;3OQJgoBfC
za1jnJ!ofv2xCjRq;ou@1T!e#*aBvY0F2cb@IJgK0>Gt@uar{4yYyP8g{II!<>-)#v
z6mlCfgMV=a{y2^=KzOytxcE4}`5(sVHO3G8e}@04UyLz<i$)vo2m53C_kJ;{f3{z&
zem=uxf&Ks4_5a;|F~;~a3+dyd(R}`(|9|$CnPMuyRN&uPfma5qF8q5FVX@>>UjM1D
zP*RRrk+Uif_fxG_8c!(&%U;Eid$+1^L_aVoGoS%)484WcJG3EV#N!M9Ko@%DoE4sP
z3?P5U(z#1sjp9sd$-*yXtkgE!fyOF+PHIrf(t>Bx1d7(9R8Vt@m-?7?>5!B*Kjm_b
zlD26Spjd{^9w}=Sr0V&WaqMyxqDrD(vKO`sQJ<_k3wJaKQJytRw{N{6L}i;ll#h)Q
zq%>v=^F9j}pj=`$(e~u<QO829m0pzdQ0#5Cj(MTn)QUiPx30?^l+NrL<0o3IRLSWF
z!aUDLpwrmCS!{AYUaVNQd-j|z_;RaTE{<r0Qqo$1Q<n9RE;y&~<bEX<vd-SM;cf})
z!#g4dC*)$~>nQGvPHFg2LhtIPCBW(E8|mxkA|N1YZU6WJ9bYGVik{fvgI5z@1hbjA
z<KMY1G1-|4Fcn}bz*OM>o&rwhDfaT)-0}Umg_>!YH@0;Kf4=SU9wB-UMsAe_B6nrl
z1KA}XQLfXV)6^dW-dc|Q-5!ZZf3`@aM)V6tS6AIy;g<;o)Xwc+oQo$*ijJjn7GmM)
zV_SQse}hDrhikfQIoOkwuT$$PF}jqlqP?XW4hNK$FIrHGsuOqT^JUZ_^oi}MH_{Cl
zIDY@Lzi=ZO#W(Jo?%W7k*`TZByGA&!)KtsRXhdCWyv9t^2Ar3O7q_G8v6Y&dwlt^~
zF7iQ_9jdD#;&NSgVQnQSH}8&-^l~WacpAf?6q1L_qQzDffjhlV%|j~>=IoW*3IejA
zu0y>mDo=ySjK+rcd7ok06lZ>2Ee?q(?6Mk~k&sk0J0fEj48L~SRafWwL+WsNd|=XR
zWT$0JV^{J(w$~hWV>72eeg9u`T+<+u;sX`S;}Kh4`C*w(g!;`VbnK3MzJ9WJ2n^CU
zUolP(M^Mm;i4$&r!hGA#im!Gt5W1A8o3SzuySK<3w*L^1i-s#BUgsxZ$~@{q=9vDT
zA(K9vb0i_%tkEKOFcDwW7OEKDN<`$kRZ*Qj3AmQpf!6c!P|n-1K(#9tNB8n&y{BUE
z(88AGV{s((ZiuVYiihLdqAxoyoD0VDOa-~V_5gUwOzzor@jVu}*^I;nzlQ#j1q&aB
zyo6+9f%W3vXMf81f3M~<J!UGvRDh`fQ-PlfEZk*()t%1?S`+VtXMAwOh3v&&l8<_#
ze8U&NKy4rFpobm5zuXTObe+lqT<O>=$erRI8Up#*(_6K#M8H*YezEf97!<dE)yr~<
zN5iLv>;#d|@a9m@Y3NQxS;VDFa}_giBR*@GqdE(c7q5R@ubYdO>5CQ~t<8to-HMY>
z9u`66<+o8~wQo4g{k~~RUm4Ciemj+uQI7U^1A<XgD#7WzIDU_M73`Wi?KkUJ!{oKb
z>RX#?u+pdYwwFLHB%Yq1ukgASa?Te=%O=#phNVa?ZfPCls;1pAnOBFlj1&i%np(`u
zmD%HcvKAWKC5BW#)u3*}%a*4V)#z`V!}qYX3M(~kA1w;4#F_aoTry27uq1FtnCP5x
zh@S`&@P1f|g&uD_qQ4Ylvm%e0GJhc+O}QMi_edVBgCrD73$tPC7yan{`AnGidhldi
zNJG8a-73-J$yhmfvF2%N0=C={S;+k|7N?sNPjEGSLhw<3R$u#YlnosoIISCkj~PEA
zRW*a~#5MWrCEq}FU#l*g>l=i*oYoigy+Xmiu39sU;}i7K+Lx63#lg;e^D^Ui$$0!u
zXZO)FnJ7Ot;&3WA4;Lc0b;rhiL!EwVzTz0qCJT-I;%xP}z3klLF5MqEc{7oBS4jto
zQbLc}^Y`L?o8J<z$U)Exc+YM+Mxzquy{HzL%tp}#c+Wd#aZoJtWFOG(a#KHKWj>`!
z@lf(g&g++5okY!<F(ZxdHb1o<uO7wQ2~f+_wsEjY2~t$F)+fDVf|PLL+Jowbf|L;Q
zquE6Tsa?-)&R?+-pd6RYUD0xopOUB%T-UvQ5_RFy3BBnZ6RFAr$~SKB=B9$La_|Rv
za!}7Wtn2(6SgDH1qC?SZXq3aM(e~$;2O%qgy|KkT@Lj#fn<us#-|U>iKN<FeGgRoI
z;;LbU<eE?EbfHnpm$gT8y0TI$q;9ZBKVqkz9cgK4J<UmVroNCCmExvk*<KDVah*V!
zJeeN#B6A{TEAXAGGl++}D70489lTWLd@9^}tRKsVmQ<p-oR<n=3ywbH%|nS~?OW6l
zF_Get<8+>Na{?8;bkMqL0yib}Nq<wPJ}31-yZ(iYH9O^5`8d(lg_T-e$X(U(ibmNO
z*A*T&9!B>|`%{lY`@x<Hsm)3~IOpJZer9h6t`_O}t>S4zOjo*n(~q%sFv|y>T~})H
zI6`X3Jf;$ST4!>W9QuyspRLW-S{I{ZiNXxO^n8>$?hU33=R!`)WUlt&Z0Kn>+!Y(h
z#7?DU3pq8$>f>XvhL1HFI6!~iY4|<^MjZEM%KBtrpM82^phpI5&t5FFGtPixgTavR
z;8)0*KPX(e_$%6!j8xkuro+#*&T^IA7wkDKef3&q3dC)eJNTD;#tM32{lnx$EGas&
zoX;;F7PcEx<2_<Asm6PM`{`($4^Nfn?uo$qP}4O#JHrsld+Jt-Q3wtn<=icHCJ41<
zOGIO<>3Dy9a|u!cAaH8vXg^l~tll+U;<)DzRxY;2Nvi%>;Jmtf28%x=PP8sP_viyI
z6i7MT86JBs9r(VkYmpx~*43L13cpABPRZ*--tTbuOO8=S-5ZpC&R%h>`ZcD@H@TJj
zdSlntHrmpip3svmerMS80y|1Z^hGtELojj6lv}k<khx4DtHy!Z-etZXnF{=k0yQn!
zv*}g#_|2*!ag&V$y7w<z_fq2tY}=Ql(8QhKTh~c{$l(II&sG0@v}ceM=g%?XafAFo
zw&H%F7l<;`&73;N16tFTuPI*n67$5kx7+G@;`<zd)#9dJxO`o@@W4}VtWJEbBpdb$
zCG$L>TKF2<_M1j%_xM1*;@#qDf^X4%;pE4qGVjpcd|PY%a$kgA>u9RNdsK!C@AHQr
zCU2LF(O&ui%R1SvM)CXOxaOG!SFZY_zr52#x7HuSX5T&~sRux`>bccH_W&&KX46Y)
z4Zxi1y|+S^(=qp2-?`pPbbMMU``q|F9X#hBY0WI5qo~F*RDu==|9X>qkHrGfwR><=
z!s0;e&$P;)uqqI08M03^*92l^u$uA{1Y*ec&h)xvf%x9I(0-3>ARJuoDmqRMgmmn2
zrIkZ;uq}QdIsBcDS<`IFtYhfNnDj&Vt{WYd0)yY_7Ic`dwyg=!rQ`gB2R9TJ)6t*R
zCs90+4!>D$##*%jU=4SPREiHkK4)u@^veLCkFMZrHTK->w6!xh9ss^6*+X$!0hlMO
zZMk(-09G2L$r&vQz-{)O%kRw|JD$*K+$R|TAF-?fnsfjr9}d`GDHDKg%7;Y0Dg~hK
zWPwTJh5)>O8a83a!2k&Hy|<LS9)R=37d!U64!~*Mk<2N%0r=h$q@>S9M?}VOw8|Pf
z?5=4RW?iS_oo2?)4JmYJ3(hji6Agrh*4&)*Gl4i=zqH%!OCV+#2C5V(1Yus-k)=u=
zLHOY!@AY_kFv@bI6nWi)(WtdcqFz1()7mxHz55!1m~6{Fkt?BS5Sf2$;i@o5Tn;i+
z5)4Pu!`!UU!El(!TvoU@`Vscgd<LnKk#N?}x?FVN6E5~G^F16Ng`r(NyBEmBAh4I$
z+JrS0*^X~?V!C6|IJL0m09!n)XLg&0FH8VOr}n{hrinP&V7I}lG6`ExO|l&}PR64>
zr*0_Areeymfq)NOX_$PsEa0-xSICzIEWSD?6REAI%^Y`UK`u%&ZSSjWsC+h6vEa?c
zcZ+#PwC?A@^p#j%{hR{ydM+?7EG$HL+>RCEUd6CGb+f+f);HvJKdkZntqglD!W}}W
za%kNs^`&jBKzqTE%XXbgL?vB~IJ2h;&Ml8x<@Q#C9*{K4O0xzoFPx;GEvm)y>Tow&
zM=cyiHV0MO)s2mBeOPyevmO_3<y^h1QIBbXLws@C^(b0Wx#)y=Jyyw!XeRj7L0Qs%
z%C*i~L^#RqsUNMuzz)%rk?d*=@@L-MU{VDq-j9cjS}L$+?ev9>Gn&xh+n_q{L@Vx0
zEDE06*N(L9iDJKLc0t!|m$kZo4;uD(C)@J(<5Jo!Q;9=^co*VRQ1)yX=WEtb3c;g@
z7Cjbin@OY8g92_uzhj}=rF`{%9Al;AWhWVJVqv3pdJV5QY{Eu~y|i<Zk7A=d{b%$`
ze`BNC6+#YgPGh4AUmTLV|A39667LJAOR`b%g8c`TZ?aN6_dZ%3|HMMoVrJ>}3>x*}
zz}uMymZO-FdB<yD!7$iEbkwE?58%zf$kb!=`(U`k*|)>38;Y`mpBBvSgd(3#+V_+;
zXkFeG`CO$16Psc=g|!=y5^40Bv$_ULhGx2T_f_IS-IJ;7AC@C|@`QED%BAp=vj|m@
zD~8<pwOMJo1+aT@Jcbg;N5<Pj+QMPE@K}52)f=B2G;6zh-i^z~O+;(|Xv)Hwu!c#m
z1+!o<bE4u#nM_F5wokF&`W1yVp|cAl)6n+t=);246y(Q>J>b~?8G|tsMfn>Nki80h
zXRgG-fJ=0_s$2}N*v@)iI~oC@6i!ZQo-p*T@!Ao2GYER_moquf2Vm#?(X6Fj@6lK)
zvf5_V8#rhE=2GSAjR$Fo)=95k;^L!@0xz`}IPYOPtodQ=x$Q;wzhUo;k@}B%zN}B-
z+`qeox6lD$PnVi4?XZVJY239f%64e6S~1a0!4|*W;M$Nbc%Ru``kQwpGs;weslfk?
z0<$;dUfHqX9`^7r)}SjtLi1?kr4r@Gpxb=uSD-%umx%Ufo*);<_~}yonq%W#+FqtE
zn_eQvUDUf`i#G}?*Ie-k^MNW~@zh`^U!2WP(!phau)D^@*0<Af=9{k7$-rP}_q?-d
zO$tK{*O2$iC6Orgm#ecIibm{yo2BYg;;{N6tH`iH0#aAM`kr5vgpI$g=PPti!KsJx
z=d7Nj!7hJB*T>`xl!zZWs_B}A%PaKigN<?^)~-L}y+Iy)Jz2P;PZwZXy4TG-yCOKf
zquVY@FF_9*Kg*68W#|vc2;x~;j&%2?3mc6pFyo<P`ctn;I2FgwC@QRi{-p1$ueocW
zt~o{Fz~WjQ=oDC|xTOx})3}{g_3ANc@#NR{wHmO(>yyDYg+{dJ^7ZHTG$N6kXUoUu
zO)#r}F5Ek_88Zh9SO=||QR&QC5*O7B)e0l~rmSWh`p_L18`%st0~wVq*P5Z-#I3n~
zdNTr9{B8EwHi28b!Rl>!Bd)e3Bq#`ul`B4IVl}$~M~vh)MM~6TbHWkplk9aUd$M!T
zxVQ$#!(Hb6@UDi@W?^sP(^c>dmZ=w9Qi+N7y8b;w<(ME{5hdgK9a<`Wx}nEQQPecm
zC}(>K&X)tB$BSV1#(&*Y{{qMu6bV0)&&Q*^4T&xdxtOVU`{b>P97ql&Ik!yBK~ek(
z71c-C5V{eO>8zd&_iWy!tt+z8IHl1v<7hT!c9gcBP0mK1LVdY|UJm5sMuP46bCJ7U
zb?*A=T=YpUc{nF657UC0+l?dgk-S1=rLlJb?i>_SD!5z-DP8Ry8cIc|y3x14KdJ}?
z?&Y5k3Krwpb;IkO(#3dK;fBu*MTlIv_x1cEMNq9y@Edtth@)>e_7|8Jz*1#iL!o#+
zI+xA<y3Qt-`5yhBd2BGPFcn}bz*K;#08;^`0)Iz=6JC)6{nOq<q2y`Yqn#nBsG`5N
z$%?^sXEhJCf@Cy(TfOhhNEU4Dcqhk>wOa*hA3sgVsz9;E7p(_->!G;JJjF`44Z22$
zRQi7F#_es@<{v@_;j1sO)8Q@))uec><9-4M^_zaQ^FGUo6#eAom9F0Wl-629wctBK
z6eVeRXxio}l-~6Qx1<qK>cCy|nxI?LsW?OT+0$dis7>E!F6oaYsC^3}#UiRDDR!mA
z#D`i^)DQN)VX*^JRGew^aRqTHs(#O$SSfKyYKwx^hRc`4DYipNH&uK}ptrhPdbFYx
zHbq)%_i&b@L6~>-OVJ8M`s6j&&#Q!{Utg@SS`{p^r!D_{yc&)Z8Xu^-*Ff;;&HEDd
zwV1j#k2hpZJuDX8J2K~01G;ro+#GqDu+MAe%BzOW7%;oV=J=uo0y~Ri#Qc7sENG>^
z$#1O))Vnw<yss4_d_1dP@wOp+enzxFUMsF@`%tR$TCv^gOU$dOKTz2|^`6JOW*jqJ
z>%6A55o(iG&D_~lkGJ7Tog!arAwo~l33*hFIbnLK{hF1K5_p!-FJ6vTImJb5CVxXo
zu-=C@st{b8cW-me&V~8bl*C2ZS=f*kKKqH)R}9ZrO~0g)3ggJthQ%RC*zn$JY1H?4
zB#fAbEOUy5ULen|BR8WVdu?h(()CE}Zc4RV-5n0rgIq&;L7})kGH+{<TQD9zzua8!
z7Krv%FTt%o0kEHIs#|vE17zk!2b4I!LkFkuqMY*Ac(k8gv2(i@*wx~17tlPQJrLs-
z+wO{236tu(Y@Q%UV)D?VLOVRoP&;4}b`SYB1Ggn5ZZX@%f9D=$9%3rMRN$|nfUC0P
zGnr#I@nkqh?WvwMJ}4hQB)0tyW`@f)U)8ih%EbkU+HxNXqsw(RsXrVWpIX1pL+la0
zgqgoB7PNy$U~k+xR(nLsIkkrWu!p1F&E-`Ek8v{LeKT*21L6)%T7B%TBQ9%PPTJ!7
z1h*Y$t<ZOR3Z4WPpIbgoXm^{hsF>>vcH1Lr8gpEsq&Aqc^7S(;_s8A&$J{XQ$WHH;
z6)$k+aNd;0vmUrVQ{Lx8$xDc*ny7EI@j`XN2fH`7Ucu5_$ntHG4;JPwdN%R#I~a)`
z;XW4a2ewr;{Z7jQu&&{|X0~o1A`9uOcP0jpeNQaP^u&WOh{TtN*oH(vXCx%vd}A~^
z2fS9By^e)3SEWf=T|DHQOJAv<OoG*M`=F*JDX1P;E%t4F8Z33iuGPQFfQ|R}N6Eff
zSl2?$GWQ?rkD1^8qu_fU?r?or+NN3vyPH!>e@iR|rN6_c@Ny~M^K8f)SzHd4=1*6Y
zKUP5hvlXxV?ke0a=Qw|fuLh2zhchD!#^#6hxh<3ls>g*N2V*5(j`dUR-n~BBy%}$P
zcQ49+JN90-xOUsvpbe|8*m>MHYRAj9Quh?DcfhWFzH;E>PV`m=JYDbIg(f<SA4^a-
z5@o6$G=Azq#PW#coiV+L2`JU(i0DIBpXXEVxBW12zg^sTZvYg_yK_GK2gl~Y<c05$
z8Up`_?-cQjAq*dLP~UfA7>cJAPP(=XgX`&Ud5`o*@cB~T`b8f`;BFD3>()2|+uQSn
z)OknIW4XXPUU(D*2XZGj(MIv=hUJdfoDp!=iT7Q+KY}h5wMl#BMi3YynNaCJ3~@7q
zu!rKq5C|B|=rSHc($*-iP`^PqW;>m@lRJP{5|y+C)&1Dq&;H_OSs${$vZcEx^kVj+
z&f?EcdoXX)bd9IGyD{$~2c3tv3%3vW+BgPuz~!sqE$OZ8plf($DkQWbY4%y_eSI_b
zX4*+h2{dB*@XWb}YIRW8F4%hfUNx@ho*LQ{Gzi*7b4{g-G|EF)Y*CXU8^t4Qx7VJ9
zlghAuV_cNTO$nZs7d>{Ahw`)0;z~-MM76yXakS3mr@l>IH1PJL02LbXChqVyK}xx}
zk(S~iNF6PGw{b$SAf<I<(i+-*L2CV6?-<{Of>cmh@fTfB0qRF+@q+#z{M6v7=PuVz
z@KG0&baX$a@KTosi)Me6<e}F1eq5ktI)Sp0we+}Bz(pM%n_RYAjgyMHcVLh4M|LVD
zLp^@I3LB*#s-<8Z%|dAk)fBGXK%??E@g4DYABKrme2zbT059CR(xbohB1h<M7CO41
zDiT`VD&K+G_hdw7OSa;Cs28<xc_R*7Fi+18uff4|R-Gk572wlf;&y(|H|Wc0)HG@4
zV-EkVbg9)@h^O(~7@nUxR-XO`$*b{jiTd=MZW#%S=c3=jFN}@r33E;g%kYJb+yW25
z7Eg#helGl2z?IoP`fEJ?nO-mzU@E{=fT_S=T!EnVPV)X8t|&V-IcU;E4?Ms5<;d_=
zFXV}NSg39HffQT5*GDs7T=;bAZQeP5bn>>{IT1m}E>@$c<9xxejFLS>(+|Z2d5yf{
z`{8KU@;Tn@5CPkUU<(6-Pbk*+r#~r+LfCyVmxr_%7`nO@Y}_A%`AtL8^qVoLyx%xW
z`g{!NiRl#zWAoBtS)^{B+Zv7Xy%sZNn?B)j%Y)=3rAQ0}mGmuD`G|<f<w4UYg+bKv
zzDsU&F!ZKfZ@z9A2&%R)XgcrMyy-|*(`v8p_?kY-ZTzDGr+0r{+NLztFZ3-_Q}Jyr
z-gHH;s^8uKf6+yH#}u0pE@!kda$76fA}X#VZEr{IMaO+tr*~qa$<&~{vEO-U#Z-5T
zma+L<g=}v&H1xqK*<X1}<p5Hv8VshT4xvn*w?OC72=0IPHi#3XQT5ApSOmgpl-k7F
zxA{w1s9?*Tn@22JC^uE*n{ymlD2J=8AreL`l*g_mb2jm@P)j9-_2NuulxwB;hu7Oi
zVOCubbwz&|1w3eYAUTMTLIVY*z&;oZnpvnW=z-R=9r1=9o#0Wv7cnB!jtN%ZJoUX>
zP&Sn<t=XUfllVH$opGvxgAkQR*Q<o{@Ll1i&EGM@mGgy4e=!{1<<tp{joZdk(+)0H
z$b;;*-fnNR97yoK8kN(@hRzq`csZvm_;!A+*dv{V+ijEDGt)BB`hLyz3s*D2qEx3^
zD3uB4UY^UNn=+s;l-epa;VZNp=grnKOM}q*yH=)dsi?AB<j!$78OxS@k($3L3G(O7
z+qHY*F=f}o!IrUk@wIE;m2HoX!n;F$7v^34h@gC~4fakUNIcTMQByn+gEIqGy%!pr
zulmZqv8C-T%<|To_1L|_iRIg-#&*4g=<@blW7-RR8?l&QujmT3kj<sVTb^R-Nt3<u
z(;ZNgmf6&E&JOjb<$|~GeTaCS625{nHc-8C>9}Z{HM5=g7k@4>-D4`iRN$|tK>1_$
zqpxl6A?Cv4)Wte>2;j=I)S)}#nZkv@C;YDH{eJXfQThuU?w)DF_t^_g-<2M_*}uV6
zu02W0Iv?;})_2$Wv3c-6VrD(pmkCAhGy(do)(A|Pudf>$7X$N2pYv13es2hJ>s@Dj
z`V4_fG7;XaX>eA#mFz$Ey}n!333_!hIf%$N*`}bFkDKj62hGOD+h(fq9p_dnMRTp~
z>zSVA;NmPA4g69G%d<EzYibQTBR{8z+^j>o#tGY+35`g8Us>exvKc-5)|PUpw<0~f
zWaH`XHay>^?s8JC1Cuj+DsSHIg!%Uuvv?!B#^ys+I;<<{M(|0m^Y-6+aKX6jnooQ$
z{7bv$gxu&upw7Hl`>Fl7z&Vwt=6XL?%v<HTHnJZAbg}vkN&T=f+tK6g-j5F}9F*p6
z?8iPurKO$;eNdkvs?j^S56wwWJu6oA!bo?)ohr2+WL)u+5na%Y+3O1pScf`s6h^m>
zg?FH=^qir{$#(3_K1+`lZi8{d!rjiMEm%B5zt8MzBX(}}`5~~h9!G_DW=RIt;N-Ec
zS)!U%xR$80Pmr$yBd107zpE>S`@P~b863snGUIYuW1bHe*{`2otj$4cih5A2K_=+k
zdD{m(zF-2cOZtGsXIS(<>=BTU$7(~pDT21q%=hPC@782`#8iN(08;^`0!#(|Z51%c
z9u}$h^uyJ6=ie4z4jud6BG(7kx3NgyapTkt;Z&%fXnxyQk&PKL=e{k>D@JmCYE^tt
zCAQwVl055a173^N&h=i<4%2RtJBn^Suz8r$7_ed(#osQfT(xGUR@sDmR8Hig_&n!2
ze3am&79UdexRWP99gb`)&`J`frVZMknwKg<aYQe$vP_spZPbvyu-Q?J+8v4f*)t@m
zch7dJ*11Sgio-j*P8v#6U;Gs!xxPzN;cVf$>1;EpF2lYXPcx;dg4tT7Yc-^)Jx%J?
zGpJ%ztrq*CdaxLKbNl=4hl^0MbwXMH^ddO_=CjK(sQ|7o%^epC=i|<{lf|p1<|2YK
zBUGd{3pYX+-Z*?H193J17cMlX;`45Lvaoj&IKwnm^?l<YufiSXtPq8d{fXxFKSJ?E
zI)Ct*GaYYsTd4N#@`dVwcoWkF-Z;^9{<bVT^LIY~ZEq|lI8y<p0!#&%3j9kbuyyy$
zcr$A^;MRk*s@bn`V%dbo`v^c<@ZCwl>%&oP@;>>vVjN^zCZ+l%rl5QA#!8;HEZF;;
zO&z*dh*`V*B1><SLsUkOby9u}ZW^uK5qP*69P`7h&&=+CrI=a0&6gg8OK^G($qZpy
z$zmPVsVtOX`o-Nn$Ji+$6V)Z!9bA;uidTv=Uhq(RtqOZi?ck>}4DJrKToR<VZ1xJ4
zoh?jdP1Z>CaGy+F^B((-V&N3(y5`LZsZk<SbyC11RSi)pxvXT?-7})p%MH_dlNX3m
zuKoURSKk$(iV8OwO}RIPisE-_q_3V#wP~Kp2niIT_-u0r;|B#Ok^O4BPHg9+*eZ<k
zLtH0Pf3_Zs)Lq8qrtUGqf3zM<@_j}==6W#y61T}e4e@8|!5EPGGhU#@86(3ZV~O$g
zivRo9gPHov6blQ(!ua|x<Kge?zx>ndBKJLFFlddVknqkQ(*L{3KduK8_6Pm{cfDg;
wXDYx{;IE{>_&T=~S%2>Lb#lqNxYYQ1a%8<bL^9UD16dcB`hA@_#`<*s0R#iZssI20

literal 0
HcmV?d00001

diff --git a/tests/core/test_core.py b/tests/core/test_core.py
index 4b72cd96..cad1bf5d 100644
--- a/tests/core/test_core.py
+++ b/tests/core/test_core.py
@@ -107,7 +107,7 @@ def test_pq_handle():
 
 
 def test_qp_handle():
-    datapath = os.path.join(RAILDIR, "rail", "examples_data", "testdata", "output_BPZ_lite.fits")
+    datapath = os.path.join(RAILDIR, "rail", "examples_data", "testdata", "output_BPZ_lite.hdf5")
     handle = do_data_handle(datapath, QPHandle)
     qpfile = handle.open()
     assert qpfile
diff --git a/tests/estimation/test_classifier.py b/tests/estimation/test_classifier.py
index 513657fa..a9a7b99d 100644
--- a/tests/estimation/test_classifier.py
+++ b/tests/estimation/test_classifier.py
@@ -12,7 +12,7 @@
 DS = RailStage.data_store
 DS.__class__.allow_overwrite = True
 
-inputdata = os.path.join(RAILDIR, 'rail/examples_data/testdata/output_BPZ_lite.fits')
+inputdata = os.path.join(RAILDIR, 'rail/examples_data/testdata/output_BPZ_lite.hdf5')
 
 @pytest.mark.parametrize(
     "input_param", 
diff --git a/tests/estimation/test_summarizers.py b/tests/estimation/test_summarizers.py
index 8b15ff34..4a59c4db 100644
--- a/tests/estimation/test_summarizers.py
+++ b/tests/estimation/test_summarizers.py
@@ -5,7 +5,7 @@
 from rail.core.utils import RAILDIR
 from rail.estimation.algos import naive_stack, point_est_hist, var_inf
 
-testdata = os.path.join(RAILDIR, "rail/examples_data/testdata/output_BPZ_lite.fits")
+testdata = os.path.join(RAILDIR, "rail/examples_data/testdata/output_BPZ_lite.hdf5")
 DS = RailStage.data_store