From 6bab7ad5a4fed82287dd54ce74fb2547f384f1d5 Mon Sep 17 00:00:00 2001 From: Mateusz Mroz Date: Mon, 13 May 2024 17:40:49 +0200 Subject: [PATCH] color constrains put outside _parse_fit_constraints() dataset for color and blending constrain can be define by index or lable --- examples/example_16/ob03235_2_full.yaml | 10 +- examples/example_16/ulens_model_fit.py | 129 +++++++++++++++++------- 2 files changed, 96 insertions(+), 43 deletions(-) diff --git a/examples/example_16/ob03235_2_full.yaml b/examples/example_16/ob03235_2_full.yaml index ed82b83c..c5bb512c 100644 --- a/examples/example_16/ob03235_2_full.yaml +++ b/examples/example_16/ob03235_2_full.yaml @@ -38,13 +38,13 @@ fit_constraints: #negative_blending_flux_sigma_mag: 20. 1 3 # Alternative sharp constraint: # no_negative_blending_flux: True - #color constrains, where color=flux_S_dataset_no._k/flux_S_dataset_no._m: gauss mu sigma k m - color : gauss 25. 5. 1 2 + #color constrains, where color=flux_S_dataset_k/flux_S_dataset_m: gauss mu sigma k m + color : gauss 25. 5. "OGLE I-band" "OB03235_MOA.txt" #Alternative for binary source models: - #color source 1 : gauss 25. 5. 1 2 + #color source 1 : gauss 25. 5. "OGLE I-band" "OB03235_MOA.txt" #and/or - #color source 2 : gauss 25. 5. 1 2 - + #color source 2 : gauss 25. 5. 0 1 + prior: t_E: Mroz et al. 2017 # Other possibility: diff --git a/examples/example_16/ulens_model_fit.py b/examples/example_16/ulens_model_fit.py index aa437c4c..66eed00b 100644 --- a/examples/example_16/ulens_model_fit.py +++ b/examples/example_16/ulens_model_fit.py @@ -9,6 +9,7 @@ import warnings import math import numpy as np +import shlex from scipy.interpolate import interp1d from matplotlib import pyplot as plt from matplotlib import gridspec, rcParams, rcParamsDefault @@ -1470,49 +1471,66 @@ def _parse_fit_constraints(self): value = self._fit_constraints[key] if key == "negative_blending_flux_sigma_mag": + self._parse_fit_constraints_soft_blending(key, value) - if isinstance(value, (float, int)): - sigma = float(value) - sets = list(range(1, len(self._datasets) + 1)) - else: - sigma = float(value.split()[0]) - sets = list(map(int, value.split()[1:])) - if len(sets) > len(self._datasets): - raise ValueError( - 'dataset number specified in negative_blending_flux_sigma_mag do not match with provided datasets') + elif key in ['color', 'color source 1', 'color source 2']: + self._parse_fit_constraints_color(key, value) - self._fit_constraints[key] = [ - mm.Utils.get_flux_from_mag(sigma), sets] + if 'prior' in self._fit_constraints: + self._parse_fit_constraints_prior() - elif key in ['color', 'color source 1', 'color source 2']: + def _parse_fit_constraints_soft_blending(self, key, value): + """ + Check if soft fit constraint on blending flux are correctly defined. + """ + if isinstance(value, (float, int)): + sigma = float(value) + sets = list(range(1, len(self._datasets) + 1)) + else: + sigma = float(value.split()[0]) + sets = list(map(int, value.split()[1:])) + if len(sets) > len(self._datasets): + raise ValueError( + 'dataset number specified in negative_blending_flux_sigma_mag do not match with provided datasets') - words = value.split() - if len(words) != 5 or words[0] != 'gauss': - msg = "Something went wrong in parsing prior for " - msg += "{:}: {:}" - if len(words) == 3 and words[0] == 'gauss': - msg += ' color priors require the specification of datasets that should be used for color calculation ' - raise ValueError(msg.format(key, value)) - try: - settings = [words[0], float(words[1]), float( - words[2]), int(words[3]), int(words[4])] - except Exception: - raise ValueError('error in parsing: ' + words[1] + " " + - words[2] + " " + words[3] + " " + words[4]) - if settings[2] < 0.: - raise ValueError('sigma cannot be negative: ' + words[2]) - if settings[3] == settings[4]: - raise ValueError( - "in " + key + " fluxes have to be from different datasets") + self._fit_constraints[key] = [ + mm.Utils.get_flux_from_mag(sigma), sets] - if (0 >= settings[3] >= len(self._datasets)+1) or (0 >= settings[4] >= len(self._datasets)+1): - raise ValueError( - 'dataset number specified in color prior do not match with provided datasets') + def _parse_fit_constraints_color(self, key, value): + """ + Check if fit constraint on color are correctly defined. + """ + words = shlex.split(value, posix=False) + if len(words) != 5 or words[0] != 'gauss': + msg = "Something went wrong in parsing prior for " + msg += "{:}: {:}" + if len(words) == 3 and words[0] == 'gauss': + msg += ' color priors require the specification of datasets that should be used for color calculation ' + raise ValueError(msg.format(key, value)) + try: + settings = [words[0], float(words[1]), float( + words[2]), words[3], words[4]] + except Exception: + raise ValueError('error in parsing: ' + words[1] + " " + + words[2] + " " + words[3] + " " + words[4]) + if settings[2] < 0.: + raise ValueError('sigma cannot be negative: ' + words[2]) + if settings[3] == settings[4]: + raise ValueError( + "in " + key + " fluxes have to be from different datasets") - self._fit_constraints[key] = settings + if isinstance(settings[3], str) and isinstance(settings[4], str): + settings[3] = settings[3].strip('"') + settings[4] = settings[4].strip('"') + settings[3] = self._get_no_of_dataset_by_lable(settings[3]) + settings[4] = self._get_no_of_dataset_by_lable(settings[4]) - if 'prior' in self._fit_constraints: - self._parse_fit_constraints_prior() + if isinstance(settings[3], int) and isinstance(settings[4], int): + if (0 >= settings[3] >= len(self._datasets)-1) or (0 >= settings[4] >= len(self._datasets)-1): + raise ValueError( + 'dataset specified in color prior do not match with provided datasets') + + self._fit_constraints[key] = settings def _parse_fit_constraints_prior(self): """ @@ -1551,6 +1569,25 @@ def _parse_fit_constraints_prior(self): if len(priors) > 0: self._priors = priors + def _get_no_of_dataset_by_lable(self, label): + """ + Returns the index of a dataset with a specific label. + + Parameters + ---------- + label : str + Label of the dataset defined by MulensData.plot_properties['label']. + + Returns + ------- + int + Sequential index of the dataset [1,2,...,n_datasets] + """ + for (i, dataset) in enumerate(self._datasets): + if dataset.plot_properties['label'] == label: + return i + return -99 + def _read_prior_t_E_data(self): """ read data that specify t_E prior and parse them appropriately @@ -2165,10 +2202,26 @@ def _get_fluxes(self): return fluxes def _sumup_inside_prior(self, fluxes, key, inside, idx_plus): + """ + Calculates the contribution to the ln_prior from specified color constraints + Parameters + ---------- + fluxes : array + Array with fluxes of the current model. + key : str + inside : float + idx_plus : int + For a single source, idx_plus=0; + for a binary source, idx_plus=0 or 1. + + Returns + ------- + inside : float + """ settings = self._fit_constraints[key] - index1 = (settings[3]-1)*self._n_fluxes_per_dataset + idx_plus - index2 = (settings[4]-1)*self._n_fluxes_per_dataset + idx_plus + index1 = (settings[3])*self._n_fluxes_per_dataset + idx_plus + index2 = (settings[4])*self._n_fluxes_per_dataset + idx_plus value = fluxes[index1]/fluxes[index2] inside += self._get_ln_prior_for_1_parameter(value, settings[:-2])