From a5354ab1f099e8e5ec876378e94a8044e6ab853c Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Andr=C3=A9=20Z=2E=20Vitorelli?= Date: Thu, 11 Aug 2022 16:11:30 +0200 Subject: [PATCH 1/6] adding simple galaxies dataset --- autometacal/python/datasets/simple.py | 134 ++++++++++++++++++++++++++ 1 file changed, 134 insertions(+) create mode 100644 autometacal/python/datasets/simple.py diff --git a/autometacal/python/datasets/simple.py b/autometacal/python/datasets/simple.py new file mode 100644 index 0000000..25a312e --- /dev/null +++ b/autometacal/python/datasets/simple.py @@ -0,0 +1,134 @@ +"""simple dataset.""" +import tensorflow_datasets as tfds +import tensorflow as tf +import numpy as np + +import galsim + +_DESCRIPTION = "This tfds generates random toy-model round galaxy stamps." +_CITATION = "{NEEDED}" +_URL = "https://github.com/CosmoStat/autometacal" + +class SimpleConfig(tfds.core.BuilderConfig): + """BuilderConfig for GalGen.""" + + def __init__(self, *,dataset_size, stamp_size, + ,**kwargs): + """BuilderConfig for Simple. + Args: + dataset_size: number of stamps + stamp_size: image stamp size in pixels. + + """ + v2 = tfds.core.Version("1.0.0") + super(SimpleConfig, self).__init__( + description=( + "%d stamps of %d x %d pixels, 0.263 arcsec/pix, flux = 1. " % + (dataset_size,stamp_size, stamp_size)), + version=v2, + **kwargs) + self.dataset_size = dataset_size + self.stamp_size = stamp_size + +class Simple(tfds.core.GeneratorBasedBuilder): + """Random galaxy image generator.""" + + MANUAL_DOWNLOAD_INSTRUCTIONS = """\ + Nothing to download. DataSet is generated at first call. + """ + + BUILDER_CONFIGS = [ + SimpleConfig(name='small_1k_noshear', dataset_size=1000, stamp_size=51), + SimpleConfig(name='large_1k_noshear', dataset_size=1000, stamp_size=101), + ] + + VERSION = tfds.core.Version('0.5.0') + RELEASE_NOTES = {'0.5.0': "Updated functionalities, simpler."} + + def _info(self): + return tfds.core.DatasetInfo( + builder=self, + # Description and homepage used for documentation + description=_DESCRIPTION, + homepage=_URL, + features=tfds.features.FeaturesDict({ + 'label': tfds.features.Tensor(shape=[2], dtype=tf.float32), + 'gal_model': tfds.features.Tensor( + shape=[self.builder_config.stamp_size,self.builder_config.stamp_size], + dtype=tf.float32 + ), + 'obs': tfds.features.Tensor( + shape=[self.builder_config.stamp_size,self.builder_config.stamp_size], + dtype=tf.float32 + ), + 'psf_image': tfds.features.Tensor( + shape=[self.builder_config.stamp_size, self.builder_config.stamp_size], + dtype=tf.float32 + ), + 'noise': tfds.features.Tensor( + shape=[self.builder_config.stamp_size, self.builder_config.stamp_size], + dtype=tf.float32 + ), + }), + supervised_keys=("obs"), + citation=_CITATION + ) + + def _split_generators(self,dl): + """Returns generators according to split.""" + return {tfds.Split.TRAIN: self._generate_examples( + self.builder_config.dataset_size, + self.builder_config.stamp_size + )} + + def _generate_examples(self, dataset_size, stamp_size,g1,g2): + """Yields examples.""" + rng = np.random.RandomState(31415) + + for i in range(dataset_size): + + #generate example + model, obs_img, psf_img, noise_img = make_data(rng,stamp_size=stamp_size) + yield '%d'%i, {'gal_model': model, #noiseless PSFless galaxy model + 'obs': obs_img, #observed image + 'psf_image': psf_img, #psf image + 'noise' : noise_img, + 'label': np.array([g1,g2],dtype='float32')} + +def make_data(rng,stamp_size = 51): + """Simple exponetial profile toy model galaxy""" + #values from the ngmix example.py + scale = 0.263 + psf_fwhm = 0.9 + gal_hlr = 0.5 + noise = 1e-6 + + psf = galsim.Moffat(beta=4.8,fwhm=psf_fwhm) + + obj = galsim.Exponential(half_light_radius=gal_hlr) + obs = galsim.Convolve([psf, obj]) + + #image drawing method, no_pixel or auto + method = 'auto' + + psf_image = psf.drawImage( + nx=stamp_size, + ny=stamp_size, + scale=scale, + method=method).array #psf image + gal_image = obj.drawImage( + nx=stamp_size, + ny=stamp_size, + scale=scale, + method=method).array #model image + obs_image = obs.drawImage( + nx=stamp_size, + ny=stamp_size, + scale=scale, + method=method).array #observed image, prenoise + + noise_image = rng.normal(scale=noise, size=obs_image.shape) + obs_image += noise_image + + return gal_image, obs_image, psf_image, noise_image.astype('float32') + From 1b112f357ca286a1112b8d9725b22fa8b3164e91 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Andr=C3=A9=20Z=2E=20Vitorelli?= Date: Thu, 11 Aug 2022 16:20:49 +0200 Subject: [PATCH 2/6] galgen reestructure --- autometacal/python/datasets/galgen.py | 181 ++++++++++++++++++++++++++ 1 file changed, 181 insertions(+) create mode 100644 autometacal/python/datasets/galgen.py diff --git a/autometacal/python/datasets/galgen.py b/autometacal/python/datasets/galgen.py new file mode 100644 index 0000000..911760e --- /dev/null +++ b/autometacal/python/datasets/galgen.py @@ -0,0 +1,181 @@ +"""GalGen dataset.""" +import tensorflow_datasets as tfds +import tensorflow as tf +import numpy as np +from scipy.stats import truncnorm, norm +import galsim + +_DESCRIPTION = "This tfds generates random toy-model galaxy stamps." +_CITATION = "{NEEDED}" +_URL = "https://github.com/CosmoStat/autometacal" + +class GalGenConfig(tfds.core.BuilderConfig): + """BuilderConfig for GalGen.""" + + def __init__(self, *,dataset_size=100000, stamp_size=51, pixel_scale=.263, **kwargs): + """BuilderConfig for GalGen. + Args: + pixel_scale: pixel_scale of the image in arcsec/pixel. + stamp_size: image stamp size in pixels. + flux: flux of the profile. + **kwargs: keyword arguments forwarded to super. + """ + v2 = tfds.core.Version("2.0.0") + super(GalGenConfig, self).__init__( + description=("%d stamps of %d x %d pixels, %.2f arcsec/pixel, flux = 1." % + (dataset_size,stamp_size, stamp_size, pixel_scale)), + version=v2, + release_notes={ + "2.0.0": "New split API (https://tensorflow.org/datasets/splits)", + }, + **kwargs) + self.dataset_size = dataset_size + self.stamp_size = stamp_size + +class GalGen(tfds.core.GeneratorBasedBuilder): + """Random galaxy image generator.""" + + MANUAL_DOWNLOAD_INSTRUCTIONS = """\ + Nothing to download. DataSet is generated at first call. + """ + + BUILDER_CONFIGS = [ + GalGenConfig(name='small_100k', dataset_size=100000, stamp_size=51), + GalGenConfig(name='large_100k', dataset_size=100000, stamp_size=101), + GalGenConfig(name='small_1k', dataset_size=1000, stamp_size=51), + GalGenConfig(name='large_1k', dataset_size=1000, stamp_size=101) + ] + + VERSION = tfds.core.Version('0.5.0') + RELEASE_NOTES = {'0.5.0': "Updated functionalities, simpler."} + + def _info(self): + return tfds.core.DatasetInfo( + builder=self, + # Description and homepage used for documentation + description=_DESCRIPTION, + homepage=_URL, + features=tfds.features.FeaturesDict( + {'label': tfds.features.Tensor( + shape=[2], + dtype=tf.float32), + 'gal_model': tfds.features.Tensor( + shape=[self.builder_config.stamp_size, self.builder_config.stamp_size], + dtype=tf.float32), + 'obs_image': tfds.features.Tensor( + shape=[self.builder_config.stamp_size, self.builder_config.stamp_size], + dtype=tf.float32), + 'psf_image': tfds.features.Tensor( + shape=[self.builder_config.stamp_size, self.builder_config.stamp_size], + dtype=tf.float32), + }), + supervised_keys=("obs_image","label"), + citation=_CITATION) + + def _split_generators(self,dl): + """Returns generators according to split.""" + return {tfds.Split.TRAIN: self._generate_examples(self.builder_config.dataset_size, + self.builder_config.stamp_size)} + + def _generate_examples(self, dataset_size, stamp_size): + """Yields examples.""" + np.random.seed(31415) + + for i in range(dataset_size): + + #generate example + label, model, obs_img, psf_img, = gs_generate_images(stamp_size = stamp_size) + yield '%d'%i, {'gal_model': model, #noiseless PSFless galaxy model + 'obs_image': obs_img, #observed image + 'psf_image': psf_img, #psf image + 'label': label} + + +def gs_generate_images(**kwargs): + + """ Random Galaxy Generator + Generates noiseless galaxy images with a simple light profile. + The resulting image is before the convolution with a PSF. + Galaxy shapes follow a bivariate normal distribution centered in zero. + Args: + g_range: galaxy shapes go from -g_range and + g_range in each g1, g2 + g_scatter: galaxy shapes scatter + flux: galaxy flux (counts) + pixel_scale: intended pixel scale in arcsec^2/pixel + stamp_size: size in pixels of the NxN resulting image + method: galsim drawing method + mean_radius: mean half-light radii of generated galaxies + scatter_radius: scatter of half-light radii of generated galaxies + mean_snr: average snr of generated stamps (snr as per galsim definition) + scatter_snr: scatter in snr of generated stamps + Returns: + g1, g2: galaxy shape parameters + mode: galaxy model without psf or noise + gal: tensor with a 2-d array representing an observed image of a galaxy (with convolved psf) + psf: tensor with a 2-d array representing the model of the psf + """ + + defaults = {'g_range' : 0.7, #elipticity + 'g_scatter' : 0.25, # + 'mean_hlr': .5, #size + 'scatter_hlr': 0.1, # + 'psf_beta': 4.8, #psf + 'psf_fwhm': 0.9, # + 'mean_snr': 200, #snr + 'scatter_snr': 20, # + 'pixel_scale' : 0.263, # + 'stamp_size' : 51, # + 'method' : "no_pixel", # + } + + defaults.update(kwargs) + + #ellipticity range + a = -defaults['g_range'] / defaults['g_scatter'] + b = defaults['g_range'] / defaults['g_scatter'] + + g1 = truncnorm.rvs(a, b, loc=0, scale=defaults['g_scatter']) + g2 = truncnorm.rvs(a, b, loc=0, scale=defaults['g_scatter']) + + hlr = norm.rvs(defaults['mean_hlr'], defaults['scatter_hlr']) + + #very simple galaxy model + gal = galsim.Exponential(half_light_radius=hlr) + + #apply 'shear' + gal = gal.shear(g1=g1,g2=g2) + + + #create constant psf + psf = galsim.Moffat(beta=defaults['psf_beta'], + fwhm=defaults['psf_fwhm']) + + #draw galaxy before convolution + model_image = gal.drawImage(nx=defaults['stamp_size'], + ny=defaults['stamp_size'], + scale=defaults['pixel_scale'], + method=defaults['method']) + + #convolve galaxy and psf + gal = galsim.Convolve([gal,psf]) + + #draw psf image + psf_image = psf.drawImage(nx=defaults['stamp_size'], + ny=defaults['stamp_size'], + scale=defaults['pixel_scale'], + method=defaults['method']) + + + #draw final observed image + obs_image = gal.drawImage(nx=defaults['stamp_size'], + ny=defaults['stamp_size'], + scale=defaults['pixel_scale'], + method=defaults['method']) + + #add noise to image with a given SNR + noise = galsim.GaussianNoise() + snr = norm.rvs(defaults['mean_snr'],defaults['scatter_snr'],) + obs_image.addNoiseSNR(noise,snr=snr) + + #output everything to tf tensors + return [g1,g2], model_image.array.astype('float32'), obs_image.array.astype('float32'), psf_image.array.astype('float32') \ No newline at end of file From 86fc9316efa0b0686a2a7154c04ffe22af251f78 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Andr=C3=A9=20Z=2E=20Vitorelli?= Date: Thu, 11 Aug 2022 17:53:40 +0200 Subject: [PATCH 3/6] reestructuring datasets for practical purposes --- autometacal/python/datasets/CFIS.py | 150 -------------- autometacal/python/datasets/__init__.py | 7 +- autometacal/python/datasets/gal_gen.py | 131 ------------ autometacal/python/datasets/galaxies.py | 257 ------------------------ autometacal/python/datasets/simple.py | 3 +- notebooks/CFIS_dataset.ipynb | 199 +++++++----------- 6 files changed, 76 insertions(+), 671 deletions(-) delete mode 100644 autometacal/python/datasets/CFIS.py delete mode 100644 autometacal/python/datasets/gal_gen.py delete mode 100644 autometacal/python/datasets/galaxies.py diff --git a/autometacal/python/datasets/CFIS.py b/autometacal/python/datasets/CFIS.py deleted file mode 100644 index c6a4feb..0000000 --- a/autometacal/python/datasets/CFIS.py +++ /dev/null @@ -1,150 +0,0 @@ -""" TensorFlow Dataset of simulated CFIS images.""" -import tensorflow_datasets as tfds -import tensorflow as tf -import numpy as np -import galsim as gs - -_CITATION = """ -""" - -_DESCRIPTION = """ -""" - -class CFISConfig(tfds.core.BuilderConfig): - """BuilderConfig for CFIS Galaxies.""" - - def __init__(self, *, galaxy_type="parametric", - data_set_size=1000, stamp_size=51, - shear_g1=0.0, shear_g2=0.0, **kwargs): - """BuilderConfig for CFIS. - Args: - size: size of sample. - stamp_size: image stamp size in pixels. - pixel_scale: pixel scale of stamps in arcsec. - **kwargs: keyword arguments forwarded to super. - """ - v1 = tfds.core.Version("0.0.1") - super(CFISConfig, self).__init__( - description=("Galaxy stamps"), - version=v1, - **kwargs) - - # Adjustable parameters - self.galaxy_type = galaxy_type - self.data_set_size = data_set_size - self.stamp_size = stamp_size - self.kstamp_size = stamp_size - self.shear_g1 = shear_g1 - self.shear_g2 = shear_g2 - - # Fixed parameters - self.pixel_scale = 0.187 - self.psf_fwhm = 0.65 # TODO: add varying PSF - self.psf_e1 = 0.0 - self.psf_e2 = 0.025 - - -class CFIS(tfds.core.GeneratorBasedBuilder): - """DatasetBuilder for simulated CFIS dataset.""" - - VERSION = tfds.core.Version('0.0.1') - RELEASE_NOTES = { - '0.0.1': 'Initial release.', - } - - BUILDER_CONFIGS = [CFISConfig(name="parametric_1k", galaxy_type="parametric", data_set_size=1000), - CFISConfig(name="parametric_shear_1k", galaxy_type="parametric", data_set_size=1000, shear_g1=0.02)] - - def _info(self) -> tfds.core.DatasetInfo: - """Returns the dataset metadata.""" - # TODO: Specifies the tfds.core.DatasetInfo object - return tfds.core.DatasetInfo( - builder=self, - description=_DESCRIPTION, - features=tfds.features.FeaturesDict({ - 'obs': tfds.features.Tensor(shape=[self.builder_config.stamp_size, - self.builder_config.stamp_size], - dtype=tf.float32), - 'psf': tfds.features.Tensor(shape=[self.builder_config.stamp_size, - self.builder_config.stamp_size], - dtype=tf.float32), - # 'gal_kimage': tfds.features.Tensor(shape=[2, self.builder_config.kstamp_size, - # self.builder_config.kstamp_size], - # dtype=tf.float32), - # 'psf_kimage': tfds.features.Tensor(shape=[2, self.builder_config.kstamp_size, - # self.builder_config.kstamp_size], - # dtype=tf.float32), - "noise_std": tfds.features.Tensor(shape=[1], dtype=tf.float32), - "mag": tfds.features.Tensor(shape=[1], dtype=tf.float32), - }), - # If there's a common (input, target) tuple from the - # features, specify them here. They'll be used if - # `as_supervised=True` in `builder.as_dataset`. - supervised_keys=("obs", "obs"), - homepage='https://dataset-homepage/', - citation=_CITATION, - ) - - def _split_generators(self, dl_manager: tfds.download.DownloadManager): - """Returns SplitGenerators.""" - if self.builder_config.data_set_size: - size = self.builder_config.data_set_size - else: - size = 81499 - - return [ - tfds.core.SplitGenerator( - name=tfds.Split.TRAIN, - gen_kwargs={ - "size": size, - },), - ] - - def _generate_examples(self, size): - """Yields examples.""" - # Loads the galsim COSMOS catalog - cat = gs.COSMOSCatalog(sample="25.2") - mag_zp = 32 - sky_level = 400 # ADU (~variance) - psf = gs.Kolmogorov(fwhm=self.builder_config.psf_fwhm, flux=1.0) - psf = psf.shear(g1=self.builder_config.psf_e1, g2=self.builder_config.psf_e2) - - # Prepare borders for kimage - Nk = self.builder_config.kstamp_size - bounds = gs._BoundsI(-Nk//2, Nk//2-1, -Nk//2, Nk//2-1) - - for i in range(size): - # retrieving galaxy and magnitude - gal = cat.makeGalaxy(i, gal_type='parametric') - gal_mag = cat.param_cat['mag_auto'][cat.orig_index[i]] - gal_flux = 10**(-(gal_mag-mag_zp)/2.5) - - gal = gal.withFlux(gal_flux) - gal = gal.shear(g1=self.builder_config.shear_g1, g2=self.builder_config.shear_g2) - - gal_conv = gs.Convolve(gal, psf) - - gal_stamp = gal_conv.drawImage(nx=self.builder_config.stamp_size, - ny=self.builder_config.stamp_size, - scale=self.builder_config.pixel_scale - ).array.astype('float32') - - psf_stamp = psf.drawImage(nx=self.builder_config.stamp_size, - ny=self.builder_config.stamp_size, - scale=self.builder_config.pixel_scale - ).array.astype('float32') - - # gal_kimage = gal.drawKImage(bounds=bounds, - # scale=2.*np.pi/(self.builder_config.stamp_size*self.builder_config.pixel_scale), - # recenter=False).array.astype('complex64') - - # psf_kimage = psf.drawKImage(bounds=bounds, - # scale=2.*np.pi/(self.builder_config.stamp_size*self.builder_config.pixel_scale), - # recenter=False).array.astype('complex64') - - yield '%d'%i, {"obs": gal_stamp, - "psf": psf_stamp, - # "gal_kimage": np.stack([gal_kimage.real, gal_kimage.imag]), - # "psf_kimage": np.stack([psf_kimage.real, psf_kimage.imag]), - "noise_std": np.array([np.sqrt(sky_level)]).astype('float32'), - "mag": np.array([gal_mag]).astype('float32')} \ No newline at end of file diff --git a/autometacal/python/datasets/__init__.py b/autometacal/python/datasets/__init__.py index 12a30af..0a3318d 100644 --- a/autometacal/python/datasets/__init__.py +++ b/autometacal/python/datasets/__init__.py @@ -1,4 +1,5 @@ -"""gal_gen dataset.""" +"""galaxy datasets""" -from .gal_gen import GalGen -from .CFIS import CFIS \ No newline at end of file +from .galgen import GalGen +from .cfis import CFIS +from .simple import Simple \ No newline at end of file diff --git a/autometacal/python/datasets/gal_gen.py b/autometacal/python/datasets/gal_gen.py deleted file mode 100644 index 2995511..0000000 --- a/autometacal/python/datasets/gal_gen.py +++ /dev/null @@ -1,131 +0,0 @@ -"""gal_gen dataset.""" -import tensorflow_datasets as tfds -import tensorflow as tf -import numpy as np - -from .galaxies import gs_generate_images, gs_drawKimage - -_DESCRIPTION = "This tfds generates random toy-model galaxy stamps." -_CITATION = "{NEEDED}" -_URL = "https://github.com/CosmoStat/autometacal" - -class GalGenConfig(tfds.core.BuilderConfig): - """BuilderConfig for GalGen.""" - - def __init__(self, *,dataset_size=None, stamp_size=None, pixel_scale=None, flux=None, interp_factor=None, padding_factor=None, **kwargs): - """BuilderConfig for SQUAD. - Args: - pixel_scale: pixel_scale of the image in arcsec/pixel. - stamp_size: image stamp size in pixels. - flux: flux of the profile. - **kwargs: keyword arguments forwarded to super. - """ - v2 = tfds.core.Version("2.0.0") - super(GalGenConfig, self).__init__( - description=("GalGen %d stamps in %d x %d resolution, %.2f arcsec/pixel, with flux of %.2f ." % - (dataset_size,stamp_size, stamp_size, pixel_scale, flux)), - version=v2, - release_notes={ - "2.0.0": "New split API (https://tensorflow.org/datasets/splits)", - }, - **kwargs) - self.dataset_size = dataset_size - self.stamp_size = stamp_size - self.pixel_scale = pixel_scale - self.flux = flux - self.interp_factor = 2 - self.padding_factor = 1 - self.kstamp_size = self.interp_factor*self.padding_factor*self.stamp_size - - - -class GalGen(tfds.core.GeneratorBasedBuilder): - """Random galaxy image generator.""" - - MANUAL_DOWNLOAD_INSTRUCTIONS = """\ - Nothing to download. DataSet is generated at first call. - """ - - BUILDER_CONFIGS = [ - GalGenConfig(name='small_stamp_100k', dataset_size=100000, stamp_size=50, pixel_scale=.2, flux=1e5), - GalGenConfig(name='large_stamp_100k', dataset_size=100000, stamp_size=100, pixel_scale=.2, flux=1e5), - GalGenConfig(name='small_stamp_100', dataset_size=100, stamp_size=50, pixel_scale=.2, flux=1e5), - GalGenConfig(name='large_stamp_100', dataset_size=100, stamp_size=100, pixel_scale=.2, flux=1e5) - ] - - VERSION = tfds.core.Version('0.1.0') - RELEASE_NOTES = {'0.1.0': "Basic functionalities, that work."} - - def _info(self): - return tfds.core.DatasetInfo( - builder=self, - # Description and homepage used for documentation - description=_DESCRIPTION, - homepage=_URL, - features=tfds.features.FeaturesDict({'label': tfds.features.Tensor(shape=[2], dtype=tf.float32), - 'gal_model': tfds.features.Tensor(shape=[self.builder_config.stamp_size, - self.builder_config.stamp_size], - dtype=tf.float32), - 'obs_image': tfds.features.Tensor(shape=[self.builder_config.stamp_size, - self.builder_config.stamp_size], - dtype=tf.float32), - 'psf_image': tfds.features.Tensor(shape=[self.builder_config.stamp_size, - self.builder_config.stamp_size], - dtype=tf.float32), - 'obs_kimage': tfds.features.Tensor(shape=[2,self.builder_config.kstamp_size, - self.builder_config.kstamp_size], - dtype=tf.float32), - 'psf_kimage': tfds.features.Tensor(shape=[2,self.builder_config.kstamp_size, - self.builder_config.kstamp_size], - dtype=tf.float32), - 'psf_deconv': tfds.features.Tensor(shape=[2,self.builder_config.kstamp_size, - self.builder_config.kstamp_size], - dtype=tf.float32) - }), - supervised_keys=("obs_image","label"), - citation=_CITATION) - - def _split_generators(self,dl): - """Returns generators according to split.""" - return {tfds.Split.TRAIN: self._generate_examples(self.builder_config.dataset_size, - self.builder_config.stamp_size, - self.builder_config.pixel_scale, - self.builder_config.flux, - self.builder_config.interp_factor, - self.builder_config.padding_factor)} - - def _generate_examples(self, dataset_size, stamp_size, pixel_scale, flux, interp_factor, padding_factor): - """Yields examples.""" - np.random.seed(31415) - - for i in range(dataset_size): - - #generate example - label, model, obs_img, psf_img, obs_kimg, psf_kimg, psf_deconv = gs_generate_images(stamp_size = stamp_size, - pixel_scale = pixel_scale, - flux = flux, - interp_factor = interp_factor, - padding_factor = padding_factor) - - - #store complex arrays in 2,N,N - obs_kimg = decomplexify(obs_kimg.numpy()) - psf_kimg = decomplexify(psf_kimg.numpy()) - psf_deconv = decomplexify(psf_deconv.numpy()) - - - yield '%d'%i, {'gal_model': model.numpy(), #noiseless PSFless galaxy model - 'obs_image': obs_img.numpy(), #observed image - 'psf_image': psf_img.numpy(), #psf image - 'obs_kimage': obs_kimg, #obs k image - 'psf_kimage': psf_kimg, #psf k image - 'psf_deconv': psf_deconv, #psf deconv kernel - 'label': label.numpy() } - - -def decomplexify(arr): - arr=np.array([arr.real,arr.imag]) - return arr - -def recomplexify(arl): - return arl[0]+1j*arl[1] diff --git a/autometacal/python/datasets/galaxies.py b/autometacal/python/datasets/galaxies.py deleted file mode 100644 index 63c781f..0000000 --- a/autometacal/python/datasets/galaxies.py +++ /dev/null @@ -1,257 +0,0 @@ -import galsim -import tensorflow as tf -import numpy as np -from scipy.stats import truncnorm, norm - -""" -Simple wrappers for toy galaxy generation -Functions with gs_ prefix depend on galsim -""" - - -def gs_generate_images(**kwargs): - - """ Random Galaxy Generator - Generates noiseless galaxy images with a simple light profile. The resulting image is before the convolution with a PSF. - Galaxy shapes follow a bivariate normal distribution centered in zero. - - Args: - g_range: galaxy shapes go from -g_range and + g_range in each g1, g2 - g_scatter: galaxy shapes scatter - flux: galaxy flux (counts) - pixel_scale: intended pixel scale in arcsec^2/pixel - stamp_size: size in pixels of the NxN resulting image - method: galsim drawing method - mean_radius: mean half-light radii of generated galaxies - scatter_radius: scatter of half-light radii of generated galaxies - mean_snr: average snr of generated stamps (snr as per galsim definition) - scatter_snr: scatter in snr of generated stamps - interp_factor: interpolation factor for drawing k space images - padding_factor: padding factor for drawing k space images - Returns: - g1, g2: galaxy shape parameters - gal: tensor with a 2-d array representing an observed image of a galaxy (with convolved psf) - psf: tensor with a 2-d array representing the model of the psf - gal_k: k space image of gal - psf_k: k space image of psf - """ - - defaults = {'g_range' : 0.8, #elipticity - 'g_scatter' : 0.25, # - 'mean_radius': 1.0, #size - 'scatter_radius': 0.1, # - 'psf_beta': 5, #psf - 'psf_fwhm': 0.7, # - 'mean_snr': 200, #snr - 'scatter_snr': 20, # - 'flux' : 1e5, #flux - 'pixel_scale' : 0.2, # - 'stamp_size' : 50, # - 'method' : "no_pixel", # - 'interp_factor': 2, #kimage interpolation - 'padding_factor': 1 #kimage padding - } - - defaults.update(kwargs) - - #ellipticity range - a, b = (-defaults['g_range'] - 0) / defaults['g_scatter'], (defaults['g_range'] - 0) / defaults['g_scatter'] - g1=g2=1 - - #select ellipticities, ensure g<=1 - while g1**2+g2**2>1: - g1 = truncnorm.rvs(a, b, loc=0, scale=defaults['g_scatter']) - g2 = truncnorm.rvs(a, b, loc=0, scale=defaults['g_scatter']) - - re = norm.rvs(defaults['mean_radius'], defaults['scatter_radius']) - - #very simple galaxy model - gal = galsim.Exponential(flux=defaults['flux'] , - half_light_radius=re) - - #apply shear - gal = gal.shear(g1=g1,g2=g2) - - - #create constant psf - psf = galsim.Moffat(beta=defaults['psf_beta'], - fwhm=defaults['psf_fwhm']) - - #draw galaxy before convolution - model_Image = gal.drawImage(nx=defaults['stamp_size'], - ny=defaults['stamp_size'], - scale=defaults['pixel_scale'], - method=defaults['method']) - - #convolve galaxy and psf - gal = galsim.Convolve([gal,psf]) - - #draw psf image - psf_Image = psf.drawImage(nx=defaults['stamp_size'], - ny=defaults['stamp_size'], - scale=defaults['pixel_scale'], - method=defaults['method']) - - - #draw final observed image - obs_Image = gal.drawImage(nx=defaults['stamp_size'], - ny=defaults['stamp_size'], - scale=defaults['pixel_scale'], - method=defaults['method']) - - #add noise to image with a given SNR - noise = galsim.GaussianNoise() - snr = norm.rvs(defaults['mean_snr'],defaults['scatter_snr'],) - obs_Image.addNoiseSNR(noise,snr=snr) - - #draw kimage of galaxy - obs_k = gs_drawKimage(obs_Image.array, - defaults['pixel_scale'], - interp_factor=defaults['interp_factor'], - padding_factor=defaults['padding_factor']) - - # draw kimage of the psf - psf_k = gs_drawKimage(psf_Image.array, - defaults['pixel_scale'], - interp_factor=defaults['interp_factor'], - padding_factor=defaults['padding_factor']) - - #draw psf deconvolution kernel - psf_deconv = gs_Deconvolve(psf_Image.array, - defaults['pixel_scale'], - interp_factor=defaults['interp_factor'], - padding_factor=defaults['padding_factor']) - - #output everything to tf tensors - # tfds names: - g = tf.convert_to_tensor([g1,g2],dtype=tf.float32) # label - model = tf.convert_to_tensor(model_Image.array,dtype=tf.float32) # model_image - obs = tf.convert_to_tensor(obs_Image.array,dtype=tf.float32) # obs_image - psf = tf.convert_to_tensor(psf_Image.array,dtype=tf.float32) # psf_image - - return g, model, obs, psf, obs_k, psf_k, psf_deconv - - -def gs_drawKimage(image, - pixel_scale=0.2, - interp_factor=2, - padding_factor=1): - """ - Args: - image: numpy array - pixel_scale: telescope pixel scale in arcsec/pixel - interp_factor: interpolation factor for superresolution - padding_factor: padding added by fraction - - returns: - tensor with k-image of object - """ - - #prepare borders - N = len(image) - Nk = N*interp_factor*padding_factor - bounds = galsim._BoundsI(-Nk//2, Nk//2-1, -Nk//2, Nk//2-1) - - #interpolated galsim object from input image - img_galsim = galsim.InterpolatedImage(galsim.Image(image,scale=pixel_scale)) - - #draw galsim output image - result = img_galsim.drawKImage(bounds=bounds, - scale=2.*np.pi/(N*padding_factor*pixel_scale), - recenter=False) - - return tf.convert_to_tensor(result.array,dtype=tf.complex64) - - -def gs_Deconvolve(psf_img, - pixel_scale=0.2, - interp_factor=2, - padding_factor=1): - """ - Returns a deconvolution kernel of a psf image. - - Args: - psf_img: numpy array representing the psf model - pixel_scale: the pixel scale of the image, in arcsec/pixel - interp_factor: the interpolation factor for super-resolution - padding_factor: a factor to add side pads to the image - Returns: - A complex tensorflow tensor that is a deconvolution kernel. - - """ - - N = len(psf_img) - - psf_galsim=galsim.InterpolatedImage(galsim.Image(psf_img,scale=pixel_scale)) - ipsf=galsim.Deconvolve(psf_galsim) - Nk = N*interp_factor*padding_factor - bounds = galsim._BoundsI(-Nk//2, - Nk//2-1, - -Nk//2, - Nk//2-1) - imipsf = ipsf.drawKImage(bounds=bounds, - scale=2.*np.pi/(N*padding_factor* pixel_scale), - recenter=False) - return tf.convert_to_tensor(imipsf.array,dtype=tf.complex64) - -def gs_noise_generator(stamp_size=50,variance=5,pixel_scale=.2,interp_factor=2,padding_factor=1): - """ Generate a noise k space image using GalSim. - - Args: - stamp_size: in pixels - variance: noise variance - pixel_scale: in arcsec/pixel - interp_factor: interpolation factor for k space - padding_factor: padding wrap factor - Returns: - A complex64 numpy array. - - """ - noise = galsim.GaussianNoise().withVariance(variance) - noise_image = galsim.Image(stamp_size,stamp_size, scale=pixel_scale) - noise.applyTo(noise_image) - noise_image = galsim.InterpolatedImage(noise_image) - Nk = stamp_size*padding_factor*interp_factor - from galsim.bounds import _BoundsI - - bounds = _BoundsI(-Nk//2, Nk//2-1, -Nk//2, Nk//2-1) - imnos = noise_image.drawKImage(bounds=bounds, - scale=2.*np.pi/(stamp_size*padding_factor*pixel_scale), - recenter=False) - return imnos.array.astype('complex64') - -def make_data(Ngals=1, - snr = 200, - scale = 0.2, - stamp_size = 51, - psf_fwhm = 0.9, - gal_hlr = 0.7, - gal_g1 = [0], - gal_g2 = [0], - flux=1.e5): - """Simple exponetial profile toy model galaxy""" - - gal_list = [] - psf_list = [] - - for n in range(Ngals): - psf = galsim.Moffat(beta=2.5, - fwhm=psf_fwhm) - - obj0 = galsim.Exponential(half_light_radius=gal_hlr,flux=flux).shear(g1=gal_g1[n],g2=gal_g2[n]) - obj = galsim.Convolve(psf, obj0) - - psf_image = psf.drawImage(nx=stamp_size, ny=stamp_size, scale=scale).array - gal_image = obj.drawImage(nx=stamp_size, ny=stamp_size, scale=scale) - noise = galsim.GaussianNoise() - gal_image.addNoiseSNR(noise,snr=snr) - - gal_image = tf.convert_to_tensor(gal_image.array) - psf_image = tf.convert_to_tensor(psf_image) - gal_list.append(gal_image) - psf_list.append(psf_image) - - gal_image_stack = tf.stack(gal_list) - psf_image_stack = tf.stack(psf_list) - - return gal_image_stack, psf_image_stack \ No newline at end of file diff --git a/autometacal/python/datasets/simple.py b/autometacal/python/datasets/simple.py index 25a312e..65f55eb 100644 --- a/autometacal/python/datasets/simple.py +++ b/autometacal/python/datasets/simple.py @@ -12,8 +12,7 @@ class SimpleConfig(tfds.core.BuilderConfig): """BuilderConfig for GalGen.""" - def __init__(self, *,dataset_size, stamp_size, - ,**kwargs): + def __init__(self, *,dataset_size, stamp_size,**kwargs): """BuilderConfig for Simple. Args: dataset_size: number of stamps diff --git a/notebooks/CFIS_dataset.ipynb b/notebooks/CFIS_dataset.ipynb index 8965016..ad14be5 100644 --- a/notebooks/CFIS_dataset.ipynb +++ b/notebooks/CFIS_dataset.ipynb @@ -3,6 +3,7 @@ { "cell_type": "code", "execution_count": 1, + "id": "54f52a71", "metadata": {}, "outputs": [ { @@ -11,6 +12,13 @@ "text": [ "Populating the interactive namespace from numpy and matplotlib\n" ] + }, + { + "name": "stderr", + "output_type": "stream", + "text": [ + ":219: RuntimeWarning: scipy._lib.messagestream.MessageStream size changed, may indicate binary incompatibility. Expected 56 from C header, got 64 from PyObject\n" + ] } ], "source": [ @@ -22,6 +30,7 @@ }, { "cell_type": "markdown", + "id": "6b8582dd", "metadata": {}, "source": [ "## Loading the data\n", @@ -32,80 +41,36 @@ }, { "cell_type": "code", - "execution_count": 2, + "execution_count": 3, + "id": "a085afdf", "metadata": {}, - "outputs": [ - { - "name": "stdout", - "output_type": "stream", - "text": [ - "\u001b[1mDownloading and preparing dataset Unknown size (download: Unknown size, generated: Unknown size, total: Unknown size) to /local/home/flanusse/tensorflow_datasets/cfis/parametric_shear_1k/0.0.1...\u001b[0m\n" - ] - }, - { - "data": { - "application/vnd.jupyter.widget-view+json": { - "model_id": "", - "version_major": 2, - "version_minor": 0 - }, - "text/plain": [ - "HBox(children=(FloatProgress(value=0.0, description='Generating splits...', max=1.0, style=ProgressStyle(descr…" - ] - }, - "metadata": {}, - "output_type": "display_data" - }, - { - "data": { - "application/vnd.jupyter.widget-view+json": { - "model_id": "", - "version_major": 2, - "version_minor": 0 - }, - "text/plain": [ - "HBox(children=(FloatProgress(value=1.0, bar_style='info', description='Generating train examples...', layout=L…" - ] - }, - "metadata": {}, - "output_type": "display_data" - }, - { - "data": { - "application/vnd.jupyter.widget-view+json": { - "model_id": "", - "version_major": 2, - "version_minor": 0 - }, - "text/plain": [ - "HBox(children=(FloatProgress(value=0.0, description='Shuffling cfis-train.tfrecord...', max=1000.0, style=Prog…" - ] - }, - "metadata": {}, - "output_type": "display_data" - }, - { - "name": "stdout", - "output_type": "stream", - "text": [ - "\u001b[1mDataset cfis downloaded and prepared to /local/home/flanusse/tensorflow_datasets/cfis/parametric_shear_1k/0.0.1. Subsequent calls will reuse this data.\u001b[0m\n" - ] - } - ], + "outputs": [], + "source": [ + "\n", + "\n", + "dset = tfds.load('CFIS/parametric_shear_1k', split='train')\n" + ] + }, + { + "cell_type": "code", + "execution_count": 4, + "id": "c5f6270c-07f1-4fc3-8d92-8b50b191583e", + "metadata": {}, + "outputs": [], "source": [ "# We create a function to add noise on the fly for data augmentation\n", "@tf.function\n", "def add_noise(example):\n", - " im_noise = example['obs'] + example['noise_std'] * tf.random.normal([51,51])\n", + " im_noise = example['gal_image'] + example['noise_std'] * tf.random.normal([51,51])\n", " return im_noise, example\n", - "\n", - "dset = tfds.load('CFIS/parametric_shear_1k', split='train')\n", + " \n", "dset = dset.map(add_noise)" ] }, { "cell_type": "code", - "execution_count": 3, + "execution_count": 5, + "id": "e21a9069", "metadata": {}, "outputs": [], "source": [ @@ -130,14 +95,15 @@ }, { "cell_type": "code", - "execution_count": 4, + "execution_count": 6, + "id": "5f0b507a", "metadata": { "scrolled": true }, "outputs": [ { "data": { - "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYsAAAEHCAYAAABfkmooAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/d3fzzAAAACXBIWXMAAAsTAAALEwEAmpwYAAAm4klEQVR4nO3df5RUZ5kn8O/TRYHVRKlG0CEVCCTkwASJ9KQ1KM4eQUcSI6QncYIxmR1HN5ido7swEQeMa8CJE2ZYTcYdZ5Rdc8xuMDYa7KWXuLArOHpQsgG7kbQBEzUhKaMhQqFJF6GofvaPqlt96/b9VdV1695b9f2ckxO6uurWyz1NPf2+z/s8r6gqiIiI3HSEPQAiIoo+BgsiIvLEYEFERJ4YLIiIyBODBREReWKwICIiT5PCHsBEzJgxQ+fOnRv2MIgo4n5x6hUAwGUzp4Y8kmg4cuTIS6o6s5bXxDpYzJ07F4cPHw57GEQUcWu+8iMAQN9H3xbySKJBRJ6t9TWRWYYSkT8UkS+LyLdE5N+HPR4iIhoTaLAQkQdE5EURecLy+LUickJEnhaRjQCgqk+q6h0AbgawLMhxERFRbYKeWXwNwLXmB0QkAeBLAK4DcCWAW0TkyvL3VgPYA+DRgMdFREQ1CDRYqOr3AZy2PPxWAE+r6i9U9TyAbwC4ofz83ap6HYBbgxwXERHVJowEdwbAc6avnwdwjYi8E8CNAKbAZWYhImsBrAWAOXPmBDZIIiIaE5ndUKr6PQDf8/G87QC2A0BPTw9b5hIRNUEYu6GyAGabvr6k/BgREUVUGMHicQBXiMg8EZkM4AMAdtdyARFZJSLbz549G8gAiYioWtBbZx8G8CMAC0TkeRH5iKpeAPAxAHsBPAlgp6oO13JdVR1Q1bXTpk1r/KCJiGicQHMWqnqLw+OPgttjiYhiIzIV3EREFF2xDBbMWRARNVcsgwVzFkREzRXLYEFERM3FYEFERJ5iGSyYsyAiaq5YBgvmLIiImiuWwYKIiJqLwYKIiDzFMlgwZ0FE1FyxDBbMWRARNVcsgwURETUXgwUREXlisCAiIk8MFkRE5CmWwYK7oYiImiuWwYK7oYiImiuWwYKIiJqLwYKIiDwxWBARkScGCyIi8sRgQUREnmIZLLh1loiouWIZLLh1loiouWIZLIiIqLliHSxyI4Wwh0BE1BZiHSx+/btzYQ+BiKgtxDpYFIqjYQ+BiKgtxDpYAMCyrfvRP5gNexhERC0t9sEim8tjXd8QlmzZx6BBRBSQ2AcLQy5fwKZdxxgwiIgCEMtgYRTljZ57perxfKGIbXtPhDQqIqLWFctgYRTldbxm6rjv/SqXD2FEREStLZbBws3F6VTYQyAiajktFSxSyQQ2rFwQ9jCIiFpOrIPF7K5OZNIpCIBMOoV7b1yM3u5M2MMiImo5k8IewESkO5M4uHFF2MMgImp5sZ5ZEBFRc8R6ZmHWP5jFtr0n8KtcHhenU9iwcgGXpIiIGqQlgkX/YBabdh1DvlAEUKrq3rTrGAAwYBARNUCsg8XxX/8e8zbuQYcIiqpV3zMK9BgsiIgmLtbBolAchQLjAoUhyAI9LnsRUTuJdbDwElSBHpe9iKjdxHI3lFNvKLMgC/S27T1RCRQG9qUiolYWy5mFqg4AGJgy64rbzY8nRDCq6rksNNElJKflLfalIqJWFctgYSeVTPiq4G7EEtLF6RSyNoGBfamIqFXFchnKkEx01NzqoxFLSBtWLkAqmah6zLzs1T+YxbKt+zFv4x6e5EdELSHWM4uFf/BaHN56fU2vacQSkhGUtu09gWwuj4RIJeAcfvY0HjmS9Zy5cDcVEcVJrGcW9XBaKqp1Cam3O1OZYRhbd7O5PHYcOuk5czGWwrK5PBRjAYUzECKKqrYLFl5LSLWwW9Kyr/ionrlwNxURxU2sl6HqYV5CmugSUC1LV+aZC3dTEVHctFSw8JsH6O3ONCQ/4LQryso6c3HbTWX+O6Q7k1AFzuYLzGsQUahaZhnKLg+wrm8IS7bsCywXYLekZcfYqWXsksrm8hDLc1LJBJYvnFn1dzgzUkAuX/DMa3D3FREFrWWChV0eAABy+UJgyePe7gzuvXExMi7J8YSUwoI5mAGl3IYRMIytvweOn7L9Oxjs8hpMlhNRM7RMsHBb7w8yedzbncHBjStw/5oltrOMoio27TqGzbuHbZPhmXQKBzeuQG93xlfOwvocJsuJqBlaJlh4bX0NOnlszDKMmYRZvlBELl/wHFe6M+n5Psbf07yk5XVdIqKJaplg4ZU/aEYrjt7uDEYd2qU7MX/4v3zugutzBaVlpu7P7sOGbx51Ta6z9QgRNVLL7IYydgltGRjGmZHq3+Kd6igaXUXdP5i1PYgJADoEGLU8bB7Xtr0nULA+wcL4rvXvZxVkx10iak8tEyyAsS2xfoJAo8+kMK7ndBCTNQ6kU0lsXr2o8l6NWjbKuPx9ne4JW48QkZeWChYGP3UUbonhej4onXZjOZk6ZVLV+/it2XBjBIpte09gfd9Q5YMfwLjAuK5vCFsGhnH9VbN89bIiovYWqWAhIr0ArgfwOgBfVdV9Qb1Xo6uoa32d9fkbVi6o+kCvlblOw/rBP2VSh+11z4wUsOPQyXEtSvKFIu7ceRQAAwYRlQSe4BaRB0TkRRF5wvL4tSJyQkSeFpGNAKCq/ap6O4A7AKwJclyNaijo9Tq73VF2z7fWbNi/akwyIUinklUt2u3qNNx2YgHOvayMLb+s1yAioDkzi68B+CcA/914QEQSAL4E4E8APA/gcRHZrao/LT/l0+XvB8buN/mJJIadrnfT1ZmqZR639zEvn1nzCMsXzsSB46dc8wrr+4bqGrsTc70GcxpE7S3wYKGq3xeRuZaH3wrgaVX9BQCIyDcA3CAiTwLYCuA7qvrjIMfl1FAQAJZt3V/zB6Nbg8KeS6fX/GFbT/8qp7xHV2cS5wqjjktcAucZhrGUxZwGUXsLK2eRAfCc6evnAVwD4OMA3g1gmojMV9UvW18oImsBrAWAOXPmTGgQ1g/kie6QcvqAb1TjQi9Os5u7Vy0CAGzePTxuScqY/Tz82HO2O7mMg53MmNMgaj+RSnCr6hcBfNHjOdsBbAeAnp6e2irgPDR6h1SzebVfd9tW3HPpdGz45tGqWo9khzjWfhg5DfP7ElHrCitYZAHMNn19SfmxULXCORNesxjX71uz6lJawnIqAoxTICWiiQkrWDwO4AoRmYdSkPgAgA/6fbGIrAKwav78+Q0dlNs5E61u294TKBSrZxGFokK1tFTllO+oJZCaZzXTUkmIlLbvJspV704FhUQUvmZsnX0YwI8ALBCR50XkI6p6AcDHAOwF8CSAnao67PeaqjqgqmunTZvW0LE28sjVuHH60D+bLzg2SAT8B1JrK/VcvlCZsZjPMOd2XaJoasZuqFscHn8UwKNBv38tGnnkaty4zaqMv7/frcZ2eRG/Fe5c2iKKpkgluKOgWTuXosar7sQukC5fONN3a5FaxClHRNQuYhksgspZtDM/sypr0WAtrUVq0Q45IqK4iWWwUNUBAAM9PT23hz2WVmI3q3Laauu0zXiigSKZkLbIERHFTSyDBTWH3exhwzeP2p4ZUqt0KolCcRSvnLcEl4ZWzhBRo7TMSXnUeHazh8KougaKrs6k64mFQKmcY+ju9yDdOXnc9wqjik/t+kld4yWi4MQyWIjIKhHZfvbs2bCH0tJqTTQbrUXuvXExulzOEzdyEk7XHymM4tP9x2p6byIKViyDRVB1FlStlkRzV2cS9964uJL3GPzMe3Db0jnjisLNO6zcrv/wY2Otw/oHs1i2dT/mbdyDZVv3sw6DKASxDBbUHHZFik7OFUbHPXZP72Lct2YJMukUBKU8xWuSHVjfN4RlW/dj+cKZjtczCvWsxXx+C/cYYIgai8GCHBkHMrktKRnMZ19Yr3Fw4wrct2YJXr0wijMjhcqH/iNHso6HPBkV427NHZ3UG2CIyBmDBbnq7c6gc7K/TXNuOQ6nD/3XJO1/BG+5ZrbrNbO5vOOHfz0BhojcxXLrLIvymstvotuagzDXaDjtiD1XGMVtS+dUztMQAVKTOvDQoZN4+LHnXHfSmluk+3kvVoYT1S+WMwsmuJvLT6Lb2ifKuhTkZFoqiXt6F+Pn974X969Zgg4RjJTzH3aHMZkZswW/78XKcKL6xTJYUHPZJbqTCUE6lYQAyKRTlZ1QBr+NA83NbO/69jEUHQ5bcpLN5X29V7t0DyYKSiyXoai56unG63fJJ2cq8BtXze2DALbdcs3ft56vns3leYYGUY0YLMiXWrvxOrU8t3veRChQ+eC3yqRTOLhxBYDSspj52FjrGRoAqpoktmObeiI3sVyGYgV39Pmp0bAuDTltozWe66So6nlo1ebdw47niZt3SnHbLZG9WAYLJrijz6jRMAryMukUbls6p+pra57j1qVzbK+VSnbgpqszjqf1Gddyu3Yu79740Fg247ZbIntchqLA1Lp0dU/vYgCobKNNiOCWa2aj59Lp2LTrmO1SkzGDmOihVV79qrjtltpd3cFCROao6slGDobont7FlaBhWLZ1v+1up4TIuBmEk67OpGO3XGu/KqfjZe2Y8xvpziRUS+eWM9dBrcZzGUpE3iYi7xeRN5S/vkpEvg7gYOCjI4Lzb/WjqlUfxm79oO5etQiJjvHLWALgpqvHZiV2uRa3s8bN+Y0zIwXk8gXmOqgluc4sRGQbgPcBGALwNyKyF8C/A3AvgA8HPjoi+Ptt37rbyTioCQAOP3saOx47CbsaPwWw41BpgnxP7+Katgl71XfkC0XcubM0Bs4wKO68lqGuB9CtqudEpAvAcwDepKrPBD4yorINKxdUndgH+NvtVBhV/HXfEMb3w61mBIyeS6dXch9+Ptz95DGKquO25hLFkdcy1DlVPQcAqnoGwFNRCBTcOtte7HZW+d3t5BUoDArY7nhyW9qalvLuxgvY76ZiC3WKG6+ZxWUistv09Tzz16q6OphhuVPVAQADPT09t4fx/tR8E93t5Id1pmB3Brl5luCwk9dWNpfHki37sHn1IgBwvS5RFHkFixssX38+qIEQTYTbbie/OkTQP5ityls41Vz0dmeqWpX4kcsXsK5vCFMnJ1yvSxRFrsFCVf+1WQMhmojrr5qFhw5NbCe3Nb/g1K7EeNxvSxMrpx5YrOWgKPPaDXUAcOz6rKr6rsYPiah2B46fash1zL/hO/WcAoC5G/c05P3MjN1d7E1FUeS1DPUJm8eWAvgkgBcbPxyi+jTyt3LjFD6v8zQarXNyB+Zt2lO1xdfYArxlYBi5ERb7UXhcd0Op6hHjPwAXAfh7ALcAuENV39KMARL54VRhnU4lK7uoujqT6HQ4xtVq065jSPvc7dQoT734im0tSGFUq84uZ7EfhcFPBfdKEfkBgP8E4HOq+g5V/U7wQyPyz6nyevPqRTi4cQXuW7ME5wqjlVP4vOQLxdIRrx6dc8NgFPtx2y01k2uwEJHHAXwFwDdQWno6KyJ/ZPzXjAES+eFVi+H35D6z3Eihcs2oKapWZhrr+oawZMs+Bg0KlFfO4hUALwN4f/k/6yR5RRCD8iIiqwCsmj9/fhhvTxHlVotRT07j4nSqck3jhL2oyuULrNWgQHktQ30SwAdVdbmqLgfwIErB4wmUgkcoeJ4F1coppzF1csL20CVrOxE/hzmFzVopzipxaiSvYPFlAK8CgIj8G5QaCD4I4CyA7cEOjahx7D7skwnB+Quj46bLXZ3Jce1ErMtcTgcxAaVAc5vDQU5BM3Zy8cQ/ajSvZaiEqp4u/3kNgO2q+giAR0RkKNCRETWQXTfZV169YNtTqnPyJNulHPMyl7UViCGdSmLz6kXo7c5MuEiwXpt2HcOUSR22VeLr+oawbe8Jbr+lmnkGCxGZpKoXALwLwNoaXksUKdacxjyHwjo/+Y1aWpk3W75QdE3msxcV1cPrA/9hAP8qIi8ByAP4AQCIyHyUlqKIYqvWU/GsnBLqRgV2lOULRWzePRzJYEfR5NUb6nMi8l0AswDsU62UDHUA+HjQgyMKkp9zMvwyAkQ2l4fAuUdOlOTyhcoyHGcb5MVzKUlVD9k89rNghkPUPI1aSrLmL+IQKOwYxX7r+oYqfbEynHFQGfMO1NYacU5GPQV/VlMnJ6CqvivMg2L0wzL+n83l8dc7hwBwxtHu/DXKISJHtRT8JURw29I5VZXmty2dg/NF/61Imm1UgfV9Q6zXaHOcWRBNkN9zLVLJxLj6DQDo/uw+FIrRXrwyRme0F9m8e7iyRRhgW/V2wGBBNEF2iXIjye1n7d/thL+pkxOVw5IEwNsvn47/98yZ0INLLl/A+r4h3PXtY+MOc2KyvDXFMliwNxRFSZA1F8OfvXbcY+bf4jtcDmgKmsL51D8eE9t6YhksVHUAwEBPT8/tYY+FCJhYojydStpWkvs5TyOsQOEHj4ltLUxwE4Vs8+pFSHZU95pKdgg2r1407rnWnk9OOpMdobdWT3c29/AoClYsZxZErcTvMlb/YBZ37jzqazYxUhjFSMi/2b987gL6B7Po7c6gfzCLzbuHKzOors4k7l61iMtUMcJgQRQB1oBhtAuxNi6M8rKTVWFUK3+PDd88isLo2NjPjBSw4VtHATAJHhcMFkQRYK0Ct+4oakThXxiMrbZ2CkXFur4h3LnzKJZe1oVnfpvn1tsIY86CKALsgoHRUtzrlD7nkzXioaiKgz8/XXX2xvq+IXy6/1jYQyMTBguiCHDbOWQ0J7STEMGtS+dE/hS/WimAhw6dZMCIEAYLogjwaouuGD+DSCUT+PzNb8Y9vYsrp/jB5nlx9tChk1iyZR9bjEQAcxZEEWBXBW6lKPWSslvXt57iZ91ZBcDz+oZEiIV+dnJ552Q424w0D4MFUQSYd0M55Scy6RQOblzh61pOH5hu1zdEKVAYjGT4tr0nsHzhTOz5yQvj2qRYNwWYA0ky0YHZ08fP3hhs/BON4A+GXz09PXr48OGwh0HUUHbnezs1IayHV8LcIALE9eOhqzOJl89dqNqu2yHAF25e4nqOeiPvc5SJyBFV7anlNcxZEEVMb3emkoMw2pg38gNsw8oF/hLiCty/ZgkSHfHLgpwZKVQFCqDUat183K3TDrSoH4kbFi5DEUVQIw5lcrs2AM9mhBenU5XnbhkYriz7pFNJvHqhiHxEz99wk83lK1XlTjvQ2NPKHoMFURuwW5s38h9OyzHms8jPmQJDLl9AMiHoABC/cIFKXsPpHBKvnWntisGCqAWZg0Pasn5vTQR79aayW64pFBVd5UaBbudxRFG+UMSWgWHcvWqRZ5CkMQwWRC3GOlOw+zC3njfhtuzltCxzZqQAQSmZ/GqhGNljYe0Y9+SmqzN4+LHnUFRFQgQ3XR3c8l/cRSZYiMhlAO4CME1V3x/2eIjiym8fKb9r827HxipKH7ypZAL3r1mCw8+exkOHTtYy3NCs7xuqavNeVMWOQyfx0KGT6OpMQhU4my9wS21ZoLuhROQBEXlRRJ6wPH6tiJwQkadFZCMAqOovVPUjQY6HqB3UEgT88LN7yljaOXD8lK9rRoHdrmDjsTMjBeTyhUqvqk27jk2oirx/MItlW/dj3sY9WLZ1fywr0oPeOvs1AFXnQopIAsCXAFwH4EoAt4jIlQGPg6ht+AkCtazNW7fyOjkzUvBVvxFH+UIRm3cP1/Va64FVjQg+YQi8KE9E5gL4X6r6pvLXbwOwWVVXlr/eBACqem/562/5XYaafukf6p986oFAxk0UVy+9/Cp++dIrMJcZCICODkFxVDG5XM0846IpdV1/8GQO54vxyU80UqJDMPf1nTXdO6f7NTnRge456QaOzr+dd7w9FkV5GQDPmb5+HkBGRF4vIl8G0G0EEDsislZEDovI4UIhXrswiJphxkVTMG/GVExOlP55T0504LKZU9FzaReumTcd3XPSdQcKAJg9PYVa6vQunzm1pudHWXFU8cuXXsFLL7/q+zVOgTVuATcyCW5V/S2AO3w8bzuA7UCp3UffR98W9NCIqMzYkjvqc0Eik07hu3e+s/K6VlimGtVS3Ynfzx6n9iqZdMr3NRptp+cn7XhhBIssgNmmry8pP0ZEEWZXvOfGnBcxtubWeo2oMirBzbUpyxfOxIHjp5DN5SudezPlxx85ko19PUcYOYtJAH4G4F0oBYnHAXxQVX1nj0RkFYBV8+fPv/2pp55q/KCJaBzfDQgB1+2m1g/Zua9P4Yc/P227O6kVJBOCNW+ZjQPHT0Wmu209jQQDDRYi8jCAdwKYAeA3AO5W1a+KyHsB3A8gAeABVf1cPddn11mi2kykJfe8jXs8P9D9tlE3j6cVZhpeujqTGPzMe8IeRkU9wSLQZShVvcXh8UcBPBrkexNRNesHs7Xthxe34jzA/9KKOWA5NTEE4LiEE0fmKnpz/sbcBj6dSmLz6kWRLf6LTIKbiILl1pLbzweU3Wl+grET/PzMUqwByylQCFA1Q9lx6GRNy1SdyY7ItR/5dP8xfP2xk1WbA8x//Vy+gHV9Qzj87Gnc07u4+QP0EMtgYcpZhD0UotiYaEtur4aDfvhtRWIuLDxw/FTN+Yy/u/GqSLUeSQh8j+WhQyfRc+l01yNk0zb9uLo6k7h7VXAzk1gGC1UdADDQ09Nze9hjIYqLRrTktjvre33fkO/A4ScwWZez6jlforc7E5lDjDoAFGuMdlsGhqvu8+bdw8jlx5ay7JpDnhlxPqu8EXhSHlGbsOvxVO8WznpbWDgFpoSI46mA9Z4vEZVDjOpZDDOCgXGfzYHCTaGogQVJnsFN1EYmshvKzK3QzG03VD3nXtezYyqTTmHk/IXYnbVhNqlDcMFv9aPFM1uvd/1+onPaL4sjZy+raTx1jSRkzFkQ1adRx7XWm//wynu4BTPrUoybbC6PZIcgmRAUal0Dioh6AwUALPrM/8bI+aLtLwT9g1lMet0b5tV6Tc4siKhm9c4s3LjNOgBUtpsmXLbbWqVTSUydMqkl2ozUy9ixZr5vLzy4Dq++8FRNHbuYsyCimjUy/2Fw2tq7ZWC4kh8BnLfb2jmbL+DgxhWVportyLhbtdw3O7FchiKicDViG62V2/Gt9TKS414dXo3fvskZgwUR1aVR+Q+DV4V4rcwzncmJDteAwUDhLZZzMxFZJSLbz549G/ZQiKhBnJa20qmk7fPTqSQSYr/snhCp2mE1e3p9229pTCyDhaoOqOraadOmhT0UImoQ6/GtRs3F5tWLbIPI5tWLMOqwDj+qWjXrmXHRFExqlROYQsJlKCKKDLelLbv8iNOBSnaFfJe+vhMnfzuCwgS2pLYzBgsiijynIGLX3NBpV9aMi6bg4yuuwLq+oSCH2rJiuQxFRAQ4L105zU56uzPI1Nk+pN3FcmbBCm4iMnjtyuofzGLwZA7ni6NYtnU/li+cib7Hn6upsjuVTOCmqzPjWoy3k1jOLJjgJiI/jKpwY9tsNpfHI0eyWPOW2ejqtN9lBZSOQk2nklWzlXt6F+MLNy9BMtGeifJYziyIiMyceko5VYUfOH6q6phTvw0WjcfaMe/BYEFEseZ2XKzfhoe1FBj2dmdw586jE26fETexXIYiIjK4HRfrdBZGvWdkGNotUAAMFkQUc26zhyAaHgKoaUdVq+Q4GCyIKNbcZg/G1lqj66zX1lq/NqxcAKcQkE4lq5LnUye3xmp/LIMFe0MRkcFr9tDbnUH3nDSumTcdBzeuaEjzw97uDG5dOmdcwEglE3jfm2fhXGGsaaHfA5uiLpbBgltnichQa2Feo9zTuxj3rVky7n0PHD9V0xGwcdEa8yMiamuNbpc+kfddX+e22q7OJDonR/dUv1jOLIiIosoph9LVmXTMcwiAu1ct8jzDPEwMFkREDWSXQxEA1181y/GQJePxiW7pDRKDBRFRA/V2Z3DT1ZmqWYQCeORI1vEgJwDYtOsYli+cOS7QRAWDBRFRgx04fmrcLCJfKEIEjsHAaENiJOujhsGCiKjBnHIPuZECbrraORH/q1wevd0ZHNy4InIBg8GCiKjB3AoFDxw/5et1UUt2xzJYsCiPiMLQP5jFsq37MW/jHizbuh/9g1nb57kVCroFAXMbkqglu2MZLFiUR0TNZnS3zebyUIx1t7ULGG6Fgm5bawFUgtHI+QtIdkSnrxSL8oiIfHDrbut09kUt54Zff9WsqsfPjBSQTAgEcNxy20yxnFkQETWb37MxvDjNOuzahBSKinRnMhKdaxksiIh8aOTZGMaOp/vWLAFQahHi1OYjN1LAtve/uaqTrbWzbTNwGYqIyAen5aN6z8awnvDnxGi1bl3S6h/MNvV4V84siIh8aHR3W7sciJVbMOrtzrhWhDcaZxZERD41srutW65DUJpRbFi5wPX9Nq9eNG52kkwIoEBhtLFpcQYLIqIQXJxO2eYpMukUDm5c4esaRiDZtvcEfpXLVwKM+bFGhQwGCyKiEPjJgfQPZscFAj8zG/MM6PJNj6KoEw8ZDBZE1Lbq/TBuBLtZwfKFM7Ft7wms7xtCujOJl89dqCwnGUWA5tdak+TGcw4/exoHjp/izIKIaKKcPmgBNDVgOH3wnxkZf3a3tQjQqVBwx6GTDS/ki+VuKPaGIqKJcqvIjsp47Pwql6/0qHKqzQii4juWwYK9oYhoohpVkd0oft93WipZ6VHVTLEMFkREE9XIiuxG8PO+qWQCInCdgQTVGITBgojaklsb8aiMJ5kQpFPJqiLAnE0uw5BJp3Dr0jmBHM3KBDcRtSWnGoVmJbfrHc+2vSc86zN6Lp2OLQPDtknyejFYEFHbamRFdiP4GY+f+gzjOv2DWfugoTpa69i4DEVEFCO19Kjq7c5g8DPvwf1rllQ9/8LvTj1b6/tyZkFEFDMTnRHJ5NRFtb6GwYKIqIXZFR8mUq+dWet1uAxFRNTC/Bb7eWGwICJqYY0qMmSwICJqYY0qMmSwICJqYXbFfvVgsCAiamF2W22L+d+fqvU63A1FRNTirFttZdNvTtZ6Dc4siIjIE4MFERF5YrAgIiJPkclZiMhUAP8M4DyA76nqjpCHREREZYHOLETkARF5UUSesDx+rYicEJGnRWRj+eEbAXxLVW8HsDrIcRERUW2CXob6GoBrzQ+ISALAlwBcB+BKALeIyJUALgHwXPlpE69NJyKihgk0WKjq9wGctjz8VgBPq+ovVPU8gG8AuAHA8ygFjMDHRUREtQkjZ5HB2AwCKAWJawB8EcA/icj1AAacXiwiawGsLX/5qnWJq43NAPBS2IOICN6LMbwXY2bsvIP3oqzms2Mjk+BW1VcA/KWP520HsB0AROSwqvYEPbY44L0Yw3sxhvdiDO/FGBE5XOtrwljuyQKYbfr6kvJjREQUUWEEi8cBXCEi80RkMoAPANgdwjiIiMinoLfOPgzgRwAWiMjzIvIRVb0A4GMA9gJ4EsBOVR2u8y22N2iorYD3YgzvxRjeizG8F2NqvheiqkEMhIiIWgi3qBIRkScGCyIi8hSLYCEis0XkgIj8VESGReQ/lh+fLiL/R0SeKv+/K+yxBs3lXvxZ+etREWmL7YEu92KbiBwXkZ+IyLdFJB3yUAPnci/+tnwfhkRkn4hcHPZYg+Z0L0zfv1NEVERmhDXGZnH5udgsItnyz8WQiLzX81pxyFmIyCwAs1T1xyLyWgBHAPQC+BCA06q6tdxjqktV/ya8kQbP5V4ogFEAXwHwCVWteR913Ljci0sA7FfVCyLy9wDQxj8Xz6vq78rP+Q8ArlTVO8IbafCc7oWq/lREZgP4bwAWArhaVVu6SM/l5+JmAC+r6n/2e61YzCxU9QVV/XH5z79HaRdVBqU2IQ+Wn/YgSjehpTndC1V9UlVPhDu65nK5F/vKu+4A4BDG2si0LJd78TvT06ai9EtFS3P5vACA+wB8Em1wHwDPe1GTWAQLMxGZC6AbwGMA3qiqL5S/9WsAbwxrXGGw3Iu25nIvPgzgO00fUIis90JEPicizwG4FcBnQhxa05nvhYjcACCrqkfDHVU4bP6NfKy8RPmAnyX8WAULEbkIwCMA1ll+Y4KW1tPa4rcFwP1etBuneyEidwG4AKBtzkaxuxeqepeqzkbpPnwszPE1k/leoPRz8Cm0WbA02Pxc/AuAywEsAfACgM97XSM2wUJEkij9ZXeo6q7yw78pr8kZa3MvhjW+ZnK4F23J6V6IyIcAvA/ArRqHxFwD+Pi52AHgpuaOKhw29+JyAPMAHBWRZ1BamvyxiPxBeKNsDrufC1X9jaoWVXUUwH9FqRu4q1gECxERAF8F8KSqfsH0rd0A/qL8578A8D+bPbZmc7kXbcfpXojItSitS69W1ZGwxtdMLvfiCtPTbgBwvNljaza7e6Gqx1T1Dao6V1XnotTt+o9U9dchDjVwLj8Xs0xP+1MAnt2747Ib6h0AfgDgGEo7foDSlPIxADsBzAHwLICbVdV6fkZLcbkXUwD8FwAzAeQADKnqyjDG2Cwu9+KLKN2P35YfO9QGO4Cc7sVHUGpHPYrSv5E7VLWlG3c63QtVfdT0nGcA9LTBbiinn4tbUFqCUgDPAPioKf9rf604BAsiIgpXLJahiIgoXAwWRETkicGCiIg8MVgQEZEnBguiCBKRd4rI28MeB5GBwYIomt4JgMGCIoPBgtqWiMwttzL/moj8TER2iMi7ReRgue39W8v//UhEBkXkhyKyoPzaThHZWW79/G0ReUxcWsOLyL+IyOFym+gtpsefMVpli0iPiHyv3MPnDgDry+2j/7g81v3lXj7fFZE5Ad8eoiqTwh4AUcjmA/gzlBoOPg7ggwDeAWA1SsVL/xbAH5fbnb8bwN+h1DLjrwCcUdUrReRNAIY83ucuVT0tIgkA3xWRq1T1J3ZPVNVnROTLMLWQFpEBAA+q6oMi8mGUCg97J/IXJ6oFgwW1u1+q6jEAEJFhAN9VVRWRYwDmApgG4MFy2wwFkCy/7h0A/hEAVPUJEbH94De5WUTWovRvbhaAKwF4vcbsbQBuLP/5fwD4hxpeSzRhXIaidveq6c+jpq9HUfpg/1sAB1T1TQBWAXhNrW8gIvMAfALAu1T1KgB7TNe5gLF/hzVfm6hZGCyI3E0DYPRS+pDp8YMonTYGEbkSwGKXa7wOwCsAzorIGwFcZ/reMwCuLv/Z3BH29wBea/r6hwA+UP7zrSj1+yFqGgYLInf/AOBeERlE9bLtPwOYKSI/BXAPgGEAZ+0uUD5sZxCljq9fRynQGLYA+EcROQygaHp8AMCfGgluAB8H8Jfl5a4/B1B1rjRR0NhIkKgO5UR1UlXPicjlAP4vgAWqej7koREFggluovp0AjhQPlhGAPwVAwW1Ms4siBpIRB5D6SwNsz83dlwRxRWDBREReWKCm4iIPDFYEBGRJwYLIiLyxGBBRESeGCyIiMgTgwUREXn6/8tBTjQNz5DgAAAAAElFTkSuQmCC\n", + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZEAAAEHCAYAAABvHnsJAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAiY0lEQVR4nO3df5RcZZ3n8fe3K5XQHTSdDKxCk5gInDhIkJY+gEZ3wV/BxUCLCkac1dElsnNwF3/EDSM7gIOTjFl/7jhqRjngksUgMj2JwRN2BQcPQxiC3RgCRJCfaXGIJh003ZBK93f/qLrdt6vr5+2uulW3Pq9zOHbduvfW46Xobz/P93m+j7k7IiIiUbTF3QAREWleCiIiIhKZgoiIiESmICIiIpEpiIiISGQKIiIiEtmsuBtQC8ccc4wvXrw47maISIN7ct8hAF577NyYWxK/Bx988Hfufmy11yUyiCxevJidO3fG3QwRaXCXfOc+ADZ/4k0xtyR+ZvZMlOsSNZxlZivNbOPBgwfjboqISEtIVBBx963uvnrevHlxN0VEpCUkKoiIiEh9KYiIiEhkCiIiIhJZooKIEusiIvWVqCCixLqISH0lKoiIiEh9KYiIiEhkCiIiIhKZgoiIiESmICIiIpElKohoiq+ISH0lKohoiq+ISH0lKoiIiEh9KYiIiEhkCiIiIhKZgoiIiESmICIiIpEpiIiISGQNH0TM7E/N7NtmdpuZ/Ze42yMiIhNiCSJmdoOZvWBmD+cdP8/M9pjZE2a2FsDdH3X3y4GLgeVxtFdERAqLqydyI3Be+ICZpYBvAu8GTgFWmdkpufcuALYBd9S3mSIiUkosQcTd7wH25x0+E3jC3Z9098PAD4ALc+dvcfd3A5fWt6UiIlLKrLgbENIFPBd6vRc4y8zOAS4C5lCiJ2Jmq4HVAIsWLapZI0VEZEIjBZGC3P1nwM8qOG8jsBGgp6fHa9sqERGBxpqdNQgsDL0+IXesYqriKyJSX40URB4ATjazJWY2G/ggsKWaG6iKr4hIfcU1xfcW4D5gqZntNbOPu/sR4ApgO/AocKu7767yvuqJiIjUUSw5EXdfVeT4HUxjGq+7bwW29vT0XBb1HiIiUrlGGs6aNvVERETqK1FBRDkREZH6SlQQERGR+kpUENFwlohIfSUqiGg4S0SkvhIVREREpL4SFUQ0nCUiUl8NXzurGsXWifT1D3Ltlt0MjWQAmN+R5pqVr6e3uyuOZoqIJEaigkghff2DrPnhQ2TGJmoyHhjOsOa2hwAUSEREpiFRw1mBXYMH6f7CneM9kHAACWRGnQ3b98TQOhGR5EhUT8TMVgIrZ3Uex4HhDFduHih5/m+GRurSLhGRpEpUTySY4tt21NyKzu/sSNe4RSIiyZaonki1DgxnWLx22/jrubNTjLkzkhkDsgn48087jrsf28dvhkY4vrOdNSuWKo8iIpLT0kEk36HDo5NeHxjOcPOOZ8dfDw6NsOaHSsiLiAQURKqUGXOu3DzAVbf/kpHMGCkzRt3pUi9FRFpQooJIOLFea8GQ16hnZ36plyIirShRQSRYbDjnuJNj2ZQq6KVAZYGkr3+QDdv3KN8iIk0rUUGkUVy5eWDS9GIDHGgzCJasdKTbeCkzxljuHPVkRKQZJWqKb6MKljqG1zwOhwJIIDPmXLulqm3lRURipZ5IgxkayU47LpeoVz0wEWkE5j61JEizm3PcyX7cR74WdzNqZn5HmoPDmSk9mVSb8Yo5szg4klGORaQCl3znPgA2f+JNMbckfmb2oLv3VHudeiJN6MBwpuDx0TEf75kMDo1w1e27AOVYRKR2EpUTCfYTGXvpUNxNaQgjmVGu3DzA8vV30dc/GHdzRCSBEhVEqq2d1SoGh0b49OYBBRIRmXEazmoRY0ydejyrzTiSmzIWTEPWynsRqYaCSAs7EppzHPyUn0vRgkgRKUVBRKYYyYyyYfsedj6zn007ni0aYEREEpUTkZkzODTCzaEAEhjJjHLdVi2IFJEs9USkauF9WJRLEWltDR9EzKwXOB94JfA9d78z3hZJWHioKz9xr1X0IskXy3CWmd1gZi+Y2cN5x88zsz1m9oSZrQVw9z53vwy4HLgkjvZKNME+94vXbuP06+4cn2Lc1z/I8vV3sWTtNq1hEWlycfVEbgT+Dvh+cMDMUsA3gXcCe4EHzGyLuz+SO+Xq3PvShIZGMlN6KqBkvUiziyWIuPs9ZrY47/CZwBPu/iSAmf0AuNDMHgXWAz9x91/Ut6VSD8HK+p3P7Of63mWA9loRaRaNlBPpAp4Lvd4LnAV8EngHMM/MTnL3bxe62MxWA6sBUq88tsZNlVq4ecezk/a0DwwOjbDmNu21ItKIGimIFOTu3wC+UcF5G4GNkK3iC5BOGbNTbRw6PFrbRkrNZUZ90nDYh89eRM9rFqi3IhKzRlonMggsDL0+IXesYuECjF2d7Wx4/xsYVgBJpJt3PMuVmwcYHBrBUX0wkbjEtp9ILifyY3c/Nfd6FvAr4O1kg8cDwIfcveqVbT09Pb5z504Alq+/i8GhkZlqtjSJjnQbF51xAnc/tk89FSlK+4lMaKr9RMzsFuAc4Bgz2wtc4+7fM7MrgO1ACrih2gBiZiuBlSeddNL4sTUrlnLV7bsYyVTWI0m3GaPuk7ayleYznBmblF8ZHBrhU3mzw1JmrDpr4XgyX0Sql8idDcM9EZg806fNskGikGDVNcBnbn2o6HmSLHNnp/jie5epl9KC1BOZELUnkqggEuqJXPb4448XPKevf3BKz6Q9nWLdRZN/iSxZu21K3ShpDeqhtA4FkQlRg0gjJdanLdiUat68eUXP6e3uYt1Fy+jqbMfI9j7yAwjA8Z3tNW6tNKpRd27e8SyL125j8dptXPoP98XdJJGG1fBTfGuht7ur7NDFmhVLWfPDh8goOdLy7v31/vGCk23AvI40Q8MZJetFSFgQKZRYjyr4xXDtlt0MjWSAbEFByNaEKtsWKDkc1mYw5pOr3y69+ie8fGRsuk2XGhpj4t9/UHTyM7cOsOqsRZoJJi0pUUHE3bcCW3t6ei6bifsV6rEUyqkUbAvZADE4NDIeMAA629Nce8HUyrZ9/YMKIE1q1JkyEyy/jItIUiUqiNRD8Mu/3Gyvrs527l37torvu2H7nhlrozSGcBmX5ScuYNNlSt5K8iQqiMzkcFYp4R5KsdlewVThSv1GCyITLZxXKaRYD1Wk0bXc7KyZVulsr3KKzQbrSLeN52KC1+3piX9t8zvSfPjsRaTMIrVfGkNQKv+dX/lZ3E0RqUqieiJxqWS2VzmFVta3p1P8TYUBaVOB6rfSfB5/4dCkHsvyExfw9O9HlLCXhpWonkgzm26PRutakuneX++fVGQy2Cmy+wt3qtikNISWW7GeVJXMGutsT49PV5bkMAPPmy4uldGK9QkqexKSXzurVfT1Dxat+RXMFlNV49aiwFKagsgElT0Reru7+PLFb6A9nZp0PDxbbM2KpVPez9fV2c7XLjmd5ScuqFlbpT6CIbAla7dp+EtqQon1hMlfx5KfjA3+97qtu6esvM8vRNnb3TVlr/NzX3csP3pwsOLS+tIYHCbtDBmm6cUyHRrOamH5AaLUsEf43HntaQ6+lCGBX52WNmdWGx/oyW7kVWmlhWan4awJTbUpVa3Ua7FhUlQ6NTk/aT80kiGdMnAqKlBZro6YNIaXj0zeyCv8r3ZoJMOnc72YpAUSmZ5EBZGZrp0lWRu275kyfJUZdeZ3pHlx5EjBRH7KjC9f/IbxXzhX9+2a9AtKms8YU4fElLgXJdalrGIlWYaGM4wVGdMac5/0i+X63mV8+OxFVLOuvk2L8BteeO3K1X274m6OxCBRPRGpjeNz1YgLHQdKvhd2fe+ySVVtT7zqjqJbEM/vSPPHl44UDVIaIms84YKT+eZ3pLlmZfJyKqKeiFSg0LTgYNpwqffKKbWHfcfsWUXzLV2d7Vx69qIKWi6N4sBwRj2WhFJPRMoqN2243HvFdBXp4XR1thcdQjMYL7G/7ZfPV7RBmDSWQj2WcG6lmlmDEj9N8ZXYFCujv+6iZWzYvqdogAmCSKHr020Glk38S3LUaoqxpvhO0BRfNMW32ZTr4ZTbp6XY9cGx/LUOxczvSKtH0+CCUvlXbh4Y77VAtB6wzCz1RKRhzeSwRrGaYcHQWfL+K2g9Blx69qKqtiRWT2SCCjCGKIhIvihDZ2HF8jfS+AqtZQn+QBkcGmF2qo0vvf+0lu/FKIiEKIi0rlK9l2LvlSujHwQbgE9vHmAs9F4b8JVLTgcmhlZA048bXf4U8fZ0ived0cXdj+1r2eExBZEQBZHWVKq3Ue6XQfgv00KChH4lQ2xL1m5TEEmASr87SaEgEqIg0ppK5T2CGV3lFAsABjy1/vxptUOaVyssltTsLGl5xdaWFDteSLnV+ZVYs2Jp2V0mpbkEiyXzS+lrfcs0goiZLXJ3VdSThlGrAFDpCvxAeOqxeiTJFtQOyw8ug0MjXHX7LnY+sz/xeZayw1lm9iagC7jH3V8ws9OAtcBb3X1hzRto9lrg88A8d39/JddoOKs1TScnkn+fmfqLsq9/sOBGUIFgf3RpHcG+LT9+6HmGRrLrk+Z3pDn/tONiDTg1yYmY2QbgPcAAcBKwHfjPwDrgO+7+UsTG3pC77wvufmro+HnA14EU8F13Xx967zYFESmnUYYUwu2gRKD48NmLtFOkFFTvxH6tgsgjwBvd/SUzmw88B5zq7k9Hbmn2vv8e+CPw/SCImFkK+BXwTmAv8ACwyt0fyb2vICJNodyUYZi8MC4ccNrMShamlNbWnm7jqHSKoeHMjP+RVKsg8gt3f2Podb+7d0dsY/69FwM/DgWRNwHXuvuK3OurANx9Xe51xUFkwWv+1N/5lzfMRDNFqtb/7BCHR8eKvj871cbCBe0cc/ScKe/d/9T+WjZNEqzU96oSt17+5prMznqtmW0JvV4Sfu3uF1T7gSV0ke3pBPYCZ5nZnwBfBLrN7KogqOQzs9XAaoCjjztxBpslUp1SAeSsJQtKXjs71Vb0+lltxpEKtiOW1nR4dIxf7zvEr/cdAhivG5dqMww4MubTDjSFlAsiF+a9/vKMfXKF3P33wOUVnLfRzJ4HVr56rp2hWjgSl1LrVcp9LwtWJq5iP3uRQPB1GQ19b/IDTbByP2VGat6rIm3SUzKIuPs/R7lpRINAeLbXCbljFdMe69IIpjNNuFBl4kMvHxmfxSMyk4LwMupOqv0Vx0a5R8kgYmZ3U7wMkLv726N8aBEPACeb2RKyweODwIequYFKwUsjqGQTr3LXh89dsnZbTdopMhPKDWd9tsCxs4HPAS9E/VAzuwU4BzjGzPYC17j798zsCrLTiFPADe6+u5r7qicijSI/EExHsUWUIo2g3HDWg8HPZvYfgP8BHAVc7u4/ifqh7r6qyPE7gDui3lckKcLTfjs70qTbbFJOJBjLDlejnTs7xeEjY8qdSF21lTvBzFaY2c/JBpAvuvtbphNAasnMVprZxoMHD8bdFJHIguT6YG6zrAPDGbDsFrFGNkF/6dmLaE+nJo01jzlccuZCuqoo8yIyXSWDiJk9AHwH+AHZIayDZvbG4J96NLAa7r7V3VfPmzcv7qaIRLZh+54pCxUzo87cObN4av353Lv2bdz92L4p54xkRrn7sX3cu/ZtWBWf15Fuq+p8kbByOZFDZFeWvz/3T34/ubL62nWixLokQSXViMudU00eZTgzxofPXsSmHc9qHxSpWrnhrM8BH3L3c939XOAmskHlYbJBpaGoJyJJUKzqcPh4uXPWrFhKezpV8Wde37uMr15y+vhQWMrUN5HKlAsi3wZehvF6V+vIBpKDwMbaNk2kNRUKAPnrTMqd09vdxbqLltHV2T6eR5k7u3hQ6f7CnZPuq/pdUqlytbMecvc35H7+JrDP3a/NvR5w99Pr0chKhYazLnv88cfjbo5IZJVUI662YnG5svTt6RRzZrVVtbCxq7N9fAbZgWEtiGxmz990JS8//3jVXdByQeRh4HR3P2JmjwGr3f2e4L1wGfdGoiq+IoUtnsGFi/nbDs/kvaX+ogaRcsNZtwD/bGb/BIwAPwcws5PIDmmJSBOZqem/hcq4aGpxayoZRNz9i8BngBuBt/hEt6UN+GRtmyYiM61cwn1+R7psQr6rs73gZknVJvMlGcruse7uOwoc+1VtmjM9muIrUlrwi//aLbun5D7a0ymuWfl6AK7buntKjqPcTnvB8ULXSnKVXbHeTDTFV6S83u4uBq55F1/LTekNZm8FAaK3u4v+v8q+P78jPX7dnFnlf12Er9XwVmso2xMRkWSqpEjkS5mJDbKGRjJcdfuu8WsruXdf/yCf2jygRYwJlqieiIjMnELlV0Yyo2zYvqfie/R2dymAJJx6IiItqtw6k0rKr1SiK2Ip+3CFYmlcieqJqIqvSGXyKwUPDo1w1e276Ouf2Ey0kvIrlYg6a8vJBqCvXXL6lPzN8hNL71Uv9ZOoIKLEukhlKhmqqqT8SiXCJViqFQQ3gHvXvo2n1p/Pua87lnt/vb/qe0ltJCqIiEhlKhmqKlR/q9QU31J6u7uqLlEfyA9ut9z/XIS7SK0oJyLSgoqVis8fqprJbX5LfW65/Ec4uKk4ZGNRT0SkBc3UUNVMfe6lZy8qOdwVDm61LFNvULLasUylnohICwp6F9VUAa7H5wYJ/3C+Jj+4rTprITfveLYm7XOy2wx3tqerqmbcykpW8W1WquIr0rwqKXF/dd8uNt3/LPm/vtrTKY5Kt1VUdmX5iQuKJujnd6R5KTM2ZfJBktWkFHyz0X4iIq2lUMABpvRmUm3G2JjjZIfDVp21kOt7lxUtX2/ApS22ZXDUIJKo4Sx33wps7enpuSzutohI7ZVK/FcyVDc71cbh0bEpx81g878+NyWAmIG7FkKGJSqIiIhA5bPKFi5o56nfHWIsLyKMOYwVGKU5fl52I66r+3Zxy/3PaaYYCiIi0sKOOXoOAE//briigPCboRH6+gf50YODCiA5CiIi0hKK5U+OOXoOT+47VNE9ju9sL7jav5VpnYiIJF6xWmG/++PLQGX1wIKpxpUWoJydqt16lkaiICIiiVesVthz+7MBodAiyHTK6GxPTyn5Mq89TTnpNsiMOl2d7XRWcH4z03CWiCResd5DMDOr0sWXff2DHDp8pOznBXt5DQ6NkE4Z6TYjk5+9T4iGDyJmNhf4e+Aw8DN33xRzk0SkyRSr2TU7NTEYU8mMrg3b95AZrS4YZEadjnQbR3LrVJImluEsM7vBzF4ws4fzjp9nZnvM7AkzW5s7fBFwm7tfBlxQ98aKSNMrVrNr4YLqytNXuyFXYDgzlsgAAvHlRG4EzgsfMLMU8E3g3cApwCozOwU4AQhqP2tKhIhUrVhZ+2CKb6VKJeC7OtvpSFf+KzUpafdYhrPc/R4zW5x3+EzgCXd/EsDMfgBcCOwlG0gG0EQAEYmo0HDVLf9aXSHHNSuWTimpkk4Zc2fPqnoL4KT0TBrpl3IXEz0OyAaPLuB24H1m9i1ga7GLzWy1me00s5379u2rbUtFpCXl92jmd6TBaemKvw2fWHf3Q8CfV3DeRmAjZKv41rpdItKawj2a5evvqqhicJI1Uk9kEFgYen1C7ljFzGylmW08ePDgjDZMRKSQqIn2KNrTbQ2ZR2mkIPIAcLKZLTGz2cAHgS3V3MDdt7r76nnz5tWkgSIiYeVWuqfbLDvkxfR2ZGxPp1h30Wl89ZLTS+4AGYdYhrPM7BbgHOAYM9sLXOPu3zOzK4DtQAq4wd13V3nfYD+RmW6yiMgUhRLtQZn4rgILFgvt3FhMyoxR9yn36e3uKroPShzimp21qsjxO4A7pnFf7SciInVT7TbD+ed3dqR5OTPKcG6Je2d7mmsveH3ZRY9dRRZPxqHhE+vVUE9EROqt0r1Lop5fSKEeUFwaKScybcqJiEgrCKYaTyfPMlMSFUQ0O0tEWkVvd1fB3RfrLVFBRD0REYlbX/8gy9ffxZK121i+/i76+qtaqVCVSvZBqbVE5UREROKUP/sq2PwKJpLqhXZYrCZHEr6+syMde5n5RPVERETiVGzzqw3b9wDFd1istLdydd8uPrV5YPz6A8MZMCZtnhWsS6mXRPVENDtLROJUbAV7cLxUkMlfT1JoP/hNO56dUrgxM+rMnTOLgWveBVD3NSSJCiJaJyIicSq2+VWQuygXZKD4kNhR6bailX+D6/v6B8cXO9aLhrNERGZIsc2vgp5EsUR4+Hix3kqpQo/B9Ru276kqgMzEBOFEBRFN8RWROBXb/CoYqioXZKD6oo6Wu2+Ua2eix6LhLBGRGVRqRXolZVKKDYl1tqd5+cjYlDpdbz5xARu27+FTmwdoy9XbqqdEBRERkUZXruxJoZIm7ekU117wemByADr3dcfyowcHx8+tdwABBRERkYZSrrcSDkDL199VsH5Wyowx97r0TBREREQaTKVFGovlQMbceWr9+QVLz7enU8yZ1TZjW/oqsS4i0qTKzfYqlui/9oLXT0nwR5WonogS6yLSSorlT8KzvUr1ajZs38Pg0AjTKQacqCAiItJKqt0Uq9C1092XREFERKSJTWeTq0ILG6uVqJyIiIhUrtrFiYUoiIiItKiZ2I9EQUREpEWtWbF02vWzEpUTUSl4EZHKN77q7e5i5zP72bTj2ciflaieiLbHFZFWV2jjqys3D9D9hTsLbn51fe8yvnrJ6ZGrMSaqJyIi0uqKzbg6MJwZ36o3OC/cUzny4gtPRfk8BRERkQQpNeNqJDPKdVt381JmrOA+8FEkajhLRKTVlZtxdWA4U3DTq9TRCyItNlEQERFJkEIbX1XCUrNmR/k8BRERkQQJii52tqenvNeeThU8DuCjRw5H+TwFERGRhOnt7mLgmnfxtUtOr6iCb3s6xegf90+dulWBhk+sm9lrgc8D89z9/XG3R0SkWZSr4BuenfXe61/cH+UzahpEzOwG4D3AC+5+auj4ecDXgRTwXXdfX+we7v4k8HEzu62WbRURaRXTKdqYr9Y9kRuBvwO+HxwwsxTwTeCdwF7gATPbQjagrMu7/mPu/kKN2ygi0vLa2l+5IMp1NQ0i7n6PmS3OO3wm8ESuh4GZ/QC40N3Xke21iIhIHfX1DzLrlce+Jsq1cSTWu4DnQq/35o4VZGZ/YmbfBrrN7KoS5602s51mtnPfvn0z11oRkYTbsH0PmEWKBw2fWHf33wOXV3DeRmAjQE9PT8QqMCIirWc6+4rE0RMZBBaGXp+QOzZtZrbSzDYePHhwJm4nItISprOvSBxB5AHgZDNbYmazgQ8CW2bixqriKyJSvTUrloL7WJRraxpEzOwW4D5gqZntNbOPu/sR4ApgO/AocKu7756hz1NPRESkSr3dXRx5cd8zUa6t9eysVUWO3wHcUYPP2wps7enpuWym7y0ikmRjI9EWGyaq7Il6IiIi9ZWoIKKciIhIfSUqiIiISDRRV6wnKohoOEtEpHrNtmK9ZjScJSJSvemsWE9UEBERkeo124r1mtFwlohI9ZptxXrNaDhLRKR6DbtiXUREGt90VqwriIiIiFasg3IiIiL1lqggopyIiEh9JSqIiIhIfSmIiIhIZAoiIiISWaKCiBLrIiL1laggosS6iEh1+voHWb7+Lma/+qQzolxf050NRUSkcfX1D3LV7bsYyYxGvkeieiIiIlK5Ddv3TCuAgIKIiEjLmk713oCCiIhIi5pO9d5AooKIZmeJiFRuzYqltKdT07pHooKIZmeJiFSut7uLdRcto0v7iYiISBS93V3cu/ZtHP7tEw9GuV5BREREIlMQERGRyBREREQkMgURERGJzNw97jbMODP7A7An7nY0iGOA38XdiAahZzFBz2KCnkXWUnd/RbUXJbV21h5374m7EY3AzHbqWWTpWUzQs5igZ5FlZjujXKfhLBERiUxBREREIktqENkYdwMaiJ7FBD2LCXoWE/QssiI9h0Qm1kVEpD6S2hMREZE6UBAREZHImjqImNlCM7vbzB4xs91m9t9yxxeY2f81s8dz/zs/7rbWWoln8YHc6zEza4lpjCWexQYze8zMfmlm/2hmnTE3teZKPIu/zj2HATO708yOj7uttVbsWYTe/4yZuZkdE1cb66XE9+JaMxvMfS8GzOw/lr1XM+dEzOw44Dh3/4WZvQJ4EOgFPgrsd/f1ZrYWmO/u/z2+ltZeiWfhwBjwHeCz7h5pLngzKfEsTgDucvcjZva3AC38vdjr7i/mzvmvwCnufnl8La29Ys/C3R8xs4XAd4HXAWe4e6IXH5b4XlwM/NHd/2el92rqnoi7P+/uv8j9/AfgUaALuBC4KXfaTWQfTqIVexbu/qi7t9Tq/RLP4k53P5I7bQfZoJJoJZ7Fi6HT5pL9YyPRSvy+APgq8Dla4DlA2WdRlaYOImFmthjoBu4HXuXuz+fe+i3wqrjaFYe8Z9HSSjyLjwE/qXuDYpT/LMzsi2b2HHAp8FcxNq3uws/CzC4EBt39oXhbFY8C/41ckRvqvKGSVEAigoiZHQ38CLgy7y8sPDte1xJ/XUDpZ9Fqij0LM/s8cATYFFfb6q3Qs3D3z7v7QrLP4Yo421dP4WdB9nvwl7RYEA0U+F58CzgROB14HvhyuXs0fRAxszTZh7DJ3W/PHf633JhfMPb3Qlztq6ciz6IlFXsWZvZR4D3Apd7MCcEqVPC92AS8r76tikeBZ3EisAR4yMyeJjvE+Qsze3V8rayPQt8Ld/83dx919zHgH4Azy92nqYOImRnwPeBRd/9K6K0twEdyP38E+Kd6t63eSjyLllPsWZjZeWTHvS9w9+G42ldPJZ7FyaHTLgQeq3fb6q3Qs3D3Xe7+79x9sbsvBvYCb3T338bY1Jor8b04LnTae4GHy96rmf8YM7O3AD8HdpGdgQTZrun9wK3AIuAZ4GJ33x9LI+ukxLOYA/wv4FhgCBhw9xVxtLFeSjyLb5B9Hr/PHdvRAjOSij2LjwNLc8eeAS5398FYGlknxZ6Fu98ROudpoKcFZmcV+16sIjuU5cDTwCdC+eXC92rmICIiIvFq6uEsERGJl4KIiIhEpiAiIiKRKYiIiEhkCiIiTcLMzjGzN8fdDpEwBRGR5nEOoCAiDUVBRCTEzBbnysXfaGa/MrNNZvYOM7s3t7XAmbl/7jOzfjP7FzNbmru2w8xuzZXX/kczu99KlN83s2+Z2c5cKe7rQsefDsqRm1mPmf0sV9/ocuBTuRLdb8219a5cnaOfmtmiGj8ekSlmxd0AkQZ0EvABskUaHwA+BLwFuIDsgqz/BLw1V1L+HcDfkC0b8hfAAXc/xcxOBQbKfM7n3X2/maWAn5rZae7+y0InuvvTZvZtQmW6zWwrcJO732RmHyO7mLJ3Ov/HRaqlICIy1VPuvgvAzHYDP3V3N7NdwGJgHnBTrnSIA+ncdW8Bvg7g7g+bWcGAEHKxma0m+9/hccApQLlrwt4EXJT7+X8DX6riWpEZoeEskaleDv08Fno9RvYX/l8Dd7v7qcBK4KhqP8DMlgCfBd7u7qcB20L3OcLEf5tV31uknhRERKo3DwjqTH00dPxesjvDYWanAMtK3OOVwCHgoJm9Cnh36L2ngTNyP4er6/4BeEXo9b8AH8z9fCnZWkgidaUgIlK9LwHrzKyfyUPCfw8ca2aPANcDu4GDhW6Q2wCpn2z13P9DNgAFrgO+bmY7gdHQ8a3Ae4PEOvBJ4M9zw2Z/BkzaM1ykHlSAUWSG5BLkaXd/ycxOBP4fsNTdD8fcNJGaUWJdZOZ0AHfnNvsx4C8UQCTp1BMRqTEzu5/sPiZhfxbMABNpZgoiIiISmRLrIiISmYKIiIhEpiAiIiKRKYiIiEhkCiIiIhKZgoiIiET2/wHNduEQPF9pVAAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] @@ -160,12 +126,13 @@ }, { "cell_type": "code", - "execution_count": 5, + "execution_count": 7, + "id": "8fa76f9d", "metadata": {}, "outputs": [ { "data": { - "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXAAAAEHCAYAAAC3Ph1GAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/d3fzzAAAACXBIWXMAAAsTAAALEwEAmpwYAAAQGElEQVR4nO3db4xldX3H8fenC1artgiMmw0LXRoJhlhd7AQxGqugditE9gHZSCxZ2zXbB9pgaqOrPiA1NVmb1D9tapMNWOcBChuULkFjJSsEmjTorFBBkIpkibvdZYcKUftAsvLtg3umTmdnd+7MvXdmf3fer2Rzz/ndc2a+v3D3wy/fPefcVBWSpPb8xmoXIElaHgNckhplgEtSowxwSWqUAS5JjTLAJalRZ/RzUJKzgJuA1wAF/BnwOHAbsAk4CGyrqmdP9XPOPffc2rRp07KLlaS16MCBA89U1cT88fRzHXiSKeD+qropyYuA3wI+Dvy0qnYn2QW8oqo+eqqfMzk5WdPT08ubgSStUUkOVNXk/PFFWyhJfgd4C3AzQFU9X1XPAdcAU91hU8DWYRUrSVpcPz3wC4EZ4J+TPJjkpiQvBdZX1ZHumKPA+lEVKUk6UT8BfgbweuCfqupS4H+AXXMPqF4fZsFeTJKdSaaTTM/MzAxarySp00+AHwIOVdUD3f7t9AL96SQbALrXYwudXFV7qmqyqiYnJk7owUuSlmnRAK+qo8BPklzcDV0JPArcCWzvxrYD+0ZSoSRpQX1dRgj8BXBLdwXKk8Cf0gv/vUl2AE8B20ZToiRpIX0FeFU9BJxwCQu91bgkaRV4J6YkNcoAl6RG9dsDX1M27fr6Sd87uPuqBY+bOy5JK8EVuCQ1ygCXpEbZQumcqm0iSacjV+CS1CgDXJIaZYBLUqMMcElqlAEuSY0ywCWpUQa4JDXKAJekRhngktQoA1ySGmWAS1KjDHBJapQBLkmNMsAlqVFr+nGyPkJWUstcgUtSowxwSWqUAS5JjTLAJalRBrgkNcoAl6RG9XUZYZKDwM+BXwHHq2oyydnAbcAm4CCwraqeHU2ZkqT5lrICf1tVba6qyW5/F7C/qi4C9nf7kqQVMkgL5RpgqtueArYOXI0kqW/9BngB30pyIMnObmx9VR3pto8C6xc6McnOJNNJpmdmZgYsV5I0q99b6d9cVYeTvBK4O8kP575ZVZWkFjqxqvYAewAmJycXPEaStHR9rcCr6nD3egy4A7gMeDrJBoDu9dioipQknWjRAE/y0iQvn90G3gk8AtwJbO8O2w7sG1WRkqQT9dNCWQ/ckWT2+C9X1TeTfBfYm2QH8BSwbXRlSpLmWzTAq+pJ4HULjP83cOUoipIkLc47MSWpUWv6Cx2Wwy+BkHS6cAUuSY0ywCWpUQa4JDXKAJekRhngktQoA1ySGmWAS1KjDHBJapQBLkmNMsAlqVEGuCQ1ygCXpEYZ4JLUKANckhq15h4n6+NgJY0LV+CS1CgDXJIaZYBLUqMMcElqlAEuSY0ywCWpUQa4JDXKAJekRhngktQoA1ySGtV3gCdZl+TBJHd1+xcmeSDJE0luS/Ki0ZUpSZpvKSvwG4DH5ux/GvhsVb0KeBbYMczCJEmn1leAJ9kIXAXc1O0HuAK4vTtkCtg6gvokSSfR7wr8c8BHgBe6/XOA56rqeLd/CDhvoROT7EwynWR6ZmZmkFolSXMsGuBJrgaOVdWB5fyCqtpTVZNVNTkxMbGcHyFJWkA/zwN/E/DuJO8CXgz8NvB54KwkZ3Sr8I3A4dGVKUmab9EVeFV9rKo2VtUm4D3At6vqvcA9wLXdYduBfSOrUpJ0gkGuA/8o8JdJnqDXE795OCVJkvqxpK9Uq6p7gXu77SeBy4ZfkiSpH96JKUmNMsAlqVEGuCQ1ygCXpEYZ4JLUKANckhplgEtSowxwSWqUAS5JjTLAJalRBrgkNcoAl6RGLelhVurPpl1f/7/tg7uvWsVKJI0zV+CS1CgDXJIaZQtlSOa2TSRpJbgCl6RGGeCS1CgDXJIaZYBLUqMMcElqlAEuSY0ywCWpUQa4JDXKAJekRhngktSoRW+lT/Ji4D7gN7vjb6+qG5NcCNwKnAMcAK6vqudHWWyLTnaLvU8plDSoflbgvwSuqKrXAZuBLUkuBz4NfLaqXgU8C+wYWZWSpBMsGuDV84tu98zuTwFXALd341PA1lEUKElaWF898CTrkjwEHAPuBn4MPFdVx7tDDgHnjaRCSdKC+grwqvpVVW0GNgKXAa/u9xck2ZlkOsn0zMzM8qqUJJ1gSVehVNVzwD3AG4Gzksz+I+hG4PBJztlTVZNVNTkxMTFIrZKkORYN8CQTSc7qtl8CvAN4jF6QX9sdth3YN6IaJUkL6OcbeTYAU0nW0Qv8vVV1V5JHgVuT/A3wIHDzCOscyDh8W45flCxpvkUDvKq+D1y6wPiT9PrhkqRV4J2YktQoA1ySGmWAS1KjDHBJapQBLkmNMsAlqVEGuCQ1ygCXpEYZ4JLUKANckhplgEtSowxwSWqUAS5JjTLAJalRBrgkNcoAl6RGGeCS1CgDXJIaZYBLUqMMcElqVD/fSq8R8xvnJS2HK3BJapQBLkmNMsAlqVEGuCQ1ygCXpEYZ4JLUqEUDPMn5Se5J8miSHyS5oRs/O8ndSX7Uvb5i9OVKkmb1swI/Dny4qi4BLgc+kOQSYBewv6ouAvZ3+5KkFbJogFfVkar6Xrf9c+Ax4DzgGmCqO2wK2DqiGiVJC1hSDzzJJuBS4AFgfVUd6d46Cqw/yTk7k0wnmZ6ZmRmkVknSHH0HeJKXAV8FPlRVP5v7XlUVUAudV1V7qmqyqiYnJiYGKlaS9Gt9BXiSM+mF9y1V9bVu+OkkG7r3NwDHRlOiJGkh/VyFEuBm4LGq+syct+4Etnfb24F9wy9PknQy/TyN8E3A9cDDSR7qxj4O7Ab2JtkBPAVsG0mFY2ruEwglaTkWDfCq+jcgJ3n7yuGWI0nql3diSlKjDHBJatTYfiPPOPSYx2EOkkbHFbgkNcoAl6RGjW0LpVW2TST1yxW4JDXKAJekRhngktQoA1ySGmWAS1KjDHBJapQBLkmNMsAlqVEGuCQ1ygCXpEYZ4JLUKANckhplgEtSo3waYePmPr3w4O6rln2MpPa4ApekRhngktQoA1ySGmUPvEF+a48kcAUuSc0ywCWpUQa4JDVq0R54ki8CVwPHquo13djZwG3AJuAgsK2qnh1dmf2xNzwcXjcutaGfFfiXgC3zxnYB+6vqImB/ty9JWkGLBnhV3Qf8dN7wNcBUtz0FbB1uWZKkxSz3MsL1VXWk2z4KrD/ZgUl2AjsBLrjggmX+Omlxtn601gz8j5hVVUCd4v09VTVZVZMTExOD/jpJUme5Af50kg0A3eux4ZUkSerHcgP8TmB7t70d2DecciRJ/ernMsKvAG8Fzk1yCLgR2A3sTbIDeArYNsoiNTwn6xN7CabUnkUDvKquO8lbVw65FknSEngnpiQ1yqcRjpGlXkZ3OrZNvBRQ6p8rcElqlAEuSY0ywCWpUfbA1bdT9cyX2q+21y0NzhW4JDXKAJekRtlCGVOnyyWCp0sdy2WrR6czV+CS1CgDXJIaZYBLUqOa74G33mMdR6P+bzLOfekW59ZizePCFbgkNcoAl6RGGeCS1Kjme+AarX772aPoe6/kz7R3qxa5ApekRhngktQoWyhadS1eCjropXP9zNl2jxbjClySGmWAS1KjDHBJapQ9cDVtkF7yyY4Z5NuF5pv7s1ay17/U/vnJ5j+q2+S9/X44XIFLUqMMcElqVKpq+ScnW4DPA+uAm6pq96mOn5ycrOnp6WX/voW0eAmaNCyjaNEM8jOXc24/5yy19bNUg7SQVuJyzyQHqmpy/viyV+BJ1gH/CPwxcAlwXZJLll+iJGkpBmmhXAY8UVVPVtXzwK3ANcMpS5K0mEEC/DzgJ3P2D3VjkqQVsOweeJJrgS1V9f5u/3rgDVX1wXnH7QR2drsXA48vs9ZzgWeWeW6rnPPa4JzH36Dz/d2qmpg/OMh14IeB8+fsb+zG/p+q2gPsGeD3AJBkeqEm/jhzzmuDcx5/o5rvIC2U7wIXJbkwyYuA9wB3DqcsSdJilr0Cr6rjST4I/Cu9ywi/WFU/GFplkqRTGuhW+qr6BvCNIdWymIHbMA1yzmuDcx5/I5nvQDfySJJWj7fSS1KjmgjwJFuSPJ7kiSS7VrueUUjyxSTHkjwyZ+zsJHcn+VH3+orVrHGYkpyf5J4kjyb5QZIbuvFxnvOLk3wnyX90c/7rbvzCJA90n+/buosCxkqSdUkeTHJXtz/Wc05yMMnDSR5KMt2NDf2zfdoH+Bq6Zf9LwJZ5Y7uA/VV1EbC/2x8Xx4EPV9UlwOXAB7r/ruM8518CV1TV64DNwJYklwOfBj5bVa8CngV2rF6JI3MD8Nic/bUw57dV1eY5lw8O/bN92gc4a+SW/aq6D/jpvOFrgKluewrYupI1jVJVHamq73XbP6f3l/s8xnvOVVW/6HbP7P4UcAVwezc+VnMGSLIRuAq4qdsPYz7nkxj6Z7uFAF/Lt+yvr6oj3fZRYP1qFjMqSTYBlwIPMOZz7loJDwHHgLuBHwPPVdXx7pBx/Hx/DvgI8EK3fw7jP+cCvpXkQHc3Oozgs+038jSiqirJ2F0ylORlwFeBD1XVz3qLs55xnHNV/QrYnOQs4A7g1atb0WgluRo4VlUHkrx1lctZSW+uqsNJXgncneSHc98c1me7hRV4X7fsj6mnk2wA6F6PrXI9Q5XkTHrhfUtVfa0bHus5z6qq54B7gDcCZyWZXUyN2+f7TcC7kxyk1/68gt53CIzznKmqw93rMXr/o76MEXy2WwjwtXzL/p3A9m57O7BvFWsZqq4PejPwWFV9Zs5b4zzniW7lTZKXAO+g1/u/B7i2O2ys5lxVH6uqjVW1id7f3W9X1XsZ4zkneWmSl89uA+8EHmEEn+0mbuRJ8i56fbTZW/Y/tboVDV+SrwBvpffUsqeBG4F/AfYCFwBPAduqav4/dDYpyZuB+4GH+XVv9OP0+uDjOufX0vvHq3X0Fk97q+qTSX6P3ur0bOBB4E+q6perV+lodC2Uv6qqq8d5zt3c7uh2zwC+XFWfSnIOQ/5sNxHgkqQTtdBCkSQtwACXpEYZ4JLUKANckhplgEtSowxwrQlJPtE9AfD73RPi3pDklu4pl490T4M8szv2fUle6C77mz3/ke6Wf+m0YYBr7CV5I3A18Pqqei3wdnrP17mF3q3svw+8BHj/nNMOAZ9Y4VKlJfFZKFoLNgDPzN4oUlXPdOP/NXtAku/Qu6V71l3AW5JcXFWPr1il0hK4Atda8C3g/CT/meQLSf5w7ptd6+R64Jtzhl8A/pbe3aHSackA19jrnsH9B8BOYAa4Lcn75hzyBeC+qrp/3qlfBi5PcuGKFCotkS0UrQndY1zvBe5N8jC9hwl9KcmNwATw5wucczzJ3wEfXclapX65AtfYS3JxkovmDG0GnkryfuCPgOuq6oUFT+591d3b6YW8dFpxBa614GXAP3SPcj0OPEGvnXKU3lPh/r37IomvVdUn555YVc8n+Xt6z7CWTis+jVCSGmULRZIaZYBLUqMMcElqlAEuSY0ywCWpUQa4JDXKAJekRhngktSo/wVUK5r5Q4ZQXwAAAABJRU5ErkJggg==\n", + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAX0AAAEGCAYAAACJnEVTAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAASUklEQVR4nO3dbaxlVX3H8e+viA9RIyDjhM5MOzRONZgq2hvAaFpFBXyIwwtrMVZHg5m+oImmNjrqCyJqgk3qU1NtiBBHgyLxoUzUqlPEaJMK3BHkUcpVIcwUmNEB1JhiRv59cdalx+Fe7rnccx/mrO8nuTl7r733OWuFw2+vWXvvdVJVSJL68AerXQFJ0sox9CWpI4a+JHXE0Jekjhj6ktSRx612BR7N8ccfX5s3b17takjSEWXPnj0/r6p1c21b06G/efNmpqenV7saknRESXLnfNsc3pGkjhj6ktQRQ1+SOmLoS1JHRgr9JHckuTHJ9UmmW9lxSXYnub29HtvKk+QTSWaS3JDkBUPvs63tf3uSbcvTJEnSfBbT039pVZ1cVVNtfQdwZVVtAa5s6wCvBLa0v+3Ap2BwkgDOB04FTgHOnz1RSJJWxlKGd7YCO9vyTuDsofLP1sAPgGOSnACcCeyuqoNVdR+wGzhrCZ8vSVqkUUO/gG8n2ZNkeytbX1V3t+V7gPVteQNw19Cxe1vZfOW/J8n2JNNJpg8cODBi9SRJoxj14awXV9W+JM8Adif58fDGqqokY5mYv6ouAi4CmJqacrJ/SRqjkUK/qva11/1JvspgTP7eJCdU1d1t+GZ/230fsGno8I2tbB/wksPKv7uk2q8hm3d8/eHlOy589SrWRJLmt+DwTpInJ3nq7DJwBnATsAuYvQNnG3BFW94FvLndxXMa8EAbBvoWcEaSY9sF3DNamSRphYzS018PfDXJ7P6fr6pvJrkWuDzJucCdwOvb/t8AXgXMAL8B3gpQVQeTfAC4tu13QVUdHFtLJEkLWjD0q+qnwPPmKP8F8LI5ygs4b573ugS4ZPHVlCSNw5qeZXOtGx7Hl6QjgdMwSFJHDH1J6oihL0kdMfQlqSOGviR1xNCXpI4Y+pLUEUNfkjriw1nLwMnXJK1V9vQlqSOGviR1xNCXpI4Y+pLUEUNfkjpi6EtSRwx9SeqIoS9JHfHhrEXy17IkHcns6UtSRwx9SeqIoS9JHTH0Jakjhr4kdcTQl6SOGPqS1BFDX5I6YuhLUkcMfUnqiKEvSR0x9CWpI4a+JHVk5NBPclSS65J8ra2fmOTqJDNJvpjk8a38CW19pm3fPPQe72nltyU5c+ytkSQ9qsX09N8O3Dq0/mHgo1X1TOA+4NxWfi5wXyv/aNuPJCcB5wDPAc4CPpnkqKVVX5K0GCOFfpKNwKuBT7f1AKcDX2q77ATObstb2zpt+8va/luBy6rqwar6GTADnDKGNkiSRjRqT/9jwLuAh9r604H7q+pQW98LbGjLG4C7ANr2B9r+D5fPcYwkaQUsGPpJXgPsr6o9K1AfkmxPMp1k+sCBAyvxkZLUjVF6+i8CXpvkDuAyBsM6HweOSTL7c4sbgX1teR+wCaBtfxrwi+HyOY55WFVdVFVTVTW1bt26RTdIkjS/BUO/qt5TVRurajODC7Hfqao3AlcBr2u7bQOuaMu72jpt+3eqqlr5Oe3unhOBLcA1Y2uJJGlBS/lh9HcDlyX5IHAdcHErvxj4XJIZ4CCDEwVVdXOSy4FbgEPAeVX1uyV8viRpkTLohK9NU1NTNT09vdrV+D2bd3z9MR97x4WvHmNNJGluSfZU1dRc23wiV5I6YuhLUkcMfUnqiKEvSR0x9CWpI4a+JHXE0Jekjhj6ktQRQ1+SOmLoS1JHDH1J6oihL0kdMfQlqSOGviR1ZCnz6XdjKdMpz/c+TrMsaTXY05ekjhj6ktQRQ1+SOmLoS1JHDH1J6oihL0kdMfQlqSOGviR1xNCXpI4Y+pLUEUNfkjpi6EtSRwx9SeqIoS9JHTH0Jakjhr4kdcTQl6SOGPqS1JEFQz/JE5Nck+RHSW5O8v5WfmKSq5PMJPlikse38ie09Zm2ffPQe72nld+W5Mxla5UkaU6j/Ebug8DpVfXrJEcD/5nk34G/Bz5aVZcl+VfgXOBT7fW+qnpmknOADwN/neQk4BzgOcAfAv+R5E+r6nfL0K41z9/LlbQaFuzp18Cv2+rR7a+A04EvtfKdwNlteWtbp21/WZK08suq6sGq+hkwA5wyjkZIkkYz0ph+kqOSXA/sB3YDPwHur6pDbZe9wIa2vAG4C6BtfwB4+nD5HMcMf9b2JNNJpg8cOLDoBkmS5jdS6FfV76rqZGAjg975s5erQlV1UVVNVdXUunXrlutjJKlLi7p7p6ruB64CXggck2T2msBGYF9b3gdsAmjbnwb8Yrh8jmMkSStglLt31iU5pi0/CXgFcCuD8H9d220bcEVb3tXWadu/U1XVys9pd/ecCGwBrhlTOyRJIxjl7p0TgJ1JjmJwkri8qr6W5BbgsiQfBK4DLm77Xwx8LskMcJDBHTtU1c1JLgduAQ4B5/V6544krZYFQ7+qbgCeP0f5T5nj7puq+l/gr+Z5rw8BH1p8NSVJ4+ATuZLUEUNfkjpi6EtSRwx9SeqIoS9JHTH0Jakjhr4kdWSUh7O0zJxmWdJKsacvSR0x9CWpI4a+JHXE0Jekjhj6ktQRQ1+SOmLoS1JHDH1J6oihL0kdMfQlqSNOwzCP4akRJGlS2NOXpI4Y+pLUEYd31hhn3JS0nOzpS1JHDH1J6oihL0kdMfQlqSOGviR1xNCXpI4Y+pLUEUNfkjpi6EtSRwx9SeqIoS9JHVkw9JNsSnJVkluS3Jzk7a38uCS7k9zeXo9t5UnyiSQzSW5I8oKh99rW9r89ybbla9Zk2Lzj6w//SdI4jNLTPwS8s6pOAk4DzktyErADuLKqtgBXtnWAVwJb2t924FMwOEkA5wOnAqcA58+eKCRJK2PB0K+qu6vqh235V8CtwAZgK7Cz7bYTOLstbwU+WwM/AI5JcgJwJrC7qg5W1X3AbuCscTZGkvToFjWmn2Qz8HzgamB9Vd3dNt0DrG/LG4C7hg7b28rmKz/8M7YnmU4yfeDAgcVUT5K0gJFDP8lTgC8D76iqXw5vq6oCahwVqqqLqmqqqqbWrVs3jreUJDUjhX6SoxkE/qVV9ZVWfG8btqG97m/l+4BNQ4dvbGXzlUuSVsgod+8EuBi4tao+MrRpFzB7B8424Iqh8je3u3hOAx5ow0DfAs5Icmy7gHtGK5MkrZBRfi7xRcCbgBuTXN/K3gtcCFye5FzgTuD1bds3gFcBM8BvgLcCVNXBJB8Arm37XVBVB8fRCEnSaDIYjl+bpqamanp6elU+ey3fG+9v50p6NEn2VNXUXNt8IleSOmLoS1JHDH1J6oihL0kdMfQlqSOGviR1xNCXpI4Y+pLUEUNfkjpi6EtSRwx9SerIKBOuaY0ZnhfIeXgkLYY9fUnqiKEvSR0x9CWpI4a+JHXE0Jekjnj3zhHOO3kkLYahP2Qt/0SiJI2DwzuS1BFDX5I64vDOBHF8X9JC7OlLUkcMfUnqiKEvSR1xTH9CHX77qWP8ksCeviR1xdCXpI4Y+pLUEUNfkjpi6EtSR7x7pxM+rSsJ7OlLUlcWDP0klyTZn+SmobLjkuxOcnt7PbaVJ8knkswkuSHJC4aO2db2vz3JtuVpjiTp0YzS0/8McNZhZTuAK6tqC3BlWwd4JbCl/W0HPgWDkwRwPnAqcApw/uyJQitv846vP/wnqS8Lhn5VfQ84eFjxVmBnW94JnD1U/tka+AFwTJITgDOB3VV1sKruA3bzyBOJJGmZPdYx/fVVdXdbvgdY35Y3AHcN7be3lc1X/ghJtieZTjJ94MCBx1g9SdJclnwht6oKqDHUZfb9LqqqqaqaWrdu3bjeVpLEYw/9e9uwDe11fyvfB2wa2m9jK5uvXJK0gh5r6O8CZu/A2QZcMVT+5nYXz2nAA20Y6FvAGUmObRdwz2hlWmVe1JX6suDDWUm+ALwEOD7JXgZ34VwIXJ7kXOBO4PVt928ArwJmgN8AbwWoqoNJPgBc2/a7oKoOvzgsSVpmGQzJr01TU1M1PT29Yp9nb3duPsErHVmS7Kmqqbm2+USuJHXE0Jekjhj6ktQRZ9nUgpyhU5ochr4WxROAdGRzeEeSOtJ9T9/bNCX1xJ6+JHWk+56+HjvH96Ujjz19SeqIPX2Nhb1+6chg6GvsPAFIa5fDO5LUEXv6WjH+C0BafYa+ltV8z0F4ApBWh6GvVecJQFo5hr7WLE8G0vgZ+lpTnBZDWl6Gvo4I9vql8TD0dcSZ718DngykhXUZ+g4hTCb/NSAtrMvQ1+TzBCDNzSdyJakj9vQ18R6t1++/CNQbQ19debTrOaNc6xk+MXjC0JHI0JcWYZQTgycDrWWGvjQGo8wxNMyTgVaLoS+tgsUOJUnjYuhLa9RSrjHMt4/UTej7QJYm0WK/14/lxOA1isnSTehLvVrsxedxvq93O609qarVrsO8pqamanp6eizvZU9fWptGGaIa5dj5LOVkc3h9jpSTVZI9VTU11zZ7+pJW1VI6ZOMa3lqq+U4sa/E6y4r39JOcBXwcOAr4dFVdON++9vQl9WopJ4Y109NPchTwL8ArgL3AtUl2VdUty/F5Br0k/b6VnnDtFGCmqn5aVb8FLgO2rnAdJKlbKz2mvwG4a2h9L3Dq8A5JtgPb2+qvk9y2hM87Hvj5Eo4/0vTWXrDNveiuzfnwktr8x/NtWHMXcqvqIuCicbxXkun5xrUmUW/tBdvcC9s8Pis9vLMP2DS0vrGVSZJWwEqH/rXAliQnJnk8cA6wa4XrIEndWtHhnao6lOTvgG8xuGXzkqq6eRk/cizDREeQ3toLtrkXtnlM1vQTuZKk8fI3ciWpI4a+JHVkIkM/yVlJbksyk2THatdnOSS5JMn+JDcNlR2XZHeS29vrsatZx3FLsinJVUluSXJzkre38oltd5InJrkmyY9am9/fyk9McnX7jn+x3RgxMZIcleS6JF9r65Pe3juS3Jjk+iTTrWxZvtcTF/pDUz28EjgJeEOSk1a3VsviM8BZh5XtAK6sqi3AlW19khwC3llVJwGnAee1/7aT3O4HgdOr6nnAycBZSU4DPgx8tKqeCdwHnLt6VVwWbwduHVqf9PYCvLSqTh66N39ZvtcTF/p0MtVDVX0POHhY8VZgZ1veCZy9knVablV1d1X9sC3/ikEobGCC210Dv26rR7e/Ak4HvtTKJ6rNSTYCrwY+3dbDBLf3USzL93oSQ3+uqR42rFJdVtr6qrq7Ld8DrF/NyiynJJuB5wNXM+HtbkMd1wP7gd3AT4D7q+pQ22XSvuMfA94FPNTWn85ktxcGJ/JvJ9nTpqKBZfper7lpGDQeVVVJJvJ+3CRPAb4MvKOqfjnoCA5MYrur6nfAyUmOAb4KPHt1a7R8krwG2F9Ve5K8ZJWrs5JeXFX7kjwD2J3kx8Mbx/m9nsSefs9TPdyb5ASA9rp/leszdkmOZhD4l1bVV1rxxLcboKruB64CXggck2S20zZJ3/EXAa9NcgeDodnTGfz+xqS2F4Cq2tde9zM4sZ/CMn2vJzH0e57qYRewrS1vA65YxbqMXRvbvRi4tao+MrRpYtudZF3r4ZPkSQx+i+JWBuH/urbbxLS5qt5TVRurajOD/3e/U1VvZELbC5DkyUmeOrsMnAHcxDJ9ryfyidwkr2IwLjg71cOHVrdG45fkC8BLGEw5ey9wPvBvwOXAHwF3Aq+vqsMv9h6xkrwY+D5wI/8/3vteBuP6E9nuJM9lcBHvKAadtMur6oIkf8KgJ3wccB3wN1X14OrVdPza8M4/VNVrJrm9rW1fbauPAz5fVR9K8nSW4Xs9kaEvSZrbJA7vSJLmYehLUkcMfUnqiKEvSR0x9CWpI4a+NIck72uzWt7QZj48NcmlbfbWm9osp0e3fd+S5KF2e+Xs8Te1qSKkNcXQlw6T5IXAa4AXVNVzgZczmM/pUgZTIPwZ8CTgbUOH7QXet8JVlRbNuXekRzoB+Pnswz9V9fNW/j+zOyS5hsF0ALO+BvxFkmdV1W0rVlNpkezpS4/0bWBTkv9O8skkfzm8sQ3rvAn45lDxQ8A/MnhCWFqzDH3pMG3++j8HtgMHgC8mecvQLp8EvldV3z/s0M8DpyU5cUUqKj0GDu9Ic2jTGX8X+G6SGxlMePWZJOcD64C/neOYQ0n+CXj3StZVWgx7+tJhkjwryZahopOBO5O8DTgTeENVPTTnwYOfsXw5gxODtObY05ce6SnAP7cpjQ8BMwyGeu5hMNvhf7UfbvlKVV0wfGBV/TbJJxjMAS+tOc6yKUkdcXhHkjpi6EtSRwx9SeqIoS9JHTH0Jakjhr4kdcTQl6SO/B+suZAgMpEHdAAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] @@ -183,6 +150,7 @@ }, { "cell_type": "markdown", + "id": "a45cceaf", "metadata": {}, "source": [ "## Testing calibration bias with ngmix and autometacal" @@ -190,7 +158,8 @@ }, { "cell_type": "code", - "execution_count": 6, + "execution_count": 8, + "id": "6cd00176", "metadata": {}, "outputs": [], "source": [ @@ -283,14 +252,15 @@ }, { "cell_type": "code", - "execution_count": 7, + "execution_count": null, + "id": "f2b6dbf6", "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ - "100%|██████████| 1000/1000 [00:41<00:00, 24.05it/s]\n" + " 31%|███ | 25297/81499 [16:41<37:00, 25.31it/s] " ] } ], @@ -304,7 +274,7 @@ " wt = np.zeros_like(im) + 1.0/example['noise_std']**2\n", "\n", " psf_obs = ngmix.Observation(\n", - " example['psf'].numpy(),\n", + " example['psf_image'].numpy(),\n", " jacobian=jacobian,\n", " )\n", " \n", @@ -324,7 +294,8 @@ }, { "cell_type": "code", - "execution_count": 8, + "execution_count": null, + "id": "6fd72de1", "metadata": {}, "outputs": [], "source": [ @@ -347,21 +318,10 @@ }, { "cell_type": "code", - "execution_count": 9, + "execution_count": null, + "id": "a9369099", "metadata": {}, - "outputs": [ - { - "name": "stdout", - "output_type": "stream", - "text": [ - "ngmix results:\n", - "S/N: 17.1766\n", - "R11: 0.326126\n", - "m: -1.41607 +/- 7.39404 (99.7% conf)\n", - "c: 0.064798 +/- 0.135984 (99.7% conf)\n" - ] - } - ], + "outputs": [], "source": [ "print('ngmix results:')\n", "print('S/N: %g' % s2n)\n", @@ -372,7 +332,8 @@ }, { "cell_type": "code", - "execution_count": 46, + "execution_count": null, + "id": "9f6d2555", "metadata": {}, "outputs": [], "source": [ @@ -382,6 +343,7 @@ }, { "cell_type": "markdown", + "id": "66912f6b", "metadata": {}, "source": [ "#### Now doing the same thing with autometacal" @@ -389,7 +351,8 @@ }, { "cell_type": "code", - "execution_count": 47, + "execution_count": null, + "id": "d3c7ae30", "metadata": {}, "outputs": [], "source": [ @@ -403,7 +366,8 @@ }, { "cell_type": "code", - "execution_count": 59, + "execution_count": null, + "id": "6eb6d850", "metadata": {}, "outputs": [], "source": [ @@ -416,19 +380,12 @@ }, { "cell_type": "code", - "execution_count": 60, + "execution_count": null, + "id": "3c4416d7", "metadata": { "scrolled": true }, - "outputs": [ - { - "name": "stderr", - "output_type": "stream", - "text": [ - "100%|██████████| 500/500 [01:05<00:00, 7.66it/s]\n" - ] - } - ], + "outputs": [], "source": [ "res_e = []\n", "res_R = []\n", @@ -442,27 +399,18 @@ }, { "cell_type": "code", - "execution_count": 61, + "execution_count": null, + "id": "4853c86d", "metadata": {}, - "outputs": [ - { - "data": { - "text/plain": [ - "(1538.4615384615386, 'measurements per sec')" - ] - }, - "execution_count": 61, - "metadata": {}, - "output_type": "execute_result" - } - ], + "outputs": [], "source": [ "100000/65, 'measurements per sec' # Batch of 200" ] }, { "cell_type": "code", - "execution_count": 62, + "execution_count": null, + "id": "eb9db1dc", "metadata": {}, "outputs": [], "source": [ @@ -473,7 +421,8 @@ }, { "cell_type": "code", - "execution_count": 63, + "execution_count": null, + "id": "fb782b7d", "metadata": {}, "outputs": [], "source": [ @@ -486,20 +435,10 @@ }, { "cell_type": "code", - "execution_count": 64, + "execution_count": null, + "id": "38f28ee0", "metadata": {}, - "outputs": [ - { - "name": "stdout", - "output_type": "stream", - "text": [ - "autometacal results:\n", - "R11: 0.188488\n", - "m: -0.829483 +/- 0.373999 (99.7% conf)\n", - "c: 0.00575358 +/- 0.00748643 (99.7% conf)\n" - ] - } - ], + "outputs": [], "source": [ "print('autometacal results:')\n", "print('R11: %g' % auto_R[0,0])\n", @@ -510,6 +449,7 @@ { "cell_type": "code", "execution_count": null, + "id": "ebbea8c3", "metadata": {}, "outputs": [], "source": [] @@ -517,6 +457,7 @@ { "cell_type": "code", "execution_count": null, + "id": "87bea4ac", "metadata": {}, "outputs": [], "source": [] @@ -524,6 +465,7 @@ { "cell_type": "code", "execution_count": null, + "id": "327662da", "metadata": {}, "outputs": [], "source": [] @@ -531,6 +473,7 @@ { "cell_type": "code", "execution_count": null, + "id": "eebd68ef", "metadata": {}, "outputs": [], "source": [] @@ -552,7 +495,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.8.6" + "version": "3.8.12" } }, "nbformat": 4, From 5a56f8074e94d809174d2a2f0fe0d798557b7dbc Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Andr=C3=A9=20Z=2E=20Vitorelli?= Date: Thu, 11 Aug 2022 18:38:16 +0200 Subject: [PATCH 4/6] reestructuring datasets for practical purposes --- notebooks/CFIS_dataset.ipynb | 107 +++++++++++++++++++++++++++-------- 1 file changed, 82 insertions(+), 25 deletions(-) diff --git a/notebooks/CFIS_dataset.ipynb b/notebooks/CFIS_dataset.ipynb index ad14be5..94bb479 100644 --- a/notebooks/CFIS_dataset.ipynb +++ b/notebooks/CFIS_dataset.ipynb @@ -41,7 +41,7 @@ }, { "cell_type": "code", - "execution_count": 3, + "execution_count": 42, "id": "a085afdf", "metadata": {}, "outputs": [], @@ -53,7 +53,7 @@ }, { "cell_type": "code", - "execution_count": 4, + "execution_count": 43, "id": "c5f6270c-07f1-4fc3-8d92-8b50b191583e", "metadata": {}, "outputs": [], @@ -252,7 +252,7 @@ }, { "cell_type": "code", - "execution_count": null, + "execution_count": 10, "id": "f2b6dbf6", "metadata": {}, "outputs": [ @@ -260,7 +260,7 @@ "name": "stderr", "output_type": "stream", "text": [ - " 31%|███ | 25297/81499 [16:41<37:00, 25.31it/s] " + "100%|██████████| 81499/81499 [53:27<00:00, 25.41it/s] \n" ] } ], @@ -294,7 +294,7 @@ }, { "cell_type": "code", - "execution_count": null, + "execution_count": 11, "id": "6fd72de1", "metadata": {}, "outputs": [], @@ -318,10 +318,22 @@ }, { "cell_type": "code", - "execution_count": null, + "execution_count": 12, "id": "a9369099", "metadata": {}, - "outputs": [], + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "ngmix results:\n", + "S/N: 16.4136\n", + "R11: 0.299841\n", + "m: 0.439048 +/- 0.887137 (99.7% conf)\n", + "c: 0.129178 +/- 0.0176289 (99.7% conf)\n" + ] + } + ], "source": [ "print('ngmix results:')\n", "print('S/N: %g' % s2n)\n", @@ -332,13 +344,13 @@ }, { "cell_type": "code", - "execution_count": null, + "execution_count": 44, "id": "9f6d2555", "metadata": {}, "outputs": [], "source": [ "# saving the reconv psf\n", - "reconv_psf = tf.convert_to_tensor(obsdict['noshear'].psf.image.reshape([1,51,51]).repeat(200,axis=0).astype('float32'))" + "reconv_psf = tf.convert_to_tensor(obsdict['noshear'].psf.image.reshape([1,51,51]).repeat(1000,axis=0).astype('float32'))" ] }, { @@ -351,7 +363,7 @@ }, { "cell_type": "code", - "execution_count": null, + "execution_count": 45, "id": "d3c7ae30", "metadata": {}, "outputs": [], @@ -361,12 +373,12 @@ "@tf.function\n", "def get_autometacal_shape(im, example):\n", " method = lambda x: autometacal.get_moment_ellipticities(x, scale=0.187, fwhm=weight_fwhm)\n", - " return autometacal.get_metacal_response(im, example['psf'], reconv_psf, method)" + " return autometacal.get_metacal_response(im, example['psf_image'], reconv_psf, method)" ] }, { "cell_type": "code", - "execution_count": null, + "execution_count": 46, "id": "6eb6d850", "metadata": {}, "outputs": [], @@ -374,18 +386,73 @@ "dset = tfds.load('CFIS/parametric_shear_1k', split='train')\n", "dset = dset.map(add_noise)\n", "dset = dset.repeat(100)\n", - "dset = dset.batch(200)\n", + "dset = dset.batch(1000)\n", "dset = dset.prefetch(buffer_size=5)" ] }, { "cell_type": "code", - "execution_count": null, + "execution_count": 47, "id": "3c4416d7", "metadata": { "scrolled": true }, - "outputs": [], + "outputs": [ + { + "name": "stderr", + "output_type": "stream", + "text": [ + " 0%| | 0/8150 [00:00ResamplerGrad\n" + ] + }, + { + "name": "stderr", + "output_type": "stream", + "text": [ + "WARNING:tensorflow:Using a while_loop for converting Addons>ResamplerGrad\n" + ] + }, + { + "name": "stdout", + "output_type": "stream", + "text": [ + "WARNING:tensorflow:Using a while_loop for converting Addons>ResamplerGrad\n" + ] + }, + { + "name": "stderr", + "output_type": "stream", + "text": [ + "WARNING:tensorflow:Using a while_loop for converting Addons>ResamplerGrad\n", + " 0%| | 34/8150 [03:12<12:45:43, 5.66s/it]\n" + ] + }, + { + "ename": "KeyboardInterrupt", + "evalue": "", + "output_type": "error", + "traceback": [ + "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", + "\u001b[0;31mKeyboardInterrupt\u001b[0m Traceback (most recent call last)", + "\u001b[0;32m\u001b[0m in \u001b[0;36m\u001b[0;34m\u001b[0m\n\u001b[1;32m 2\u001b[0m \u001b[0mres_R\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0;34m[\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 3\u001b[0m \u001b[0;32mfor\u001b[0m \u001b[0mim\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mexample\u001b[0m \u001b[0;32min\u001b[0m \u001b[0mtqdm\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mdset\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 4\u001b[0;31m \u001b[0me\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mR\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mget_autometacal_shape\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mim\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mexample\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 5\u001b[0m \u001b[0mres_e\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mappend\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0me\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 6\u001b[0m \u001b[0mres_R\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mappend\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mR\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", + "\u001b[0;32m~/.local/lib/python3.8/site-packages/tensorflow/python/util/traceback_utils.py\u001b[0m in \u001b[0;36merror_handler\u001b[0;34m(*args, **kwargs)\u001b[0m\n\u001b[1;32m 148\u001b[0m \u001b[0mfiltered_tb\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0;32mNone\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 149\u001b[0m \u001b[0;32mtry\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 150\u001b[0;31m \u001b[0;32mreturn\u001b[0m \u001b[0mfn\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m**\u001b[0m\u001b[0mkwargs\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 151\u001b[0m \u001b[0;32mexcept\u001b[0m \u001b[0mException\u001b[0m \u001b[0;32mas\u001b[0m \u001b[0me\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 152\u001b[0m \u001b[0mfiltered_tb\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0m_process_traceback_frames\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0me\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m__traceback__\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", + "\u001b[0;32m~/.local/lib/python3.8/site-packages/tensorflow/python/eager/def_function.py\u001b[0m in \u001b[0;36m__call__\u001b[0;34m(self, *args, **kwds)\u001b[0m\n\u001b[1;32m 908\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 909\u001b[0m \u001b[0;32mwith\u001b[0m \u001b[0mOptionalXlaContext\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_jit_compile\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 910\u001b[0;31m \u001b[0mresult\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_call\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m**\u001b[0m\u001b[0mkwds\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 911\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 912\u001b[0m \u001b[0mnew_tracing_count\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mexperimental_get_tracing_count\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", + "\u001b[0;32m~/.local/lib/python3.8/site-packages/tensorflow/python/eager/def_function.py\u001b[0m in \u001b[0;36m_call\u001b[0;34m(self, *args, **kwds)\u001b[0m\n\u001b[1;32m 947\u001b[0m \u001b[0;31m# In this case we have not created variables on the first call. So we can\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 948\u001b[0m \u001b[0;31m# run the first trace but we should fail if variables are created.\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 949\u001b[0;31m \u001b[0mresults\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_stateful_fn\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m**\u001b[0m\u001b[0mkwds\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 950\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_created_variables\u001b[0m \u001b[0;32mand\u001b[0m \u001b[0;32mnot\u001b[0m \u001b[0mALLOW_DYNAMIC_VARIABLE_CREATION\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 951\u001b[0m raise ValueError(\"Creating variables on a non-first call to a function\"\n", + "\u001b[0;32m~/.local/lib/python3.8/site-packages/tensorflow/python/eager/function.py\u001b[0m in \u001b[0;36m__call__\u001b[0;34m(self, *args, **kwargs)\u001b[0m\n\u001b[1;32m 3128\u001b[0m (graph_function,\n\u001b[1;32m 3129\u001b[0m filtered_flat_args) = self._maybe_define_function(args, kwargs)\n\u001b[0;32m-> 3130\u001b[0;31m return graph_function._call_flat(\n\u001b[0m\u001b[1;32m 3131\u001b[0m filtered_flat_args, captured_inputs=graph_function.captured_inputs) # pylint: disable=protected-access\n\u001b[1;32m 3132\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n", + "\u001b[0;32m~/.local/lib/python3.8/site-packages/tensorflow/python/eager/function.py\u001b[0m in \u001b[0;36m_call_flat\u001b[0;34m(self, args, captured_inputs, cancellation_manager)\u001b[0m\n\u001b[1;32m 1957\u001b[0m and executing_eagerly):\n\u001b[1;32m 1958\u001b[0m \u001b[0;31m# No tape is watching; skip to running the function.\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m-> 1959\u001b[0;31m return self._build_call_outputs(self._inference_function.call(\n\u001b[0m\u001b[1;32m 1960\u001b[0m ctx, args, cancellation_manager=cancellation_manager))\n\u001b[1;32m 1961\u001b[0m forward_backward = self._select_forward_and_backward_functions(\n", + "\u001b[0;32m~/.local/lib/python3.8/site-packages/tensorflow/python/eager/function.py\u001b[0m in \u001b[0;36mcall\u001b[0;34m(self, ctx, args, cancellation_manager)\u001b[0m\n\u001b[1;32m 596\u001b[0m \u001b[0;32mwith\u001b[0m \u001b[0m_InterpolateFunctionError\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 597\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0mcancellation_manager\u001b[0m \u001b[0;32mis\u001b[0m \u001b[0;32mNone\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 598\u001b[0;31m outputs = execute.execute(\n\u001b[0m\u001b[1;32m 599\u001b[0m \u001b[0mstr\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0msignature\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mname\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 600\u001b[0m \u001b[0mnum_outputs\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_num_outputs\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", + "\u001b[0;32m~/.local/lib/python3.8/site-packages/tensorflow/python/eager/execute.py\u001b[0m in \u001b[0;36mquick_execute\u001b[0;34m(op_name, num_outputs, inputs, attrs, ctx, name)\u001b[0m\n\u001b[1;32m 56\u001b[0m \u001b[0;32mtry\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 57\u001b[0m \u001b[0mctx\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mensure_initialized\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m---> 58\u001b[0;31m tensors = pywrap_tfe.TFE_Py_Execute(ctx._handle, device_name, op_name,\n\u001b[0m\u001b[1;32m 59\u001b[0m inputs, attrs, num_outputs)\n\u001b[1;32m 60\u001b[0m \u001b[0;32mexcept\u001b[0m \u001b[0mcore\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_NotOkStatusException\u001b[0m \u001b[0;32mas\u001b[0m \u001b[0me\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", + "\u001b[0;31mKeyboardInterrupt\u001b[0m: " + ] + } + ], "source": [ "res_e = []\n", "res_R = []\n", @@ -397,16 +464,6 @@ "res_R = np.concatenate(res_R)" ] }, - { - "cell_type": "code", - "execution_count": null, - "id": "4853c86d", - "metadata": {}, - "outputs": [], - "source": [ - "100000/65, 'measurements per sec' # Batch of 200" - ] - }, { "cell_type": "code", "execution_count": null, From bfac86a5e654d680839b710a4ce53aae8147d242 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Andr=C3=A9=20Z=2E=20Vitorelli?= Date: Thu, 11 Aug 2022 19:03:42 +0200 Subject: [PATCH 5/6] datasets testing --- autometacal/python/datasets/simple.py | 15 +++++++-------- 1 file changed, 7 insertions(+), 8 deletions(-) diff --git a/autometacal/python/datasets/simple.py b/autometacal/python/datasets/simple.py index 65f55eb..ff53aa1 100644 --- a/autometacal/python/datasets/simple.py +++ b/autometacal/python/datasets/simple.py @@ -51,12 +51,11 @@ def _info(self): description=_DESCRIPTION, homepage=_URL, features=tfds.features.FeaturesDict({ - 'label': tfds.features.Tensor(shape=[2], dtype=tf.float32), 'gal_model': tfds.features.Tensor( shape=[self.builder_config.stamp_size,self.builder_config.stamp_size], dtype=tf.float32 ), - 'obs': tfds.features.Tensor( + 'obs_image': tfds.features.Tensor( shape=[self.builder_config.stamp_size,self.builder_config.stamp_size], dtype=tf.float32 ), @@ -64,12 +63,12 @@ def _info(self): shape=[self.builder_config.stamp_size, self.builder_config.stamp_size], dtype=tf.float32 ), - 'noise': tfds.features.Tensor( + 'noise_image': tfds.features.Tensor( shape=[self.builder_config.stamp_size, self.builder_config.stamp_size], dtype=tf.float32 ), }), - supervised_keys=("obs"), + supervised_keys=("obs","obs"), citation=_CITATION ) @@ -80,7 +79,7 @@ def _split_generators(self,dl): self.builder_config.stamp_size )} - def _generate_examples(self, dataset_size, stamp_size,g1,g2): + def _generate_examples(self, dataset_size, stamp_size): """Yields examples.""" rng = np.random.RandomState(31415) @@ -89,10 +88,10 @@ def _generate_examples(self, dataset_size, stamp_size,g1,g2): #generate example model, obs_img, psf_img, noise_img = make_data(rng,stamp_size=stamp_size) yield '%d'%i, {'gal_model': model, #noiseless PSFless galaxy model - 'obs': obs_img, #observed image + 'obs_image': obs_img, #observed image 'psf_image': psf_img, #psf image - 'noise' : noise_img, - 'label': np.array([g1,g2],dtype='float32')} + 'noise_image' : noise_img, + } def make_data(rng,stamp_size = 51): """Simple exponetial profile toy model galaxy""" From e78cbbeadbbbaf9668e1a7dde6d02b8e61d83e92 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Andr=C3=A9=20Z=2E=20Vitorelli?= Date: Thu, 11 Aug 2022 19:19:00 +0200 Subject: [PATCH 6/6] changing filename --- autometacal/python/datasets/cfis.py | 137 ++++++++++++++++++++++++++++ 1 file changed, 137 insertions(+) create mode 100644 autometacal/python/datasets/cfis.py diff --git a/autometacal/python/datasets/cfis.py b/autometacal/python/datasets/cfis.py new file mode 100644 index 0000000..ea11fd3 --- /dev/null +++ b/autometacal/python/datasets/cfis.py @@ -0,0 +1,137 @@ +""" TensorFlow Dataset of simulated CFIS images.""" +import tensorflow_datasets as tfds +import tensorflow as tf +import numpy as np +import galsim as gs + +_CITATION = """{NEEDED}""" +_URL = "https://github.com/CosmoStat/autometacal" +_DESCRIPTION = """Noiseless CFHT-pixel galaxies with shapes from COSMOS""" + +class CFISConfig(tfds.core.BuilderConfig): + """BuilderConfig for CFIS Galaxies.""" + + def __init__(self, *, galaxy_type="parametric", + data_set_size=1000, stamp_size=51, + shear_g1=0.0, shear_g2=0.0, **kwargs): + """BuilderConfig for CFIS. + Args: + size: size of sample. + stamp_size: image stamp size in pixels. + pixel_scale: pixel scale of stamps in arcsec. + **kwargs: keyword arguments forwarded to super. + """ + v1 = tfds.core.Version("0.0.1") + super(CFISConfig, self).__init__( + description=("Galaxy stamps"), + version=v1, + **kwargs) + + # Adjustable parameters + self.galaxy_type = galaxy_type + self.data_set_size = data_set_size + self.stamp_size = stamp_size + self.kstamp_size = stamp_size + self.shear_g1 = shear_g1 + self.shear_g2 = shear_g2 + + # Fixed parameters + self.pixel_scale = 0.187 + self.psf_fwhm = 0.65 # TODO: add varying PSF + self.psf_e1 = 0.0 + self.psf_e2 = 0.025 + + +class CFIS(tfds.core.GeneratorBasedBuilder): + """DatasetBuilder for simulated CFIS dataset.""" + + VERSION = tfds.core.Version('0.0.1') + RELEASE_NOTES = { + '0.0.1': 'pre alpha release.', + } + + BUILDER_CONFIGS = [ + CFISConfig( + name="parametric_1k", + galaxy_type="parametric", + data_set_size=81499), + CFISConfig( + name="parametric_shear_1k", + galaxy_type="parametric", + data_set_size=81499, + shear_g1=0.02) + ] + + def _info(self) -> tfds.core.DatasetInfo: + """Returns the dataset metadata.""" + # TODO: Specifies the tfds.core.DatasetInfo object + return tfds.core.DatasetInfo( + builder=self, + description=_DESCRIPTION, + features=tfds.features.FeaturesDict({ + 'gal_image': tfds.features.Tensor( + shape=[self.builder_config.stamp_size, self.builder_config.stamp_size], + dtype=tf.float32 + ), + 'psf_image': tfds.features.Tensor( + shape=[self.builder_config.stamp_size, self.builder_config.stamp_size], + dtype=tf.float32), + "noise_std": tfds.features.Tensor(shape=[1], dtype=tf.float32), + "mag": tfds.features.Tensor(shape=[1], dtype=tf.float32), + }), + supervised_keys=("obs", "obs"), + homepage=_URL, + citation=_CITATION, + ) + + def _split_generators(self, dl_manager: tfds.download.DownloadManager): + """Returns SplitGenerators.""" + if self.builder_config.data_set_size: + size = self.builder_config.data_set_size + else: + size = 81499 + + return [ + tfds.core.SplitGenerator( + name=tfds.Split.TRAIN, + gen_kwargs={ + "size": size, + },), + ] + + def _generate_examples(self, size): + """Yields examples.""" + # Loads the galsim COSMOS catalog + cat = gs.COSMOSCatalog(sample="25.2") + psf = gs.Kolmogorov(fwhm=self.builder_config.psf_fwhm, flux=1.0) + + for i in range(size): + # retrieving galaxy and magnitude + gal = cat.makeGalaxy(i, gal_type='parametric') + gal_mag = cat.param_cat['mag_auto'][cat.orig_index[i]] + sky_level = 400 + mag_zp = 32. + gal_flux = 10**(-(gal_mag-mag_zp)/2.5) + + gal = gal.withFlux(gal_flux) + gal = gal.shear(g1=self.builder_config.shear_g1, g2=self.builder_config.shear_g2) + + gal_conv = gs.Convolve(gal, psf) + method="auto" + gal_stamp = gal_conv.drawImage(nx=self.builder_config.stamp_size, + ny=self.builder_config.stamp_size, + scale=self.builder_config.pixel_scale, + method=method + ).array.astype('float32') + + psf_stamp = psf.drawImage(nx=self.builder_config.stamp_size, + ny=self.builder_config.stamp_size, + scale=self.builder_config.pixel_scale, + method=method + ).array.astype('float32') + + + yield '%d'%i, {"gal_image": gal_stamp, + "psf_image": psf_stamp, + "noise_std": np.array([np.sqrt(sky_level)]).astype('float32'), + "mag": np.array([gal_mag]).astype('float32')} \ No newline at end of file