implementing mass pipeline for Y3

Author: vargatn

.gitmodules

@@ -1,6 +1,6 @@
 [submodule "submodules/xshear"]
 	path = submodules/xshear
-	url = https://github.com/danielgruen/xshear.git
+	url = https://github.com/vargatn/xshear.git
 	branch = redmapper
 
-# https://github.com/esheldon/xshear.git
+# https://github.com/esheldon/xshear.git

README.md

@@ -18,4 +18,4 @@ For detailed documentation see http://xpipe.readthedocs.io/
 * [API reference](https://xpipe.readthedocs.io/en/stable/api.html)
 
 
-Tamas Norbert Varga 2018
+Tamas Norbert Varga 2021

bin/redpipe/mkbins.py bin/pipe-y1/mkbins.py

@@ -0,0 +1,81 @@
+import pandas as pd
+import numpy as np
+import healpy as hp
+import h5py
+import matplotlib.pyplot as plt
+import yaml
+import fitsio
+
+
+main_file_path = "/e/ocean1/users/vargatn/DESY3/Y3_mastercat___UNBLIND___final_v1.1_12_22_20.h5"
+main_file = h5py.File(main_file_path, mode = 'r')
+fname_root = "/e/ocean1/users/vargatn/DES/DES_Y3A2_cluster/data/shearcat/"
+fname_prefix = "y3_mcal_sompz_v3_unblind"
+
+fmt = ['%d', '%d'] + 10 * ['%.6f']
+sheared_names = ["unsheared", "_1p", "_1m", "_2p", "_2m"]
+# sheared_names = ["unsheared", "sheared_1p", "sheared_1m", "sheared_2p", "sheared_2m"]
+select_names = ["select", "select_1p", "select_1m", "select_2p", "select_2m"]
+bin_names = ["bin1", "bin2", "bin3", "bin4"]
+
+for ibin, bin_name in enumerate(bin_names):
+    for ishear, sheared_name in enumerate(sheared_names):
+        if ishear == 0:
+            _sheared_name = ""
+            ztag = sheared_name
+            _name = ""
+        else:
+            _sheared_name = "_sheared" + sheared_name
+            ztag = "sheared" + sheared_name
+            _name = sheared_name
+
+        oname = fname_root + fname_prefix + "_" + bin_name + _name + ".dat"
+        print(oname)
+
+        index = main_file["index/"+select_names[ishear] + "_" + bin_name][:]
+
+        table = pd.DataFrame()
+
+        # loading part 1
+        _tmp = main_file["index/coadd_object_id"][:]
+        table["coadd_objects_id"] = _tmp[index]
+        print("loaded coadd_id")
+
+        _tmp = main_file["catalog/sompz/" + ztag + "/cell_wide"][:]
+        table["cell_id"] = _tmp[index]
+        print("loaded cell_id")
+
+        _tmp = main_file["catalog/gold/ra"][:]
+        table["ra"] = _tmp[index]
+        _tmp = main_file["catalog/gold/dec"][:]
+        table["dec"]= _tmp[index]
+        print("loaded ra, dec")
+
+        #         print("reading metacal unsheared")
+        _tmp = main_file["catalog/metacal/" + ztag + "/e_1"][:]
+        table["e1"] = _tmp[index]
+        _tmp = main_file["catalog/metacal/" + ztag + "/e_2"][:]
+        table["e2"] = _tmp[index]
+        print("loaded, e1, e2")
+        if ishear == 0:
+            _tmp = main_file["catalog/metacal/" + ztag + "/R11"][:]
+            table["R11"] = _tmp[index]
+            _tmp = main_file["catalog/metacal/" + ztag + "/R22"][:]
+            table["R22"] = _tmp[index]
+            print("loaded, R11, R22")
+            _R_12 = _tmp[index]
+            _tmp = main_file["catalog/metacal/" + ztag + "/R12"][:]
+            _R_21 = _tmp[index]
+            _tmp = main_file["catalog/metacal/" + ztag + "/R21"][:]
+            table["R12"] = (_R_12 + _R_21) / 2.
+            print("loaded, R12, R21")
+
+        _tmp = main_file["catalog/metacal/" + ztag + "/weight"][:]
+        table["weight"] = _tmp[index]
+
+        placeholder = np.ones(len(table)) * 1.5
+        table["z_mean"] = placeholder
+        table["z_mc"] = placeholder
+
+        table.to_hdf(oname.replace(".dat", ".h5"), key="data")
+        table.to_csv(oname, index=False, sep=" ", header=False)

bin/redpipe/postprocess.py bin/pipe-y1/postprocess.py

@@ -0,0 +1,33 @@
+"""
+Divides the lens and random points catalog into bins by parameters (and JK regions)
+"""
+
+import argparse
+import pickle
+import xpipe.paths as paths
+import xpipe.xhandle.parbins as parbins
+
+
+parser = argparse.ArgumentParser(description='Creates parameter bins for xshear calculations')
+parser.add_argument('--params', type=str)
+parser.add_argument('--norands', default=False, action="store_true")
+
+
+if __name__ == '__main__':
+    args = parser.parse_args()
+    paths.update_params(args.params)
+
+    lenses = parbins.prepare_lenses()
+    if not args.norands:
+        randoms = parbins.prepare_random()
+    else:
+        randoms = None
+
+    xio = parbins.XIO(lenses, randoms, force_centers=paths.params["njk_max"])
+    xio.mkdir()
+
+    logfile = xio.dpath + '/' + paths.params['tag'] + '_params.p'
+    pickle.dump(paths.params, open(logfile, 'wb'))
+
+    xio.loop_bins(norands=args.norands)
+

bin/tools/lwb_blind.py bin/pipe-y1/tools/lwb_blind.py

@@ -0,0 +1,52 @@
+"""
+Collects xshear output files into DeltaSigma profiles
+"""
+from __future__ import print_function, division
+import argparse
+import copy
+import numpy as np
+import pandas as pd
+import fitsio as fio
+
+import xpipe.paths as paths
+import xpipe.xhandle.parbins as parbins
+import xpipe.xhandle.xwrap as xwrap
+import xpipe.xhandle.shearops as shearops
+
+
+
+parser = argparse.ArgumentParser(description='postprocesses xshear output')
+parser.add_argument('--nometa', action="store_false", default=True)
+parser.add_argument('--npatch', default="auto")
+parser.add_argument('--params', type=str)
+parser.add_argument('--calibs', action="store_true", default=False)
+parser.add_argument('--norands', action="store_true", default=False)
+parser.add_argument('--lensweight', action="store_true", default=False)
+
+# TODO update with intuitive source bin detection, and check with a proper example
+if __name__ == '__main__':
+    args = parser.parse_args()
+    paths.update_params(args.params)
+
+    clusts = []
+    flist, flist_jk, rlist, rlist_jk = parbins.get_file_lists(paths.params, paths.dirpaths)
+    for i, clust_name in enumerate(flist_jk):
+        print("processing bin", i)
+
+        clust_infos = xwrap.create_infodict(clust_name)
+        clust_files = [info["outfile"] for info in clust_infos]
+
+        weights = None
+        if args.lensweight:
+            weights = shearops.load_weights(i)
+
+        bin_tag = clust_files[0].split("_" + paths.params["lens_prefix"])[1].split("_patch")[0]
+
+        clust = shearops.process_profile(clust_files, ismeta=args.nometa, weights=weights)
+
+        resroot = paths.dirpaths["results"] + "/" + \
+                  paths.params["tag"] + "/" + paths.params["tag"] + "_" + paths.params["lens_prefix"] + bin_tag
+        xwrap.write_profile(clust, resroot)
+
+        if args.calibs:
+            ccont = xwrap.append_scrit_inv(ccont, clust)

bin/tools/photoz/extract_full_PDF.py bin/pipe-y1/tools/photoz/extract_full_PDF.py

@@ -0,0 +1,125 @@
+"""
+Wrapper to run xshear with METACALIBRATION shears
+
+USAGE:
+run 5 times in total:
+ - first for shear responses and shear,
+ - 4 times to calculate inputs for selection responses
+
+"""
+from __future__ import print_function, division
+import argparse
+import numpy as np
+from collections import OrderedDict
+
+import xpipe.paths as paths
+import xpipe.xhandle.parbins as parbins
+import xpipe.xhandle.xwrap as xwrap
+import xpipe.tools.selector as selector
+
+parser = argparse.ArgumentParser(description='Runs xshear')
+parser.add_argument('--runall', action="store_true")
+parser.add_argument('--head', action="store_true")
+parser.add_argument('--nopairs', action="store_false")
+parser.add_argument('--noclust', action="store_true")
+parser.add_argument('--norands', action="store_true")
+parser.add_argument('--params', type=str)
+
+parser.add_argument("--metachunk", type=int, default=0)
+parser.add_argument("--nchunk", type=int, default=1)
+parser.add_argument("--ichunk", type=int, default=0)
+
+
+shape_path = paths.fullpaths[paths.params["shear_to_use"]]
+sheared_shape_paths = []
+for tag in xwrap.sheared_tags:
+    sheared_shape_paths.append(shape_path.replace(".dat", tag + ".dat"))
+
+if __name__ == '__main__':
+    args = parser.parse_args()
+    paths.update_params(args.params)
+
+    xpath = parbins.get_dpath(paths.params, paths.dirpaths) + "/" + \
+            paths.params["tag"] + "_xconfig.cfg"
+
+    flist, flist_jk, rlist, rlist_jk = parbins.get_file_lists(paths.params, paths.dirpaths)
+    print(flist)
+    print(rlist)
+    if args.metachunk == 0 or args.runall:
+        print(xpath)
+
+        shear_settings = OrderedDict([
+            ('shear_style', "metacal"),
+            ('sigmacrit_style', 'sample'),
+            ('shear_units', 'both'),
+            ('sourceid_style', 'somcell'),
+            ('weight_style', "uniform")
+        ])
+        xsettings=xwrap.get_main_source_settings()
+        xsettings.update({"shear": shear_settings})
+        # print(xsettings)
+        xwrap.write_custom_xconf(xpath, xsettings=xsettings)
+
+        if not args.noclust:
+            _alist = np.concatenate(flist_jk)
+            alist = selector.partition(_alist, args.nchunk)[args.ichunk]
+            clust_infos = xwrap.create_infodict(alist,
+                                                head=args.head,
+                                                pairs=args.nopairs,
+                                                src_bins=paths.params["source_bins_to_use"],
+                                                xconfig=xpath)
+            print(clust_infos[0])
+            xwrap.multi_xrun(clust_infos, nprocess=paths.params['nprocess'])
+
+        if not args.norands:
+            _blist = np.concatenate(rlist_jk)
+            blist = selector.partition(_blist, args.nchunk)[args.ichunk]
+            rands_infos = xwrap.create_infodict(blist,
+                                                head=args.head,
+                                                pairs=args.nopairs,
+                                                src_bins=paths.params["source_bins_to_use"],
+                                                xconfig=xpath)
+            xwrap.multi_xrun(rands_infos, nprocess=paths.params['nprocess'])
+
+    # this is the sheared computation with flags_select_*
+    #sheared_tags =["_sheared_1p", "_sheared_1m", "_sheared_2p", "_sheared_2m"]
+    for i, tag in enumerate(xwrap.sheared_tags):
+        if i == (args.metachunk - 1) or args.runall:
+            sheared_xconfig_fname = xpath.replace("xconfig", "xconfig" + tag)
+            print(sheared_xconfig_fname)
+
+            xsettings = {
+                "shear": OrderedDict([
+                    ('shear_style', "reduced"),
+                    ('sigmacrit_style', 'sample'),
+                    ('shear_units', 'both'),
+                    ('sourceid_style', 'somcell'),
+                ]),
+                "pairlog": xwrap.pairlog_nopairs
+            }
+            xwrap.write_custom_xconf(sheared_xconfig_fname, xsettings=xsettings)
+
+            if not args.noclust:
+                _alist = np.concatenate(flist_jk)
+                alist = selector.partition(_alist, args.nchunk)[args.ichunk]
+                sheared_clust_infos = xwrap.create_infodict(alist,
+                                                            head=args.head,
+                                                            pairs=args.nopairs,
+                                                            src_bins=paths.params["source_bins_to_use"],
+                                                            xconfig=sheared_xconfig_fname,
+                                                            metatag=tag,
+                                                            shape_path=sheared_shape_paths[i], )
+
+                xwrap.multi_xrun(sheared_clust_infos, nprocess=paths.params['nprocess'])
+
+            if not args.norands:
+                _blist = np.concatenate(rlist_jk)
+                blist = selector.partition(_blist, args.nchunk)[args.ichunk]
+                sheared_rand_infos = xwrap.create_infodict(blist,
+                                                           head=args.head,
+                                                           pairs=args.nopairs,
+                                                           src_bins=paths.params["source_bins_to_use"],
+                                                           xconfig=sheared_xconfig_fname,
+                                                           metatag=tag,
+                                                           shape_path=sheared_shape_paths[i], )
+                xwrap.multi_xrun(sheared_rand_infos, nprocess=paths.params['nprocess'])

bin/tools/photoz/extract_full_pwsum.py bin/pipe-y1/tools/photoz/extract_full_pwsum.py

@@ -2,10 +2,10 @@
 
 setup(name="xpipe",
       packages=find_packages(),
-      description="Measurement and data processing for DES Y1 cluster lensing",
+      description="Measurement and data processing for DES cluster lensing",
       install_requires=['numpy', 'scipy', 'pandas', 'astropy', 'fitsio', 'kmeans_radec', ],
       dependency_links=['http://github.com/user/repo/tarball/master#egg=package-1.0'],
       author="Tamas Norbert Varga",
       author_email="[email protected]",
-      version="0.45")
+      version="0.5")