From b50343edb4f7b4bb860897e1298ceed31b7c4e75 Mon Sep 17 00:00:00 2001 From: Naim Goksel Karacayli Date: Thu, 31 Oct 2024 17:16:29 -0400 Subject: [PATCH 01/13] function to read attributes file and modify weights of Deltas --- py/picca/utils.py | 32 +++++++++++++++++++++++++++++++- 1 file changed, 31 insertions(+), 1 deletion(-) diff --git a/py/picca/utils.py b/py/picca/utils.py index 95ce0579d..cb4ffb1e3 100644 --- a/py/picca/utils.py +++ b/py/picca/utils.py @@ -13,9 +13,10 @@ import sys import numpy as np import fitsio -import scipy.interpolate as interpolate import iminuit +from scipy import interpolate + def userprint(*args, **kwds): """Defines an extension of the print function. @@ -491,3 +492,32 @@ def unred(wave, ebv, R_V=3.1, LMC2=False, AVGLMC=False): corr = 1. / (10.**(0.4 * curve)) return corr + + +def modify_weights_with_varlss_factor(data, attributes, varlss_mod_factor): + """Recalculates the weights using eta and var_lss from attributes file. + Multiplies var_lss with varlss_mod_factor. + + Args: + data: dict(Delta) + A dictionary with the data. Changes in place. Keys are the healpix + numbers of each spectrum. Values are lists of delta instances. + attributes: str + Path to attributes file. + varlss_mod_factor: float + Multiplicative factor for var_lss + """ + varfunc = fitsio.read(attributes, ext="VAR_FUNC") + varfunc['lambda'] = np.log10(varfunc['lambda']) + interp_eta = interpolate.interp1d( + varfunc['lambda'], varfunc['eta'], kind='cubic') + interp_varlss = interpolate.interp1d( + varfunc['lambda'], varfunc['var_lss'], kind='cubic') + + for delta_list in data.values(): + for delta in delta_list: + if delta.ivar is None: + raise Exception("Need IVAR in deltas") + eta = interp_eta(delta.log_lambda) + var_lss = interp_varlss(delta.log_lambda) + delta.weights = delta.ivar / (eta + delta.ivar * varlss_mod_factor * var_lss) From 3bf3c10467488a5e01e011198c428cac2e1c5979 Mon Sep 17 00:00:00 2001 From: Naim Goksel Karacayli Date: Thu, 31 Oct 2024 17:16:49 -0400 Subject: [PATCH 02/13] read IVAR column if present in files --- py/picca/data.py | 10 ++++++++-- 1 file changed, 8 insertions(+), 2 deletions(-) diff --git a/py/picca/data.py b/py/picca/data.py index 46fa9e394..187d3b54c 100644 --- a/py/picca/data.py +++ b/py/picca/data.py @@ -432,7 +432,10 @@ def from_fitsio(cls, hdu, pk1d_type=False): weights = None cont = None else: - ivar = None + if 'IVAR' in hdu: + ivar = hdu['IVAR'][:].astype(float) + else: + ivar = None exposures_diff = None mean_snr = None mean_reso = None @@ -556,7 +559,10 @@ def from_image(cls, hdul, pk1d_type=False, z_min_qso=0, z_max_qso=10): else: raise KeyError("Did not find LOGLAM or LAMBDA in delta file") - ivar = Nones + if 'IVAR' in hdul: + ivar = hdul["IVAR"][:].astype(float) + else: + ivar = Nones exposures_diff = Nones mean_snr = Nones mean_reso = Nones From 3690cab4a574fda6d4ad7ed362136f2302465530 Mon Sep 17 00:00:00 2001 From: Naim Goksel Karacayli Date: Thu, 31 Oct 2024 17:17:29 -0400 Subject: [PATCH 03/13] two new arguments in picca_cf.py to read attributes file and modify weights --- bin/picca_cf.py | 14 ++++++++++++++ 1 file changed, 14 insertions(+) diff --git a/bin/picca_cf.py b/bin/picca_cf.py index dbf272727..038ddc9bb 100755 --- a/bin/picca_cf.py +++ b/bin/picca_cf.py @@ -231,6 +231,14 @@ def main(cmdargs): help='Rebin factor for deltas. If not None, deltas will ' 'be rebinned by that factor') + parser.add_argument('--varlss-mod-factor', type=float, default=None, + help='Modifies weights with this factor. Requires ' + 'IVAR column to be present in deltas and an input' + ' attributes file.') + parser.add_argument('--attributes', type=str, default=None, + help='Attributes file with VAR_FUNC extension with ' + 'lambda, eta, var_lss columns.') + args = parser.parse_args(cmdargs) if args.nproc is None: @@ -276,6 +284,12 @@ def main(cmdargs): rebin_factor=args.rebin_factor, z_min_qso=args.z_min_sources, z_max_qso=args.z_max_sources,) + + if args.varlss_mod_factor is not None: + assert args.attributes + utils.modify_weights_with_varlss_factor( + data, args.attributes, args.varlss_mod_factor) + del z_max cf.data = data cf.num_data = num_data From 0283ff0516597a58c5397238433ba2f9cadf0b2a Mon Sep 17 00:00:00 2001 From: Naim Goksel Karacayli Date: Fri, 1 Nov 2024 09:51:00 -0400 Subject: [PATCH 04/13] fix reading IVAR --- py/picca/data.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/py/picca/data.py b/py/picca/data.py index 187d3b54c..f38e1e32a 100644 --- a/py/picca/data.py +++ b/py/picca/data.py @@ -432,7 +432,7 @@ def from_fitsio(cls, hdu, pk1d_type=False): weights = None cont = None else: - if 'IVAR' in hdu: + if 'IVAR' in hdu.get_colnames(): ivar = hdu['IVAR'][:].astype(float) else: ivar = None From e89bc84c00f4b3aff26b8933baf0cc6e53b7795c Mon Sep 17 00:00:00 2001 From: Naim Goksel Karacayli Date: Fri, 1 Nov 2024 11:11:06 -0400 Subject: [PATCH 05/13] move rebin factor outside of multiprocess. Check with WAVE_SOLUTION is removed. Alsre recalculate weights before rebinning --- py/picca/io.py | 67 ++++++++++++++++++++++---------------------------- 1 file changed, 29 insertions(+), 38 deletions(-) diff --git a/py/picca/io.py b/py/picca/io.py index 8789483d1..e587fc4b6 100644 --- a/py/picca/io.py +++ b/py/picca/io.py @@ -21,7 +21,7 @@ import warnings from multiprocessing import Pool -from .utils import userprint +from .utils import userprint, modify_weights_with_varlss_factor from .data import Delta, QSO from .pk1d.prep_pk1d import exp_diff, spectral_resolution from .pk1d.prep_pk1d import spectral_resolution_desi @@ -245,7 +245,7 @@ def read_blinding(in_dir): return blinding -def read_delta_file(filename, z_min_qso=0, z_max_qso=10, rebin_factor=None): +def read_delta_file(filename, z_min_qso=0, z_max_qso=10): """Extracts deltas from a single file. Args: filename: str @@ -254,8 +254,6 @@ def read_delta_file(filename, z_min_qso=0, z_max_qso=10, rebin_factor=None): Specifies the minimum redshift for QSOs z_max_qso: float - default: 10 Specifies the maximum redshift for QSOs - rebin_factor: int - default: None - Factor to rebin the lambda grid by. If None, no rebinning is done. Returns: deltas: A dictionary with the data. Keys are the healpix numbers of each @@ -269,39 +267,16 @@ def read_delta_file(filename, z_min_qso=0, z_max_qso=10, rebin_factor=None): else: deltas = [Delta.from_fitsio(hdu) for hdu in hdul[1:] if z_min_qso 5 and in_dir[-5:] == '.fits': files += sorted(glob.glob(in_dir)) else: - files += sorted(glob.glob(in_dir + '/*.fits') + glob.glob(in_dir + - '/*.fits.gz')) + files += sorted( + glob.glob(in_dir + '/*.fits') + glob.glob(in_dir + '/*.fits.gz')) files = sorted(files) - if rebin_factor is not None: - userprint(f"Rebinning deltas by a factor of {rebin_factor}\n") - - arguments = [(f, z_min_qso, z_max_qso, rebin_factor) for f in files] + arguments = [(f, z_min_qso, z_max_qso) for f in files] pool = Pool(processes=nproc) results = pool.starmap(read_delta_file, arguments) pool.close() @@ -384,6 +362,19 @@ def read_deltas(in_dir, userprint("\n") + # Recalculate weights + if varlss_mod_factor is not None: + userprint(f"Recalculating weights with eta_lss={varlss_mod_factor}\n") + if not attributes: + raise ValueError('Recalculating weights require an attributes file.') + modify_weights_with_varlss_factor(deltas, attributes, varlss_mod_factor) + + # Rebin + if rebin_factor is not None: + userprint(f"Rebinning deltas by a factor of {rebin_factor}\n") + for delta_ in deltas: + delta_.rebin(rebin_factor) + # compute healpix numbers phi = [delta.ra for delta in deltas] theta = [np.pi / 2. - delta.dec for delta in deltas] From c4db61402853df8d6568816dbcdd755739931ef2 Mon Sep 17 00:00:00 2001 From: Naim Goksel Karacayli Date: Fri, 1 Nov 2024 11:11:54 -0400 Subject: [PATCH 06/13] estimate dwave from existing wave grid. Weakness is pathological input wavelength with gaps every other pixel. --- py/picca/data.py | 10 +++++++--- 1 file changed, 7 insertions(+), 3 deletions(-) diff --git a/py/picca/data.py b/py/picca/data.py index f38e1e32a..c870afd26 100644 --- a/py/picca/data.py +++ b/py/picca/data.py @@ -649,16 +649,20 @@ def project(self): self.delta -= mean_delta + res - def rebin(self, factor, dwave=0.8): + def rebin(self, factor): """Rebin deltas by an integer factor Args: factor: int Factor to rebin deltas (new_bin_size = factor * old_bin_size) - dwave: float - Delta lambda of original deltas """ wave = 10**np.array(self.log_lambda) + dwaves = np.round(np.diff(wave), decimals=3) + dwave = np.min(dwaves) + + # This check would not work if there are gaps. + # if not np.allclose(dwave, dwaves): + # raise ValueError('Delta rebinning only implemented for linear lambda bins') start = wave.min() - dwave / 2 num_bins = np.ceil(((wave[-1] - wave[0]) / dwave + 1) / factor) From 4da0f6551821166d07aa169ff409c2fb4cd95b07 Mon Sep 17 00:00:00 2001 From: Naim Goksel Karacayli Date: Fri, 1 Nov 2024 11:12:23 -0400 Subject: [PATCH 07/13] recalculating weights moved to reading deltas function --- bin/picca_cf.py | 28 +++++++++------------------- py/picca/utils.py | 15 +++++++-------- 2 files changed, 16 insertions(+), 27 deletions(-) diff --git a/bin/picca_cf.py b/bin/picca_cf.py index 038ddc9bb..467f526bc 100755 --- a/bin/picca_cf.py +++ b/bin/picca_cf.py @@ -272,24 +272,12 @@ def main(cmdargs): ### Read data 1 print('z_min_sources', args.z_min_sources) - data, num_data, z_min, z_max = io.read_deltas(args.in_dir, - cf.nside, - cf.lambda_abs, - cf.alpha, - cf.z_ref, - cosmo, - max_num_spec=args.nspec, - no_project=args.no_project, - nproc=args.nproc, - rebin_factor=args.rebin_factor, - z_min_qso=args.z_min_sources, - z_max_qso=args.z_max_sources,) - - if args.varlss_mod_factor is not None: - assert args.attributes - utils.modify_weights_with_varlss_factor( - data, args.attributes, args.varlss_mod_factor) - + data, num_data, z_min, z_max = io.read_deltas( + args.in_dir, cf.nside, cf.lambda_abs, cf.alpha, cf.z_ref, cosmo, + max_num_spec=args.nspec, no_project=args.no_project, nproc=args.nproc, + rebin_factor=args.rebin_factor, z_min_qso=args.z_min_sources, + z_max_qso=args.z_max_sources, varlss_mod_factor=args.varlss_mod_factor, + attributes=args.attributes) del z_max cf.data = data cf.num_data = num_data @@ -321,7 +309,9 @@ def main(cmdargs): nproc=args.nproc, rebin_factor=args.rebin_factor, z_min_qso=args.z_min_sources, - z_max_qso=args.z_max_sources,) + z_max_qso=args.z_max_sources, + varlss_mod_factor=args.varlss_mod_factor, + attributes=args.attributes) del z_max2 cf.data2 = data2 cf.num_data2 = num_data2 diff --git a/py/picca/utils.py b/py/picca/utils.py index cb4ffb1e3..3e4967a51 100644 --- a/py/picca/utils.py +++ b/py/picca/utils.py @@ -499,7 +499,7 @@ def modify_weights_with_varlss_factor(data, attributes, varlss_mod_factor): Multiplies var_lss with varlss_mod_factor. Args: - data: dict(Delta) + data: list(Delta) A dictionary with the data. Changes in place. Keys are the healpix numbers of each spectrum. Values are lists of delta instances. attributes: str @@ -514,10 +514,9 @@ def modify_weights_with_varlss_factor(data, attributes, varlss_mod_factor): interp_varlss = interpolate.interp1d( varfunc['lambda'], varfunc['var_lss'], kind='cubic') - for delta_list in data.values(): - for delta in delta_list: - if delta.ivar is None: - raise Exception("Need IVAR in deltas") - eta = interp_eta(delta.log_lambda) - var_lss = interp_varlss(delta.log_lambda) - delta.weights = delta.ivar / (eta + delta.ivar * varlss_mod_factor * var_lss) + for delta in data: + if delta.ivar is None: + raise Exception("Need IVAR in deltas") + eta = interp_eta(delta.log_lambda) + var_lss = varlss_mod_factor * interp_varlss(delta.log_lambda) + delta.weights = delta.ivar / (eta + delta.ivar * var_lss) From 6cc67fd51a9d55989cebb842d2ba4f98e6cc77b0 Mon Sep 17 00:00:00 2001 From: Naim Goksel Karacayli Date: Fri, 1 Nov 2024 13:38:37 -0400 Subject: [PATCH 08/13] fix: eta and varlss need to be extrapolated --- py/picca/utils.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/py/picca/utils.py b/py/picca/utils.py index 3e4967a51..6463bdd9d 100644 --- a/py/picca/utils.py +++ b/py/picca/utils.py @@ -510,9 +510,9 @@ def modify_weights_with_varlss_factor(data, attributes, varlss_mod_factor): varfunc = fitsio.read(attributes, ext="VAR_FUNC") varfunc['lambda'] = np.log10(varfunc['lambda']) interp_eta = interpolate.interp1d( - varfunc['lambda'], varfunc['eta'], kind='cubic') + varfunc['lambda'], varfunc['eta'], fill_value='extrapolate', kind='cubic') interp_varlss = interpolate.interp1d( - varfunc['lambda'], varfunc['var_lss'], kind='cubic') + varfunc['lambda'], varfunc['var_lss'], fill_value='extrapolate', kind='cubic') for delta in data: if delta.ivar is None: From 1b27e2b620ce24a1ff3feacc67bf0e65ae939523 Mon Sep 17 00:00:00 2001 From: Naim Goksel Karacayli Date: Mon, 4 Nov 2024 13:55:14 -0500 Subject: [PATCH 09/13] function to stack and renormalize deltas --- bin/picca_cf.py | 7 +++++-- py/picca/io.py | 13 +++++++++++-- py/picca/utils.py | 25 +++++++++++++++++++++++++ 3 files changed, 41 insertions(+), 4 deletions(-) diff --git a/bin/picca_cf.py b/bin/picca_cf.py index 467f526bc..f48408948 100755 --- a/bin/picca_cf.py +++ b/bin/picca_cf.py @@ -238,6 +238,8 @@ def main(cmdargs): parser.add_argument('--attributes', type=str, default=None, help='Attributes file with VAR_FUNC extension with ' 'lambda, eta, var_lss columns.') + parser.add_argument('--renormalize-deltas', action="store_true", + help="Stacks deltas and renormalizes deltas.") args = parser.parse_args(cmdargs) @@ -277,7 +279,7 @@ def main(cmdargs): max_num_spec=args.nspec, no_project=args.no_project, nproc=args.nproc, rebin_factor=args.rebin_factor, z_min_qso=args.z_min_sources, z_max_qso=args.z_max_sources, varlss_mod_factor=args.varlss_mod_factor, - attributes=args.attributes) + attributes=args.attributes, renormalize_deltas=args.renormalize_deltas) del z_max cf.data = data cf.num_data = num_data @@ -311,7 +313,8 @@ def main(cmdargs): z_min_qso=args.z_min_sources, z_max_qso=args.z_max_sources, varlss_mod_factor=args.varlss_mod_factor, - attributes=args.attributes) + attributes=args.attributes, + renormalize_deltas=args.renormalize_deltas) del z_max2 cf.data2 = data2 cf.num_data2 = num_data2 diff --git a/py/picca/io.py b/py/picca/io.py index e587fc4b6..7650e6c16 100644 --- a/py/picca/io.py +++ b/py/picca/io.py @@ -21,7 +21,8 @@ import warnings from multiprocessing import Pool -from .utils import userprint, modify_weights_with_varlss_factor +from .utils import ( + userprint, modify_weights_with_varlss_factor, stack_flux_and_normalize) from .data import Delta, QSO from .pk1d.prep_pk1d import exp_diff, spectral_resolution from .pk1d.prep_pk1d import spectral_resolution_desi @@ -275,7 +276,8 @@ def read_delta_file(filename, z_min_qso=0, z_max_qso=10): def read_deltas( in_dir, nside, lambda_abs, alpha, z_ref, cosmo, max_num_spec=None, no_project=False, nproc=None, rebin_factor=None, - z_min_qso=0, z_max_qso=10, varlss_mod_factor=None, attributes=None + z_min_qso=0, z_max_qso=10, varlss_mod_factor=None, attributes=None, + renormalize_deltas=False ): """Reads deltas and computes their redshifts. @@ -317,6 +319,8 @@ def read_deltas( attributes: str - default: None Attributes file with VAR_FUNC extension with lambda, eta, var_lss columns. + renormalize_deltas: bool - default: False + If True, stack 1 + delta and renormalize after rebinning. Returns: The following variables: @@ -402,6 +406,11 @@ def read_deltas( data[healpix] = [] data[healpix].append(delta) + if renormalize_deltas: + stack_flux_and_normalize( + deltas, lambda_abs * (1 + z_min) + 4.0, lambda_abs * (1 + z_max) - 4.0, + dlambda=8.0) + return data, num_data, z_min, z_max diff --git a/py/picca/utils.py b/py/picca/utils.py index 6463bdd9d..345644964 100644 --- a/py/picca/utils.py +++ b/py/picca/utils.py @@ -520,3 +520,28 @@ def modify_weights_with_varlss_factor(data, attributes, varlss_mod_factor): eta = interp_eta(delta.log_lambda) var_lss = varlss_mod_factor * interp_varlss(delta.log_lambda) delta.weights = delta.ivar / (eta + delta.ivar * var_lss) + + +def stack_flux_and_normalize(data, l1, l2, dlambda): + n = round((l2 - l1) / dlambda + 1.1) + stack_delta = np.zeros(n) + stack_weight = np.zeros(n) + + for delta in data: + wave = 10**delta.log_lambda + bins = ((wave - l1) / dlambda + 0.5).astype(np.int64) + stack_delta += np.bincount( + bins, weights=delta.delta * delta.weights, minlength=n) + stack_weight += np.bincount(bins, weights=delta.weights, minlength=n) + + w = stack_weight != 0 + stack_delta = stack_delta[w] + stack_weight = stack_weight[w] + stack_delta /= stack_weight + wave = np.linspace(l1, l2, n)[w] + spl = interpolate.interp1d( + wave, stack_delta + 1, fill_value='extrapolate', kind='nearest') + + for delta in data: + wave = 10**delta.log_lambda + delta.delta = (1 + delta.delta) / spl(wave) - 1 From ea5564259cf38e333b2ea2c89ad2d7ea858317fd Mon Sep 17 00:00:00 2001 From: "wynneturner04@gmail.com" Date: Mon, 13 Jan 2025 15:50:16 -0800 Subject: [PATCH 10/13] Stack and renormalize deltas before rebinning --- py/picca/io.py | 19 +++++++++++-------- py/picca/utils.py | 5 ++++- 2 files changed, 15 insertions(+), 9 deletions(-) diff --git a/py/picca/io.py b/py/picca/io.py index 7650e6c16..1229f7272 100644 --- a/py/picca/io.py +++ b/py/picca/io.py @@ -229,8 +229,8 @@ def read_blinding(in_dir): elif len(in_dir) > 5 and in_dir[-5:] == '.fits': files += glob.glob(in_dir) else: - files += glob.glob(in_dir + '/*.fits') + glob.glob(in_dir - + '/*.fits.gz') + files += glob.glob(in_dir + '/delta*.fits') + glob.glob(in_dir + + '/delta*.fits.gz') filename = files[0] hdul = fitsio.FITS(filename) if "LAMBDA" in hdul: # This is for ImageHDU format @@ -342,7 +342,7 @@ def read_deltas( files += sorted(glob.glob(in_dir)) else: files += sorted( - glob.glob(in_dir + '/*.fits') + glob.glob(in_dir + '/*.fits.gz')) + glob.glob(in_dir + '/delta*.fits') + glob.glob(in_dir + '/delta*.fits.gz')) files = sorted(files) arguments = [(f, z_min_qso, z_max_qso) for f in files] @@ -372,6 +372,14 @@ def read_deltas( if not attributes: raise ValueError('Recalculating weights require an attributes file.') modify_weights_with_varlss_factor(deltas, attributes, varlss_mod_factor) + + # Renormalize deltas + if renormalize_deltas: + userprint(f"Stacking and renormalizing deltas\n") + z_min, z_max = 0, 10 + stack_flux_and_normalize( + deltas, lambda_abs * (1 + z_min) + 4.0, lambda_abs * (1 + z_max) - 4.0, + dlambda=8.0) # Rebin if rebin_factor is not None: @@ -406,11 +414,6 @@ def read_deltas( data[healpix] = [] data[healpix].append(delta) - if renormalize_deltas: - stack_flux_and_normalize( - deltas, lambda_abs * (1 + z_min) + 4.0, lambda_abs * (1 + z_max) - 4.0, - dlambda=8.0) - return data, num_data, z_min, z_max diff --git a/py/picca/utils.py b/py/picca/utils.py index 345644964..a7b063f00 100644 --- a/py/picca/utils.py +++ b/py/picca/utils.py @@ -507,7 +507,10 @@ def modify_weights_with_varlss_factor(data, attributes, varlss_mod_factor): varlss_mod_factor: float Multiplicative factor for var_lss """ - varfunc = fitsio.read(attributes, ext="VAR_FUNC") + try: + varfunc = fitsio.read(attributes, ext="VAR_FUNC") + except OSError: + varfunc = fitsio.read(attributes, ext="VAR_FUNC-fit") varfunc['lambda'] = np.log10(varfunc['lambda']) interp_eta = interpolate.interp1d( varfunc['lambda'], varfunc['eta'], fill_value='extrapolate', kind='cubic') From 6e9ec383859fa296022841055f8bab9d3fa02056 Mon Sep 17 00:00:00 2001 From: "wynneturner04@gmail.com" Date: Fri, 17 Jan 2025 09:24:18 -0800 Subject: [PATCH 11/13] Added functionality to xcf --- bin/picca_xcf.py | 15 ++++++++++++++- 1 file changed, 14 insertions(+), 1 deletion(-) diff --git a/bin/picca_xcf.py b/bin/picca_xcf.py index 63c693356..07e17c646 100755 --- a/bin/picca_xcf.py +++ b/bin/picca_xcf.py @@ -252,6 +252,16 @@ def main(cmdargs): help='Rebin factor for deltas. If not None, deltas will ' 'be rebinned by that factor') + parser.add_argument('--varlss-mod-factor', type=float, default=None, + help='Modifies weights with this factor. Requires ' + 'IVAR column to be present in deltas and an input' + ' attributes file.') + parser.add_argument('--attributes', type=str, default=None, + help='Attributes file with VAR_FUNC extension with ' + 'lambda, eta, var_lss columns.') + parser.add_argument('--renormalize-deltas', action="store_true", + help="Stacks deltas and renormalizes deltas.") + args = parser.parse_args(cmdargs) if args.nproc is None: @@ -310,7 +320,10 @@ def main(cmdargs): nproc=args.nproc, rebin_factor=args.rebin_factor, z_min_qso=args.z_min_sources, - z_max_qso=args.z_max_sources) + z_max_qso=args.z_max_sources, + varlss_mod_factor=args.varlss_mod_factor, + attributes=args.attributes, + renormalize_deltas=args.renormalize_deltas) xcf.data = data xcf.num_data = num_data userprint("") From 16d9a251c38d03df39b8838d9eef269fe0edbba6 Mon Sep 17 00:00:00 2001 From: "wynneturner04@gmail.com" Date: Wed, 22 Jan 2025 15:55:33 -0800 Subject: [PATCH 12/13] fix for lyb region and any attributes file --- bin/picca_cf.py | 6 +++++- py/picca/utils.py | 15 ++++++++++----- 2 files changed, 15 insertions(+), 6 deletions(-) diff --git a/bin/picca_cf.py b/bin/picca_cf.py index f48408948..bacf8278b 100755 --- a/bin/picca_cf.py +++ b/bin/picca_cf.py @@ -238,6 +238,8 @@ def main(cmdargs): parser.add_argument('--attributes', type=str, default=None, help='Attributes file with VAR_FUNC extension with ' 'lambda, eta, var_lss columns.') + parser.add_argument('--attributes2', type=str, default=None, required=False, + help='Attributes file for 2nd tracer') parser.add_argument('--renormalize-deltas', action="store_true", help="Stacks deltas and renormalizes deltas.") @@ -294,6 +296,8 @@ def main(cmdargs): cf.alpha2 = args.z_evol2 if args.in_dir2 is None: args.in_dir2 = args.in_dir + if args.attributes2 is None: + args.attributes2 = args.attributes if args.lambda_abs2: cf.lambda_abs2 = constants.ABSORBER_IGM[args.lambda_abs2] else: @@ -313,7 +317,7 @@ def main(cmdargs): z_min_qso=args.z_min_sources, z_max_qso=args.z_max_sources, varlss_mod_factor=args.varlss_mod_factor, - attributes=args.attributes, + attributes=args.attributes2, renormalize_deltas=args.renormalize_deltas) del z_max2 cf.data2 = data2 diff --git a/py/picca/utils.py b/py/picca/utils.py index a7b063f00..d2609b573 100644 --- a/py/picca/utils.py +++ b/py/picca/utils.py @@ -13,6 +13,7 @@ import sys import numpy as np import fitsio +from astropy.io import fits import iminuit from scipy import interpolate @@ -507,15 +508,19 @@ def modify_weights_with_varlss_factor(data, attributes, varlss_mod_factor): varlss_mod_factor: float Multiplicative factor for var_lss """ + att = fits.open(attributes) # need astropy's case-insensitivity to work with qsonic and picca attributes files try: - varfunc = fitsio.read(attributes, ext="VAR_FUNC") + varfunc = att['VAR_FUNC'] except OSError: - varfunc = fitsio.read(attributes, ext="VAR_FUNC-fit") - varfunc['lambda'] = np.log10(varfunc['lambda']) + varfunc = att['VAR_FUNC-fit'] + try: + loglam = np.log10(varfunc.data['lambda']) + except KeyError: + loglam = varfunc.data['loglam'] interp_eta = interpolate.interp1d( - varfunc['lambda'], varfunc['eta'], fill_value='extrapolate', kind='cubic') + loglam, varfunc.data['eta'], fill_value='extrapolate', kind='cubic') interp_varlss = interpolate.interp1d( - varfunc['lambda'], varfunc['var_lss'], fill_value='extrapolate', kind='cubic') + loglam, varfunc.data['var_lss'], fill_value='extrapolate', kind='cubic') for delta in data: if delta.ivar is None: From 8b0a63ee42bac0d05d97520ab16937e095cf4aac Mon Sep 17 00:00:00 2001 From: "wynneturner04@gmail.com" Date: Wed, 22 Jan 2025 19:45:59 -0800 Subject: [PATCH 13/13] added option to update weights in dmat, xdmat --- bin/picca_dmat.py | 24 ++++++++++++++++++++++-- bin/picca_xdmat.py | 16 +++++++++++++++- 2 files changed, 37 insertions(+), 3 deletions(-) diff --git a/bin/picca_dmat.py b/bin/picca_dmat.py index 959d6bf37..8685d9391 100755 --- a/bin/picca_dmat.py +++ b/bin/picca_dmat.py @@ -249,6 +249,18 @@ def main(cmdargs): action='store_true', help='Ignore redshift evolution when computing distortion matrix') + parser.add_argument('--varlss-mod-factor', type=float, default=None, + help='Modifies weights with this factor. Requires ' + 'IVAR column to be present in deltas and an input' + ' attributes file.') + parser.add_argument('--attributes', type=str, default=None, + help='Attributes file with VAR_FUNC extension with ' + 'lambda, eta, var_lss columns.') + parser.add_argument('--attributes2', type=str, default=None, required=False, + help='Attributes file for 2nd tracer') + parser.add_argument('--renormalize-deltas', action="store_true", + help="Stacks deltas and renormalizes deltas.") + args = parser.parse_args(cmdargs) if args.nproc is None: @@ -310,7 +322,10 @@ def main(cmdargs): nproc=args.nproc, rebin_factor=args.rebin_factor, z_min_qso=args.z_min_sources, - z_max_qso=args.z_max_sources) + z_max_qso=args.z_max_sources, + varlss_mod_factor=args.varlss_mod_factor, + attributes=args.attributes, + renormalize_deltas=args.renormalize_deltas) del z_max cf.data = data cf.num_data = num_data @@ -325,6 +340,8 @@ def main(cmdargs): cf.alpha2 = args.z_evol2 if args.in_dir2 is None: args.in_dir2 = args.in_dir + if args.attributes2 is None: + args.attributes2 = args.attributes if args.lambda_abs2: cf.lambda_abs2 = constants.ABSORBER_IGM[args.lambda_abs2] else: @@ -341,7 +358,10 @@ def main(cmdargs): nproc=args.nproc, rebin_factor=args.rebin_factor, z_min_qso=args.z_min_sources, - z_max_qso=args.z_max_sources) + z_max_qso=args.z_max_sources, + varlss_mod_factor=args.varlss_mod_factor, + attributes=args.attributes2, + renormalize_deltas=args.renormalize_deltas) del z_max2 cf.data2 = data2 cf.num_data2 = num_data2 diff --git a/bin/picca_xdmat.py b/bin/picca_xdmat.py index 1da22d49b..b1daf5787 100755 --- a/bin/picca_xdmat.py +++ b/bin/picca_xdmat.py @@ -248,6 +248,17 @@ def main(cmdargs): action='store_true', help='Ignore redshift evolution when computing distortion matrix') + + parser.add_argument('--varlss-mod-factor', type=float, default=None, + help='Modifies weights with this factor. Requires ' + 'IVAR column to be present in deltas and an input' + ' attributes file.') + parser.add_argument('--attributes', type=str, default=None, + help='Attributes file with VAR_FUNC extension with ' + 'lambda, eta, var_lss columns.') + parser.add_argument('--renormalize-deltas', action="store_true", + help="Stacks deltas and renormalizes deltas.") + args = parser.parse_args(cmdargs) if args.nproc is None: args.nproc = cpu_count() // 2 @@ -303,7 +314,10 @@ def main(cmdargs): nproc=args.nproc, rebin_factor=args.rebin_factor, z_min_qso=args.z_min_sources, - z_max_qso=args.z_max_sources) + z_max_qso=args.z_max_sources, + varlss_mod_factor=args.varlss_mod_factor, + attributes=args.attributes, + renormalize_deltas=args.renormalize_deltas) xcf.data = data xcf.num_data = num_data userprint("\n")