pynkowski.data.bak_data
1import numpy as np 2import healpy as hp 3 4try: 5 from tqdm.notebook import tqdm 6except: 7 tqdm = lambda x: x 8 print('tqdm not loaded') 9 10 11 12def get_theta(nside): 13 """Define a HEALPix map with the value of θ in each pixel at the input nside""" 14 theta, _ = hp.pix2ang(nside, np.arange(12 * nside ** 2)) 15 return np.pi/2. - theta 16 17def derivatives(mapp, lmax=None, gradient=False, **kwargs): 18 """Find the derivatives d_theta, d_phi of a Healpix map. It uses the healpy alm2map_der1 function. 19 20 Parameters 21 ---------- 22 mapp : np.array 23 The Healpix map to find its derivatives. 24 25 lmax : int, optional 26 Maximum multipole to get the alm and the derivatives. It can create numerical errors if it is too high. 27 Default: 3*nside-1. 28 29 gradient : bool, optional 30 If True, return the covariant derivatives. If False, return the partial derivatives. Default: False. 31 32 **kwargs : 33 Extra keywords to pass to the map2alm function. 34 35 Returns 36 ------- 37 d_theta : np.array 38 A healpix map with the derivative with respect to theta. 39 40 d_phi : np.array 41 A healpix map with the derivatives with respect to phi, without normalizing. 42 """ 43 nside = hp.get_nside(mapp) 44 if lmax is None: 45 lmax = 3 * nside - 1 46 47 alm = hp.map2alm(mapp, lmax=lmax, **kwargs) 48 [_, d_theta, d_phiosin] = hp.alm2map_der1(alm, nside, lmax=lmax) 49 d_theta = -d_theta 50 51 if gradient: 52 return (d_theta, d_phiosin) 53 54 theta = get_theta(nside) 55 56 d_phi = d_phiosin * np.cos(theta) 57 return (d_theta, d_phi) 58 59 60def second_derivatives(d_theta, d_phi, lmax=None, **kwargs): 61 """Find the Second derivatives for every pixel of a Healpix map given the first partial derivatives. 62 63 Parameters 64 ---------- 65 d_theta : np.array 66 The partial theta derivative of the Healpix map. 67 68 d_phi : np.array 69 The partial phi derivative of the Healpix map. 70 71 lmax : int, optional 72 Maximum multipole to get the alm and the derivatives. It can create numerical errors if it is too high. 73 Default: 3*nside-1. 74 75 Returns 76 ------- 77 d_thetatheta : np.array 78 A Healpix map of the second partial derivative wrt theta. 79 80 d_phiphi : np.array 81 A Healpix map of the second partial derivative wrt phi. 82 83 d_phitheta : np.array 84 A Healpix map of the second partial derivative wrt theta and phi. 85 86 """ 87 nside = hp.get_nside(d_theta) 88 89 theta = get_theta(nside) 90 91 if lmax is None: 92 lmax = 3 * nside - 1 93 94 d_phitheta, d_phiphi = derivatives(d_phi, lmax=lmax, **kwargs) 95 96 alm_theta = hp.map2alm(d_theta, lmax=lmax, **kwargs) 97 [_, d_thetatheta, _] = hp.alm2map_der1(alm_theta, nside, lmax=lmax) 98 d_thetatheta = -d_thetatheta 99 100 101 return (d_thetatheta, d_phiphi, d_phitheta) 102 103 104class Scalar(): 105 """Class to compute Minkowski functionals (MFs) and extrema of Healpix scalar maps. It computes and stores spatial first and 106 second derivatives of Healpix scalar maps. 107 108 Parameters 109 ---------- 110 Smap : np.array 111 The input Healpix map where all the statistics will be computed. 112 113 normalise : bool, optional 114 If True, divide input Smap by its standard deviation. Default: True. 115 116 mask : np.array, optional 117 An input Healpix mask. All the statiscal quantities will be computed only within this mask. Default: None (all sky is considered). 118 119 Attributes 120 ---------- 121 Smap : np.array 122 The Healpix map where all the statistics are computed. 123 124 nside : int 125 The Nside parameter of Smap 126 127 mask : np.array 128 A Healpix mask whithin which all the statiscal quantities are computed. 129 130 grad_phi : np.array 131 The phi covariant first derivative of Smap (only if first or second MFs are computed). 132 133 grad_theta : np.array 134 The theta covariant first derivative of Smap (only if first or second MFs are computed). 135 136 der_phi : np.array 137 The phi partial first derivative of Smap (only if first or second MFs are computed). 138 139 der_theta_phi : np.array 140 The second covariant derivative wrt theta and phi of Smap (only if second MF is computed). 141 142 der_phi_phi : np.array 143 The second covariant derivative wrt phi of Smap (only if second MF is computed). 144 145 der_theta_theta : np.array 146 The second covariant derivative wrt theta of Smap (only if second MF is computed). 147 148 Notes 149 ----- 150 The derivatives are always computed full-sky regardless of input mask. 151 152 """ 153 def __init__(self, Smap, normalise=True, mask=None): 154 """Initialise the class to compute Minkowski functionals (MFs) and extrema of Healpix scalar maps. 155 156 Parameters 157 ---------- 158 Smap : np.array 159 The input Healpix map where all the statistics will be computed. 160 161 normalise : bool, optional 162 If True, divide input Smap by its standard deviation. Default: True. 163 164 mask : np.array, optional 165 An input Healpix mask. All the statiscal quantities will be computed only within this mask. Default: None (all sky is considered). 166 167 168 """ 169 self.Smap = Smap.copy() 170 self.nside = hp.get_nside(Smap) 171 if mask is None: 172 self.mask = np.ones(Smap.shape, dtype='bool') 173 else: 174 self.mask = mask 175 if hp.get_nside(mask) != self.nside: 176 raise ValueError('The map and the mask have different nside') 177 if normalise: 178 σ = self.get_variance() 179 self.Smap /= np.sqrt(σ) 180 self.grad_phi = None 181 self.der_phi_phi = None 182 183 184 def __repr__(self): 185 return(f'map = {self.Smap}') 186 187 def get_variance(self): 188 """compute the variance of the input Healpix scalar map within the input mask. 189 190 Returns 191 ------- 192 var : float 193 The variance of the input Healpix map within the input mask. 194 195 """ 196 return (np.var(self.Smap[self.mask])) 197 198 def set_pix(self, pixs): 199 """return the values of the input Healpix scalar map in pixels pixs. 200 201 Parameters 202 ---------- 203 pixs : np.array 204 The indices of the input map pixels whose values of the map are returned. 205 206 Returns 207 ------- 208 values : np.array 209 The values of the input Healpix scalar map in pixels pixs. 210 211 """ 212 return self.Smap[pixs] 213 214 def get_gradient(self): 215 """Compute the covariant and partial first derivatives of the input Healpix scalar map. 216 It stores: 217 218 - first covariant derivative wrt theta in self.grad_theta 219 - first partial derivative wrt phi in self.der_phi 220 - first covariant derivative wrt phi in self.grad_phi 221 222 """ 223 S_grad = derivatives(self.Smap, gradient=True) 224 theta = get_theta(self.nside) 225 self.grad_theta = Scalar(S_grad[0], normalise=False) 226 self.der_phi = Scalar(np.cos(theta) * S_grad[1], normalise=False) 227 self.grad_phi = Scalar(S_grad[1], normalise=False) 228 229 230 def get_hessian(self): 231 """compute the covariant second derivatives of the input Healpix scalar map. 232 It stores: 233 234 - second covariant derivative wrt theta in self.der_theta_theta 235 - second covariant derivative wrt phi in self.der_phi_phi 236 - second covariant derivative wrt theta and phi in self.der_theta_phi 237 238 """ 239 if self.grad_phi == None: 240 self.get_gradient() 241 theta = get_theta(self.nside) 242 243 S_der_der = second_derivatives(self.grad_theta.Smap, self.der_phi.Smap) 244 245 self.der_theta_theta = Scalar(S_der_der[0], normalise=False) 246 self.der_phi_phi = Scalar(S_der_der[1]/np.cos(theta)**2. + self.grad_theta.Smap*np.tan(theta), normalise=False) 247 self.der_theta_phi = Scalar((S_der_der[2]/np.cos(theta) - self.grad_phi.Smap * np.tan(theta)) , normalise=False) 248 249 250 def get_κ(self, pixs): 251 """Compute the geodesic curvature multiplied by the modulus of the gradient in pixels pixs. If not already computed, it computes 252 the first and second covariant derivatives of the input map. 253 254 Parameters 255 ---------- 256 pixs : np.array 257 The indices of the input map pixels where geodesic curvature is computed. 258 259 Returns 260 ------- 261 k : np.array 262 The geodesic curvature in pixels `pixs`. 263 264 """ 265 num = 2.*self.grad_theta.set_pix(pixs)*self.grad_phi.set_pix(pixs)*self.der_theta_phi.set_pix(pixs) - self.grad_phi.set_pix(pixs)**2. * self.der_theta_theta.set_pix(pixs) - self.grad_theta.set_pix(pixs)**2. * self.der_phi_phi.set_pix(pixs) 266 den = self.grad_theta.set_pix(pixs)**2. + self.grad_phi.set_pix(pixs)**2. 267 return num / den 268 269 270 def V0_pixel(self, u): 271 """Determine where input Healpix scalar map Smap is greater than threshold u. 272 273 Parameters 274 ---------- 275 u : float 276 The threshold considered for the computation of first Minkowski functional V0. 277 278 Returns 279 ------- 280 v0map : np.array 281 a bool array with the same shape as the input map, with False where input map values are lower than threshold u. 282 283 """ 284 return self.Smap>u 285 286 def V0_iter(self, u): 287 """Compute the normalised first Minkowski functional v0 at the threshold u within the given mask. 288 289 Parameters 290 ---------- 291 u : float 292 The threshold considered for the computation of v0. 293 294 Returns 295 ------- 296 v0 : np.array 297 First normalised Minkowski functional evaluated at threshold u within the given mask. 298 299 """ 300 return np.mean(self.V0_pixel(u)[self.mask]) 301 302 def V0(self, us): 303 """Compute the normalised first Minkowski functional v0 at the different thresholds us within the given mask. 304 305 Parameters 306 ---------- 307 us : np.array 308 The thresholds considered for the computation of v0. 309 310 Returns 311 ------- 312 v0s : np.array 313 First normalised Minkowski functional evaluated at thresholds us within the given mask. 314 315 """ 316 us = np.atleast_1d(us) 317 return np.array([self.V0_iter(u) for u in tqdm(us)]) 318 319 320 def V1_pixel(self, u, du): 321 """Compute the modulus of the gradient where the values of the input map are between `u-du/2` and `u+du/2`. 322 323 Parameters 324 ---------- 325 u : float 326 The centered value of the bin considered for the computation. 327 328 du : float 329 The width of the bin considered for the computation. 330 331 Returns 332 ------- 333 v1map : np.array 334 Modulus of the gradient where `u-du/2 < Smap < u+du/2`. 335 336 """ 337 theta = get_theta(self.nside) 338 339 areas = np.zeros_like(theta) 340 341 mask = (u + du/2. > self.Smap) & (u - du/2. <= self.Smap) & ~(np.isclose(np.cos(theta),0, atol=1.e-2)) 342 pixs = np.arange(12*self.nside**2)[mask] 343 344 areas[pixs] = np.sqrt(self.grad_theta.set_pix(pixs)**2. + self.grad_phi.set_pix(pixs)**2.) 345 return areas/du / 4. 346 347 def V1_iter(self, u, du): 348 """Compute the normalised second Minkowski functional v1 in the bin `u-du/2` and `u+du/2` within the given mask. 349 350 Parameters 351 ---------- 352 u : float 353 The centered value of the bin considered for the computation of v1. 354 355 du : float 356 The width of the bin considered for the computation of v1. 357 358 Returns 359 ------- 360 v1 : float 361 Second normalised Minkowski functional evaluated in the bin u-du/2 and u+du/2 within the given mask. 362 363 """ 364 return np.mean(self.V1_pixel(u, du)[self.mask]) 365 366 def V1(self, us, edges=False): 367 """Compute the normalised second Minkowski functional v1 at the different thresholds us within the given mask. 368 369 Parameters 370 ---------- 371 us : np.array 372 The thresholds considered for the computation of v2. See 'edges' for details. 373 374 edges : bool, optional 375 If False, us is considered as an array of uniformly distributed thresholds. 376 If True, us is considered as a monotonically increasing array of bin edges, including the rightmost edge, allowing for non-uniform distributions of thresholds. 377 In this last case, the thresholds are the central value of the given bins. 378 Default: False. 379 380 Returns 381 ------- 382 v1s : np.array 383 Second normalised Minkowski functional evaluated at thresholds us within the given mask. 384 385 """ 386 if self.grad_phi == None: 387 self.get_gradient() 388 389 us = np.atleast_1d(us) 390 391 if edges: 392 dus = (us[1:]-us[:-1]) 393 us = (us[1:]+us[:-1])/2. 394 else: 395 if us.shape == (1,): 396 dus = np.array([0.1]) 397 else: 398 dus = (us[1]-us[0])*np.ones(us.shape[0]) 399 if not (np.isclose(us[1:]-us[:-1], dus[0])).all(): 400 raise ValueError('The thresholds distributions is not uniform. Please set edges=True') 401 402 return np.array([self.V1_iter(u, du) for (u, du) in zip(tqdm(us), dus)]) 403 404 405 def V2_pixel(self, u, du): 406 """Compute the geodesic curvature multiplied by the modulus of the gradient where the values of the input map are between `u-du/2` and `u+du/2`. 407 408 Parameters 409 ---------- 410 u : float 411 The centered value of the bin considered for the computation. 412 413 du : float 414 The width of the bin considered for the computation. 415 416 Returns 417 ------- 418 v2map : np.array 419 Geodesic curvature multiplied by the modulus of the gradient where `u-du/2 < Smap < u+du/2`. 420 421 """ 422 theta = get_theta(self.nside) 423 424 mask = (u + du/2. > self.Smap) & (u -du/2. <= self.Smap) & ~(np.isclose(np.cos(theta),0, atol=1.e-2)) 425 pixs = np.arange(12*self.nside**2)[mask] 426 areas = np.zeros_like(theta) 427 areas[mask] = self.get_κ(pixs) 428 429 return (areas/du / (2.*np.pi)) 430 431 432 def V2_iter(self, u, du): 433 """Compute the normalised third Minkowski functional v2 in the bin `u-du/2` and `u+du/2` within the given mask. 434 435 Parameters 436 ---------- 437 u : float 438 The centered value of the bin considered for the computation of v2. 439 440 du : float 441 The width of the bin considered for the computation of v2. 442 443 Returns 444 ------- 445 v2 : np.array 446 Third normalised Minkowski functional evaluated in the bin u-du/2 and u+du/2 within the given mask. 447 448 """ 449 return np.mean(self.V2_pixel(u, du)[self.mask]) 450 451 def V2(self, us, edges=False): 452 """Compute the normalised third Minkowski functional v2 at the different thresholds us within the given mask. 453 454 Parameters 455 ---------- 456 us : np.array 457 The thresholds considered for the computation of v2. See 'edges' for details. 458 459 edges : bool, optional 460 If False, us is considered as an array of uniformly distributed thresholds. 461 If True, us is considered as a monotonically increasing array of bin edges, including the rightmost edge, allowing for non-uniform distributions of thresholds. 462 In this last case, the thresholds are the central value of the given bins. 463 Default: False. 464 465 Returns 466 ------- 467 v2s : np.array 468 Third normalised Minkowski functional evaluated at thresholds us within the given mask. 469 470 """ 471 if self.grad_phi == None: 472 self.get_gradient() 473 if self.der_phi_phi == None: 474 self.get_hessian() 475 476 us = np.atleast_1d(us) 477 478 if edges: 479 dus = (us[1:]-us[:-1]) 480 us = (us[1:]+us[:-1])/2. 481 else: 482 if us.shape == (1,): 483 dus = np.array([0.1]) 484 else: 485 dus = (us[1]-us[0])*np.ones(us.shape[0]) 486 if not (np.isclose(us[1:]-us[:-1], dus[0])).all(): 487 raise ValueError('The thresholds distributions is not uniform. Please set edges=True') 488 489 return np.array([self.V2_iter(u, du) for (u, du) in zip(tqdm(us), dus)]) 490 491 492 493 def get_maxima(self): 494 """Find the local maxima of the input scalar map. 495 496 Returns 497 ------- 498 pixels : np.array 499 Indices of the pixels which are local maxima. 500 501 values : np.array 502 Values of input map which are local maxima. 503 504 """ 505 neigh = hp.get_all_neighbours(self.nside, np.arange(12*self.nside**2)) 506 507 extT = np.concatenate([self.Smap, [np.min(self.Smap)-1.]]) 508 neigh_matrix = extT[neigh] 509 510 mask = np.all(self.Smap > neigh_matrix, axis=0) 511 pixels = np.argwhere(mask).flatten() 512 values = self.Smap[pixels].flatten() 513 514 return(pixels, values) 515 516 def get_minima(self): 517 """Find the local minima of the input scalar map. 518 519 Returns 520 ------- 521 pixels : np.array 522 Indices of the pixels which are local minima 523 524 values : np.array 525 Values of input map which are local minima 526 527 """ 528 self.Smap = -self.Smap 529 pixels, values = self.get_maxima() 530 self.Smap = -self.Smap 531 532 return(pixels, -values) 533 534 535__all__ = ["Scalar"] 536__docformat__ = "numpy"
105class Scalar(): 106 """Class to compute Minkowski functionals (MFs) and extrema of Healpix scalar maps. It computes and stores spatial first and 107 second derivatives of Healpix scalar maps. 108 109 Parameters 110 ---------- 111 Smap : np.array 112 The input Healpix map where all the statistics will be computed. 113 114 normalise : bool, optional 115 If True, divide input Smap by its standard deviation. Default: True. 116 117 mask : np.array, optional 118 An input Healpix mask. All the statiscal quantities will be computed only within this mask. Default: None (all sky is considered). 119 120 Attributes 121 ---------- 122 Smap : np.array 123 The Healpix map where all the statistics are computed. 124 125 nside : int 126 The Nside parameter of Smap 127 128 mask : np.array 129 A Healpix mask whithin which all the statiscal quantities are computed. 130 131 grad_phi : np.array 132 The phi covariant first derivative of Smap (only if first or second MFs are computed). 133 134 grad_theta : np.array 135 The theta covariant first derivative of Smap (only if first or second MFs are computed). 136 137 der_phi : np.array 138 The phi partial first derivative of Smap (only if first or second MFs are computed). 139 140 der_theta_phi : np.array 141 The second covariant derivative wrt theta and phi of Smap (only if second MF is computed). 142 143 der_phi_phi : np.array 144 The second covariant derivative wrt phi of Smap (only if second MF is computed). 145 146 der_theta_theta : np.array 147 The second covariant derivative wrt theta of Smap (only if second MF is computed). 148 149 Notes 150 ----- 151 The derivatives are always computed full-sky regardless of input mask. 152 153 """ 154 def __init__(self, Smap, normalise=True, mask=None): 155 """Initialise the class to compute Minkowski functionals (MFs) and extrema of Healpix scalar maps. 156 157 Parameters 158 ---------- 159 Smap : np.array 160 The input Healpix map where all the statistics will be computed. 161 162 normalise : bool, optional 163 If True, divide input Smap by its standard deviation. Default: True. 164 165 mask : np.array, optional 166 An input Healpix mask. All the statiscal quantities will be computed only within this mask. Default: None (all sky is considered). 167 168 169 """ 170 self.Smap = Smap.copy() 171 self.nside = hp.get_nside(Smap) 172 if mask is None: 173 self.mask = np.ones(Smap.shape, dtype='bool') 174 else: 175 self.mask = mask 176 if hp.get_nside(mask) != self.nside: 177 raise ValueError('The map and the mask have different nside') 178 if normalise: 179 σ = self.get_variance() 180 self.Smap /= np.sqrt(σ) 181 self.grad_phi = None 182 self.der_phi_phi = None 183 184 185 def __repr__(self): 186 return(f'map = {self.Smap}') 187 188 def get_variance(self): 189 """compute the variance of the input Healpix scalar map within the input mask. 190 191 Returns 192 ------- 193 var : float 194 The variance of the input Healpix map within the input mask. 195 196 """ 197 return (np.var(self.Smap[self.mask])) 198 199 def set_pix(self, pixs): 200 """return the values of the input Healpix scalar map in pixels pixs. 201 202 Parameters 203 ---------- 204 pixs : np.array 205 The indices of the input map pixels whose values of the map are returned. 206 207 Returns 208 ------- 209 values : np.array 210 The values of the input Healpix scalar map in pixels pixs. 211 212 """ 213 return self.Smap[pixs] 214 215 def get_gradient(self): 216 """Compute the covariant and partial first derivatives of the input Healpix scalar map. 217 It stores: 218 219 - first covariant derivative wrt theta in self.grad_theta 220 - first partial derivative wrt phi in self.der_phi 221 - first covariant derivative wrt phi in self.grad_phi 222 223 """ 224 S_grad = derivatives(self.Smap, gradient=True) 225 theta = get_theta(self.nside) 226 self.grad_theta = Scalar(S_grad[0], normalise=False) 227 self.der_phi = Scalar(np.cos(theta) * S_grad[1], normalise=False) 228 self.grad_phi = Scalar(S_grad[1], normalise=False) 229 230 231 def get_hessian(self): 232 """compute the covariant second derivatives of the input Healpix scalar map. 233 It stores: 234 235 - second covariant derivative wrt theta in self.der_theta_theta 236 - second covariant derivative wrt phi in self.der_phi_phi 237 - second covariant derivative wrt theta and phi in self.der_theta_phi 238 239 """ 240 if self.grad_phi == None: 241 self.get_gradient() 242 theta = get_theta(self.nside) 243 244 S_der_der = second_derivatives(self.grad_theta.Smap, self.der_phi.Smap) 245 246 self.der_theta_theta = Scalar(S_der_der[0], normalise=False) 247 self.der_phi_phi = Scalar(S_der_der[1]/np.cos(theta)**2. + self.grad_theta.Smap*np.tan(theta), normalise=False) 248 self.der_theta_phi = Scalar((S_der_der[2]/np.cos(theta) - self.grad_phi.Smap * np.tan(theta)) , normalise=False) 249 250 251 def get_κ(self, pixs): 252 """Compute the geodesic curvature multiplied by the modulus of the gradient in pixels pixs. If not already computed, it computes 253 the first and second covariant derivatives of the input map. 254 255 Parameters 256 ---------- 257 pixs : np.array 258 The indices of the input map pixels where geodesic curvature is computed. 259 260 Returns 261 ------- 262 k : np.array 263 The geodesic curvature in pixels `pixs`. 264 265 """ 266 num = 2.*self.grad_theta.set_pix(pixs)*self.grad_phi.set_pix(pixs)*self.der_theta_phi.set_pix(pixs) - self.grad_phi.set_pix(pixs)**2. * self.der_theta_theta.set_pix(pixs) - self.grad_theta.set_pix(pixs)**2. * self.der_phi_phi.set_pix(pixs) 267 den = self.grad_theta.set_pix(pixs)**2. + self.grad_phi.set_pix(pixs)**2. 268 return num / den 269 270 271 def V0_pixel(self, u): 272 """Determine where input Healpix scalar map Smap is greater than threshold u. 273 274 Parameters 275 ---------- 276 u : float 277 The threshold considered for the computation of first Minkowski functional V0. 278 279 Returns 280 ------- 281 v0map : np.array 282 a bool array with the same shape as the input map, with False where input map values are lower than threshold u. 283 284 """ 285 return self.Smap>u 286 287 def V0_iter(self, u): 288 """Compute the normalised first Minkowski functional v0 at the threshold u within the given mask. 289 290 Parameters 291 ---------- 292 u : float 293 The threshold considered for the computation of v0. 294 295 Returns 296 ------- 297 v0 : np.array 298 First normalised Minkowski functional evaluated at threshold u within the given mask. 299 300 """ 301 return np.mean(self.V0_pixel(u)[self.mask]) 302 303 def V0(self, us): 304 """Compute the normalised first Minkowski functional v0 at the different thresholds us within the given mask. 305 306 Parameters 307 ---------- 308 us : np.array 309 The thresholds considered for the computation of v0. 310 311 Returns 312 ------- 313 v0s : np.array 314 First normalised Minkowski functional evaluated at thresholds us within the given mask. 315 316 """ 317 us = np.atleast_1d(us) 318 return np.array([self.V0_iter(u) for u in tqdm(us)]) 319 320 321 def V1_pixel(self, u, du): 322 """Compute the modulus of the gradient where the values of the input map are between `u-du/2` and `u+du/2`. 323 324 Parameters 325 ---------- 326 u : float 327 The centered value of the bin considered for the computation. 328 329 du : float 330 The width of the bin considered for the computation. 331 332 Returns 333 ------- 334 v1map : np.array 335 Modulus of the gradient where `u-du/2 < Smap < u+du/2`. 336 337 """ 338 theta = get_theta(self.nside) 339 340 areas = np.zeros_like(theta) 341 342 mask = (u + du/2. > self.Smap) & (u - du/2. <= self.Smap) & ~(np.isclose(np.cos(theta),0, atol=1.e-2)) 343 pixs = np.arange(12*self.nside**2)[mask] 344 345 areas[pixs] = np.sqrt(self.grad_theta.set_pix(pixs)**2. + self.grad_phi.set_pix(pixs)**2.) 346 return areas/du / 4. 347 348 def V1_iter(self, u, du): 349 """Compute the normalised second Minkowski functional v1 in the bin `u-du/2` and `u+du/2` within the given mask. 350 351 Parameters 352 ---------- 353 u : float 354 The centered value of the bin considered for the computation of v1. 355 356 du : float 357 The width of the bin considered for the computation of v1. 358 359 Returns 360 ------- 361 v1 : float 362 Second normalised Minkowski functional evaluated in the bin u-du/2 and u+du/2 within the given mask. 363 364 """ 365 return np.mean(self.V1_pixel(u, du)[self.mask]) 366 367 def V1(self, us, edges=False): 368 """Compute the normalised second Minkowski functional v1 at the different thresholds us within the given mask. 369 370 Parameters 371 ---------- 372 us : np.array 373 The thresholds considered for the computation of v2. See 'edges' for details. 374 375 edges : bool, optional 376 If False, us is considered as an array of uniformly distributed thresholds. 377 If True, us is considered as a monotonically increasing array of bin edges, including the rightmost edge, allowing for non-uniform distributions of thresholds. 378 In this last case, the thresholds are the central value of the given bins. 379 Default: False. 380 381 Returns 382 ------- 383 v1s : np.array 384 Second normalised Minkowski functional evaluated at thresholds us within the given mask. 385 386 """ 387 if self.grad_phi == None: 388 self.get_gradient() 389 390 us = np.atleast_1d(us) 391 392 if edges: 393 dus = (us[1:]-us[:-1]) 394 us = (us[1:]+us[:-1])/2. 395 else: 396 if us.shape == (1,): 397 dus = np.array([0.1]) 398 else: 399 dus = (us[1]-us[0])*np.ones(us.shape[0]) 400 if not (np.isclose(us[1:]-us[:-1], dus[0])).all(): 401 raise ValueError('The thresholds distributions is not uniform. Please set edges=True') 402 403 return np.array([self.V1_iter(u, du) for (u, du) in zip(tqdm(us), dus)]) 404 405 406 def V2_pixel(self, u, du): 407 """Compute the geodesic curvature multiplied by the modulus of the gradient where the values of the input map are between `u-du/2` and `u+du/2`. 408 409 Parameters 410 ---------- 411 u : float 412 The centered value of the bin considered for the computation. 413 414 du : float 415 The width of the bin considered for the computation. 416 417 Returns 418 ------- 419 v2map : np.array 420 Geodesic curvature multiplied by the modulus of the gradient where `u-du/2 < Smap < u+du/2`. 421 422 """ 423 theta = get_theta(self.nside) 424 425 mask = (u + du/2. > self.Smap) & (u -du/2. <= self.Smap) & ~(np.isclose(np.cos(theta),0, atol=1.e-2)) 426 pixs = np.arange(12*self.nside**2)[mask] 427 areas = np.zeros_like(theta) 428 areas[mask] = self.get_κ(pixs) 429 430 return (areas/du / (2.*np.pi)) 431 432 433 def V2_iter(self, u, du): 434 """Compute the normalised third Minkowski functional v2 in the bin `u-du/2` and `u+du/2` within the given mask. 435 436 Parameters 437 ---------- 438 u : float 439 The centered value of the bin considered for the computation of v2. 440 441 du : float 442 The width of the bin considered for the computation of v2. 443 444 Returns 445 ------- 446 v2 : np.array 447 Third normalised Minkowski functional evaluated in the bin u-du/2 and u+du/2 within the given mask. 448 449 """ 450 return np.mean(self.V2_pixel(u, du)[self.mask]) 451 452 def V2(self, us, edges=False): 453 """Compute the normalised third Minkowski functional v2 at the different thresholds us within the given mask. 454 455 Parameters 456 ---------- 457 us : np.array 458 The thresholds considered for the computation of v2. See 'edges' for details. 459 460 edges : bool, optional 461 If False, us is considered as an array of uniformly distributed thresholds. 462 If True, us is considered as a monotonically increasing array of bin edges, including the rightmost edge, allowing for non-uniform distributions of thresholds. 463 In this last case, the thresholds are the central value of the given bins. 464 Default: False. 465 466 Returns 467 ------- 468 v2s : np.array 469 Third normalised Minkowski functional evaluated at thresholds us within the given mask. 470 471 """ 472 if self.grad_phi == None: 473 self.get_gradient() 474 if self.der_phi_phi == None: 475 self.get_hessian() 476 477 us = np.atleast_1d(us) 478 479 if edges: 480 dus = (us[1:]-us[:-1]) 481 us = (us[1:]+us[:-1])/2. 482 else: 483 if us.shape == (1,): 484 dus = np.array([0.1]) 485 else: 486 dus = (us[1]-us[0])*np.ones(us.shape[0]) 487 if not (np.isclose(us[1:]-us[:-1], dus[0])).all(): 488 raise ValueError('The thresholds distributions is not uniform. Please set edges=True') 489 490 return np.array([self.V2_iter(u, du) for (u, du) in zip(tqdm(us), dus)]) 491 492 493 494 def get_maxima(self): 495 """Find the local maxima of the input scalar map. 496 497 Returns 498 ------- 499 pixels : np.array 500 Indices of the pixels which are local maxima. 501 502 values : np.array 503 Values of input map which are local maxima. 504 505 """ 506 neigh = hp.get_all_neighbours(self.nside, np.arange(12*self.nside**2)) 507 508 extT = np.concatenate([self.Smap, [np.min(self.Smap)-1.]]) 509 neigh_matrix = extT[neigh] 510 511 mask = np.all(self.Smap > neigh_matrix, axis=0) 512 pixels = np.argwhere(mask).flatten() 513 values = self.Smap[pixels].flatten() 514 515 return(pixels, values) 516 517 def get_minima(self): 518 """Find the local minima of the input scalar map. 519 520 Returns 521 ------- 522 pixels : np.array 523 Indices of the pixels which are local minima 524 525 values : np.array 526 Values of input map which are local minima 527 528 """ 529 self.Smap = -self.Smap 530 pixels, values = self.get_maxima() 531 self.Smap = -self.Smap 532 533 return(pixels, -values)
Class to compute Minkowski functionals (MFs) and extrema of Healpix scalar maps. It computes and stores spatial first and second derivatives of Healpix scalar maps.
Parameters
- Smap (np.array): The input Healpix map where all the statistics will be computed.
- normalise (bool, optional): If True, divide input Smap by its standard deviation. Default: True.
- mask (np.array, optional): An input Healpix mask. All the statiscal quantities will be computed only within this mask. Default: None (all sky is considered).
Attributes
- Smap (np.array): The Healpix map where all the statistics are computed.
- nside (int): The Nside parameter of Smap
- mask (np.array): A Healpix mask whithin which all the statiscal quantities are computed.
- grad_phi (np.array): The phi covariant first derivative of Smap (only if first or second MFs are computed).
- grad_theta (np.array): The theta covariant first derivative of Smap (only if first or second MFs are computed).
- der_phi (np.array): The phi partial first derivative of Smap (only if first or second MFs are computed).
- der_theta_phi (np.array): The second covariant derivative wrt theta and phi of Smap (only if second MF is computed).
- der_phi_phi (np.array): The second covariant derivative wrt phi of Smap (only if second MF is computed).
- der_theta_theta (np.array): The second covariant derivative wrt theta of Smap (only if second MF is computed).
Notes
The derivatives are always computed full-sky regardless of input mask.
154 def __init__(self, Smap, normalise=True, mask=None): 155 """Initialise the class to compute Minkowski functionals (MFs) and extrema of Healpix scalar maps. 156 157 Parameters 158 ---------- 159 Smap : np.array 160 The input Healpix map where all the statistics will be computed. 161 162 normalise : bool, optional 163 If True, divide input Smap by its standard deviation. Default: True. 164 165 mask : np.array, optional 166 An input Healpix mask. All the statiscal quantities will be computed only within this mask. Default: None (all sky is considered). 167 168 169 """ 170 self.Smap = Smap.copy() 171 self.nside = hp.get_nside(Smap) 172 if mask is None: 173 self.mask = np.ones(Smap.shape, dtype='bool') 174 else: 175 self.mask = mask 176 if hp.get_nside(mask) != self.nside: 177 raise ValueError('The map and the mask have different nside') 178 if normalise: 179 σ = self.get_variance() 180 self.Smap /= np.sqrt(σ) 181 self.grad_phi = None 182 self.der_phi_phi = None
Initialise the class to compute Minkowski functionals (MFs) and extrema of Healpix scalar maps.
Parameters
- Smap (np.array): The input Healpix map where all the statistics will be computed.
- normalise (bool, optional): If True, divide input Smap by its standard deviation. Default: True.
- mask (np.array, optional): An input Healpix mask. All the statiscal quantities will be computed only within this mask. Default: None (all sky is considered).
188 def get_variance(self): 189 """compute the variance of the input Healpix scalar map within the input mask. 190 191 Returns 192 ------- 193 var : float 194 The variance of the input Healpix map within the input mask. 195 196 """ 197 return (np.var(self.Smap[self.mask]))
compute the variance of the input Healpix scalar map within the input mask.
Returns
- var (float): The variance of the input Healpix map within the input mask.
199 def set_pix(self, pixs): 200 """return the values of the input Healpix scalar map in pixels pixs. 201 202 Parameters 203 ---------- 204 pixs : np.array 205 The indices of the input map pixels whose values of the map are returned. 206 207 Returns 208 ------- 209 values : np.array 210 The values of the input Healpix scalar map in pixels pixs. 211 212 """ 213 return self.Smap[pixs]
return the values of the input Healpix scalar map in pixels pixs.
Parameters
- pixs (np.array): The indices of the input map pixels whose values of the map are returned.
Returns
- values (np.array): The values of the input Healpix scalar map in pixels pixs.
215 def get_gradient(self): 216 """Compute the covariant and partial first derivatives of the input Healpix scalar map. 217 It stores: 218 219 - first covariant derivative wrt theta in self.grad_theta 220 - first partial derivative wrt phi in self.der_phi 221 - first covariant derivative wrt phi in self.grad_phi 222 223 """ 224 S_grad = derivatives(self.Smap, gradient=True) 225 theta = get_theta(self.nside) 226 self.grad_theta = Scalar(S_grad[0], normalise=False) 227 self.der_phi = Scalar(np.cos(theta) * S_grad[1], normalise=False) 228 self.grad_phi = Scalar(S_grad[1], normalise=False)
Compute the covariant and partial first derivatives of the input Healpix scalar map. It stores:
- first covariant derivative wrt theta in self.grad_theta
- first partial derivative wrt phi in self.der_phi
- first covariant derivative wrt phi in self.grad_phi
231 def get_hessian(self): 232 """compute the covariant second derivatives of the input Healpix scalar map. 233 It stores: 234 235 - second covariant derivative wrt theta in self.der_theta_theta 236 - second covariant derivative wrt phi in self.der_phi_phi 237 - second covariant derivative wrt theta and phi in self.der_theta_phi 238 239 """ 240 if self.grad_phi == None: 241 self.get_gradient() 242 theta = get_theta(self.nside) 243 244 S_der_der = second_derivatives(self.grad_theta.Smap, self.der_phi.Smap) 245 246 self.der_theta_theta = Scalar(S_der_der[0], normalise=False) 247 self.der_phi_phi = Scalar(S_der_der[1]/np.cos(theta)**2. + self.grad_theta.Smap*np.tan(theta), normalise=False) 248 self.der_theta_phi = Scalar((S_der_der[2]/np.cos(theta) - self.grad_phi.Smap * np.tan(theta)) , normalise=False)
compute the covariant second derivatives of the input Healpix scalar map. It stores:
- second covariant derivative wrt theta in self.der_theta_theta
- second covariant derivative wrt phi in self.der_phi_phi
- second covariant derivative wrt theta and phi in self.der_theta_phi
251 def get_κ(self, pixs): 252 """Compute the geodesic curvature multiplied by the modulus of the gradient in pixels pixs. If not already computed, it computes 253 the first and second covariant derivatives of the input map. 254 255 Parameters 256 ---------- 257 pixs : np.array 258 The indices of the input map pixels where geodesic curvature is computed. 259 260 Returns 261 ------- 262 k : np.array 263 The geodesic curvature in pixels `pixs`. 264 265 """ 266 num = 2.*self.grad_theta.set_pix(pixs)*self.grad_phi.set_pix(pixs)*self.der_theta_phi.set_pix(pixs) - self.grad_phi.set_pix(pixs)**2. * self.der_theta_theta.set_pix(pixs) - self.grad_theta.set_pix(pixs)**2. * self.der_phi_phi.set_pix(pixs) 267 den = self.grad_theta.set_pix(pixs)**2. + self.grad_phi.set_pix(pixs)**2. 268 return num / den
Compute the geodesic curvature multiplied by the modulus of the gradient in pixels pixs. If not already computed, it computes the first and second covariant derivatives of the input map.
Parameters
- pixs (np.array): The indices of the input map pixels where geodesic curvature is computed.
Returns
- k (np.array):
The geodesic curvature in pixels
pixs
.
271 def V0_pixel(self, u): 272 """Determine where input Healpix scalar map Smap is greater than threshold u. 273 274 Parameters 275 ---------- 276 u : float 277 The threshold considered for the computation of first Minkowski functional V0. 278 279 Returns 280 ------- 281 v0map : np.array 282 a bool array with the same shape as the input map, with False where input map values are lower than threshold u. 283 284 """ 285 return self.Smap>u
Determine where input Healpix scalar map Smap is greater than threshold u.
Parameters
- u (float): The threshold considered for the computation of first Minkowski functional V0.
Returns
- v0map (np.array): a bool array with the same shape as the input map, with False where input map values are lower than threshold u.
287 def V0_iter(self, u): 288 """Compute the normalised first Minkowski functional v0 at the threshold u within the given mask. 289 290 Parameters 291 ---------- 292 u : float 293 The threshold considered for the computation of v0. 294 295 Returns 296 ------- 297 v0 : np.array 298 First normalised Minkowski functional evaluated at threshold u within the given mask. 299 300 """ 301 return np.mean(self.V0_pixel(u)[self.mask])
Compute the normalised first Minkowski functional v0 at the threshold u within the given mask.
Parameters
- u (float): The threshold considered for the computation of v0.
Returns
- v0 (np.array): First normalised Minkowski functional evaluated at threshold u within the given mask.
303 def V0(self, us): 304 """Compute the normalised first Minkowski functional v0 at the different thresholds us within the given mask. 305 306 Parameters 307 ---------- 308 us : np.array 309 The thresholds considered for the computation of v0. 310 311 Returns 312 ------- 313 v0s : np.array 314 First normalised Minkowski functional evaluated at thresholds us within the given mask. 315 316 """ 317 us = np.atleast_1d(us) 318 return np.array([self.V0_iter(u) for u in tqdm(us)])
Compute the normalised first Minkowski functional v0 at the different thresholds us within the given mask.
Parameters
- us (np.array): The thresholds considered for the computation of v0.
Returns
- v0s (np.array): First normalised Minkowski functional evaluated at thresholds us within the given mask.
321 def V1_pixel(self, u, du): 322 """Compute the modulus of the gradient where the values of the input map are between `u-du/2` and `u+du/2`. 323 324 Parameters 325 ---------- 326 u : float 327 The centered value of the bin considered for the computation. 328 329 du : float 330 The width of the bin considered for the computation. 331 332 Returns 333 ------- 334 v1map : np.array 335 Modulus of the gradient where `u-du/2 < Smap < u+du/2`. 336 337 """ 338 theta = get_theta(self.nside) 339 340 areas = np.zeros_like(theta) 341 342 mask = (u + du/2. > self.Smap) & (u - du/2. <= self.Smap) & ~(np.isclose(np.cos(theta),0, atol=1.e-2)) 343 pixs = np.arange(12*self.nside**2)[mask] 344 345 areas[pixs] = np.sqrt(self.grad_theta.set_pix(pixs)**2. + self.grad_phi.set_pix(pixs)**2.) 346 return areas/du / 4.
Compute the modulus of the gradient where the values of the input map are between u-du/2
and u+du/2
.
Parameters
- u (float): The centered value of the bin considered for the computation.
- du (float): The width of the bin considered for the computation.
Returns
- v1map (np.array):
Modulus of the gradient where
u-du/2 < Smap < u+du/2
.
348 def V1_iter(self, u, du): 349 """Compute the normalised second Minkowski functional v1 in the bin `u-du/2` and `u+du/2` within the given mask. 350 351 Parameters 352 ---------- 353 u : float 354 The centered value of the bin considered for the computation of v1. 355 356 du : float 357 The width of the bin considered for the computation of v1. 358 359 Returns 360 ------- 361 v1 : float 362 Second normalised Minkowski functional evaluated in the bin u-du/2 and u+du/2 within the given mask. 363 364 """ 365 return np.mean(self.V1_pixel(u, du)[self.mask])
Compute the normalised second Minkowski functional v1 in the bin u-du/2
and u+du/2
within the given mask.
Parameters
- u (float): The centered value of the bin considered for the computation of v1.
- du (float): The width of the bin considered for the computation of v1.
Returns
- v1 (float): Second normalised Minkowski functional evaluated in the bin u-du/2 and u+du/2 within the given mask.
367 def V1(self, us, edges=False): 368 """Compute the normalised second Minkowski functional v1 at the different thresholds us within the given mask. 369 370 Parameters 371 ---------- 372 us : np.array 373 The thresholds considered for the computation of v2. See 'edges' for details. 374 375 edges : bool, optional 376 If False, us is considered as an array of uniformly distributed thresholds. 377 If True, us is considered as a monotonically increasing array of bin edges, including the rightmost edge, allowing for non-uniform distributions of thresholds. 378 In this last case, the thresholds are the central value of the given bins. 379 Default: False. 380 381 Returns 382 ------- 383 v1s : np.array 384 Second normalised Minkowski functional evaluated at thresholds us within the given mask. 385 386 """ 387 if self.grad_phi == None: 388 self.get_gradient() 389 390 us = np.atleast_1d(us) 391 392 if edges: 393 dus = (us[1:]-us[:-1]) 394 us = (us[1:]+us[:-1])/2. 395 else: 396 if us.shape == (1,): 397 dus = np.array([0.1]) 398 else: 399 dus = (us[1]-us[0])*np.ones(us.shape[0]) 400 if not (np.isclose(us[1:]-us[:-1], dus[0])).all(): 401 raise ValueError('The thresholds distributions is not uniform. Please set edges=True') 402 403 return np.array([self.V1_iter(u, du) for (u, du) in zip(tqdm(us), dus)])
Compute the normalised second Minkowski functional v1 at the different thresholds us within the given mask.
Parameters
- us (np.array): The thresholds considered for the computation of v2. See 'edges' for details.
- edges (bool, optional): If False, us is considered as an array of uniformly distributed thresholds. If True, us is considered as a monotonically increasing array of bin edges, including the rightmost edge, allowing for non-uniform distributions of thresholds. In this last case, the thresholds are the central value of the given bins. Default: False.
Returns
- v1s (np.array): Second normalised Minkowski functional evaluated at thresholds us within the given mask.
406 def V2_pixel(self, u, du): 407 """Compute the geodesic curvature multiplied by the modulus of the gradient where the values of the input map are between `u-du/2` and `u+du/2`. 408 409 Parameters 410 ---------- 411 u : float 412 The centered value of the bin considered for the computation. 413 414 du : float 415 The width of the bin considered for the computation. 416 417 Returns 418 ------- 419 v2map : np.array 420 Geodesic curvature multiplied by the modulus of the gradient where `u-du/2 < Smap < u+du/2`. 421 422 """ 423 theta = get_theta(self.nside) 424 425 mask = (u + du/2. > self.Smap) & (u -du/2. <= self.Smap) & ~(np.isclose(np.cos(theta),0, atol=1.e-2)) 426 pixs = np.arange(12*self.nside**2)[mask] 427 areas = np.zeros_like(theta) 428 areas[mask] = self.get_κ(pixs) 429 430 return (areas/du / (2.*np.pi))
Compute the geodesic curvature multiplied by the modulus of the gradient where the values of the input map are between u-du/2
and u+du/2
.
Parameters
- u (float): The centered value of the bin considered for the computation.
- du (float): The width of the bin considered for the computation.
Returns
- v2map (np.array):
Geodesic curvature multiplied by the modulus of the gradient where
u-du/2 < Smap < u+du/2
.
433 def V2_iter(self, u, du): 434 """Compute the normalised third Minkowski functional v2 in the bin `u-du/2` and `u+du/2` within the given mask. 435 436 Parameters 437 ---------- 438 u : float 439 The centered value of the bin considered for the computation of v2. 440 441 du : float 442 The width of the bin considered for the computation of v2. 443 444 Returns 445 ------- 446 v2 : np.array 447 Third normalised Minkowski functional evaluated in the bin u-du/2 and u+du/2 within the given mask. 448 449 """ 450 return np.mean(self.V2_pixel(u, du)[self.mask])
Compute the normalised third Minkowski functional v2 in the bin u-du/2
and u+du/2
within the given mask.
Parameters
- u (float): The centered value of the bin considered for the computation of v2.
- du (float): The width of the bin considered for the computation of v2.
Returns
- v2 (np.array): Third normalised Minkowski functional evaluated in the bin u-du/2 and u+du/2 within the given mask.
452 def V2(self, us, edges=False): 453 """Compute the normalised third Minkowski functional v2 at the different thresholds us within the given mask. 454 455 Parameters 456 ---------- 457 us : np.array 458 The thresholds considered for the computation of v2. See 'edges' for details. 459 460 edges : bool, optional 461 If False, us is considered as an array of uniformly distributed thresholds. 462 If True, us is considered as a monotonically increasing array of bin edges, including the rightmost edge, allowing for non-uniform distributions of thresholds. 463 In this last case, the thresholds are the central value of the given bins. 464 Default: False. 465 466 Returns 467 ------- 468 v2s : np.array 469 Third normalised Minkowski functional evaluated at thresholds us within the given mask. 470 471 """ 472 if self.grad_phi == None: 473 self.get_gradient() 474 if self.der_phi_phi == None: 475 self.get_hessian() 476 477 us = np.atleast_1d(us) 478 479 if edges: 480 dus = (us[1:]-us[:-1]) 481 us = (us[1:]+us[:-1])/2. 482 else: 483 if us.shape == (1,): 484 dus = np.array([0.1]) 485 else: 486 dus = (us[1]-us[0])*np.ones(us.shape[0]) 487 if not (np.isclose(us[1:]-us[:-1], dus[0])).all(): 488 raise ValueError('The thresholds distributions is not uniform. Please set edges=True') 489 490 return np.array([self.V2_iter(u, du) for (u, du) in zip(tqdm(us), dus)])
Compute the normalised third Minkowski functional v2 at the different thresholds us within the given mask.
Parameters
- us (np.array): The thresholds considered for the computation of v2. See 'edges' for details.
- edges (bool, optional): If False, us is considered as an array of uniformly distributed thresholds. If True, us is considered as a monotonically increasing array of bin edges, including the rightmost edge, allowing for non-uniform distributions of thresholds. In this last case, the thresholds are the central value of the given bins. Default: False.
Returns
- v2s (np.array): Third normalised Minkowski functional evaluated at thresholds us within the given mask.
494 def get_maxima(self): 495 """Find the local maxima of the input scalar map. 496 497 Returns 498 ------- 499 pixels : np.array 500 Indices of the pixels which are local maxima. 501 502 values : np.array 503 Values of input map which are local maxima. 504 505 """ 506 neigh = hp.get_all_neighbours(self.nside, np.arange(12*self.nside**2)) 507 508 extT = np.concatenate([self.Smap, [np.min(self.Smap)-1.]]) 509 neigh_matrix = extT[neigh] 510 511 mask = np.all(self.Smap > neigh_matrix, axis=0) 512 pixels = np.argwhere(mask).flatten() 513 values = self.Smap[pixels].flatten() 514 515 return(pixels, values)
Find the local maxima of the input scalar map.
Returns
- pixels (np.array): Indices of the pixels which are local maxima.
- values (np.array): Values of input map which are local maxima.
517 def get_minima(self): 518 """Find the local minima of the input scalar map. 519 520 Returns 521 ------- 522 pixels : np.array 523 Indices of the pixels which are local minima 524 525 values : np.array 526 Values of input map which are local minima 527 528 """ 529 self.Smap = -self.Smap 530 pixels, values = self.get_maxima() 531 self.Smap = -self.Smap 532 533 return(pixels, -values)
Find the local minima of the input scalar map.
Returns
- pixels (np.array): Indices of the pixels which are local minima
- values (np.array): Values of input map which are local minima