1# Licensed under a 3-clause BSD style license - see LICENSE.rst 2 3import math 4 5import numpy as np 6 7from .core import Kernel1D, Kernel2D, Kernel 8from .utils import has_even_axis, raise_even_kernel_exception, KernelSizeError 9from astropy.modeling import models 10from astropy.modeling.core import Fittable1DModel, Fittable2DModel 11from astropy.utils.decorators import deprecated 12 13__all__ = ['Gaussian1DKernel', 'Gaussian2DKernel', 'CustomKernel', 14 'Box1DKernel', 'Box2DKernel', 'Tophat2DKernel', 15 'Trapezoid1DKernel', 'RickerWavelet1DKernel', 'RickerWavelet2DKernel', 16 'AiryDisk2DKernel', 'Moffat2DKernel', 'Model1DKernel', 17 'Model2DKernel', 'TrapezoidDisk2DKernel', 'Ring2DKernel'] 18 19 20def _round_up_to_odd_integer(value): 21 i = math.ceil(value) 22 if i % 2 == 0: 23 return i + 1 24 else: 25 return i 26 27 28class Gaussian1DKernel(Kernel1D): 29 """ 30 1D Gaussian filter kernel. 31 32 The Gaussian filter is a filter with great smoothing properties. It is 33 isotropic and does not produce artifacts. 34 35 The generated kernel is normalized so that it integrates to 1. 36 37 Parameters 38 ---------- 39 stddev : number 40 Standard deviation of the Gaussian kernel. 41 x_size : int, optional 42 Size of the kernel array. Default = ⌊8*stddev+1⌋. 43 mode : str, optional 44 One of the following discretization modes: 45 * 'center' (default) 46 Discretize model by taking the value 47 at the center of the bin. 48 * 'linear_interp' 49 Discretize model by linearly interpolating 50 between the values at the corners of the bin. 51 * 'oversample' 52 Discretize model by taking the average 53 on an oversampled grid. 54 * 'integrate' 55 Discretize model by integrating the 56 model over the bin. Very slow. 57 factor : number, optional 58 Factor of oversampling. Default factor = 10. If the factor 59 is too large, evaluation can be very slow. 60 61 62 See Also 63 -------- 64 Box1DKernel, Trapezoid1DKernel, RickerWavelet1DKernel 65 66 67 Examples 68 -------- 69 Kernel response: 70 71 .. plot:: 72 :include-source: 73 74 import matplotlib.pyplot as plt 75 from astropy.convolution import Gaussian1DKernel 76 gauss_1D_kernel = Gaussian1DKernel(10) 77 plt.plot(gauss_1D_kernel, drawstyle='steps') 78 plt.xlabel('x [pixels]') 79 plt.ylabel('value') 80 plt.show() 81 """ 82 _separable = True 83 _is_bool = False 84 85 def __init__(self, stddev, **kwargs): 86 self._model = models.Gaussian1D(1. / (np.sqrt(2 * np.pi) * stddev), 87 0, stddev) 88 self._default_size = _round_up_to_odd_integer(8 * stddev) 89 super().__init__(**kwargs) 90 self._truncation = np.abs(1. - self._array.sum()) 91 92 93class Gaussian2DKernel(Kernel2D): 94 """ 95 2D Gaussian filter kernel. 96 97 The Gaussian filter is a filter with great smoothing properties. It is 98 isotropic and does not produce artifacts. 99 100 The generated kernel is normalized so that it integrates to 1. 101 102 Parameters 103 ---------- 104 x_stddev : float 105 Standard deviation of the Gaussian in x before rotating by theta. 106 y_stddev : float 107 Standard deviation of the Gaussian in y before rotating by theta. 108 theta : float or `~astropy.units.Quantity` ['angle'] 109 Rotation angle. If passed as a float, it is assumed to be in radians. 110 The rotation angle increases counterclockwise. 111 x_size : int, optional 112 Size in x direction of the kernel array. Default = ⌊8*stddev + 1⌋. 113 y_size : int, optional 114 Size in y direction of the kernel array. Default = ⌊8*stddev + 1⌋. 115 mode : str, optional 116 One of the following discretization modes: 117 * 'center' (default) 118 Discretize model by taking the value 119 at the center of the bin. 120 * 'linear_interp' 121 Discretize model by performing a bilinear interpolation 122 between the values at the corners of the bin. 123 * 'oversample' 124 Discretize model by taking the average 125 on an oversampled grid. 126 * 'integrate' 127 Discretize model by integrating the 128 model over the bin. 129 factor : number, optional 130 Factor of oversampling. Default factor = 10. 131 132 133 See Also 134 -------- 135 Box2DKernel, Tophat2DKernel, RickerWavelet2DKernel, Ring2DKernel, 136 TrapezoidDisk2DKernel, AiryDisk2DKernel, Moffat2DKernel 137 138 Examples 139 -------- 140 Kernel response: 141 142 .. plot:: 143 :include-source: 144 145 import matplotlib.pyplot as plt 146 from astropy.convolution import Gaussian2DKernel 147 gaussian_2D_kernel = Gaussian2DKernel(10) 148 plt.imshow(gaussian_2D_kernel, interpolation='none', origin='lower') 149 plt.xlabel('x [pixels]') 150 plt.ylabel('y [pixels]') 151 plt.colorbar() 152 plt.show() 153 154 """ 155 _separable = True 156 _is_bool = False 157 158 def __init__(self, x_stddev, y_stddev=None, theta=0.0, **kwargs): 159 if y_stddev is None: 160 y_stddev = x_stddev 161 self._model = models.Gaussian2D(1. / (2 * np.pi * x_stddev * y_stddev), 162 0, 0, x_stddev=x_stddev, 163 y_stddev=y_stddev, theta=theta) 164 self._default_size = _round_up_to_odd_integer( 165 8 * np.max([x_stddev, y_stddev])) 166 super().__init__(**kwargs) 167 self._truncation = np.abs(1. - self._array.sum()) 168 169 170class Box1DKernel(Kernel1D): 171 """ 172 1D Box filter kernel. 173 174 The Box filter or running mean is a smoothing filter. It is not isotropic 175 and can produce artifacts when applied repeatedly to the same data. 176 177 The generated kernel is normalized so that it integrates to 1. 178 179 By default the Box kernel uses the ``linear_interp`` discretization mode, 180 which allows non-shifting, even-sized kernels. This is achieved by 181 weighting the edge pixels with 1/2. E.g a Box kernel with an effective 182 smoothing of 4 pixel would have the following array: [0.5, 1, 1, 1, 0.5]. 183 184 185 Parameters 186 ---------- 187 width : number 188 Width of the filter kernel. 189 mode : str, optional 190 One of the following discretization modes: 191 * 'center' 192 Discretize model by taking the value 193 at the center of the bin. 194 * 'linear_interp' (default) 195 Discretize model by linearly interpolating 196 between the values at the corners of the bin. 197 * 'oversample' 198 Discretize model by taking the average 199 on an oversampled grid. 200 * 'integrate' 201 Discretize model by integrating the 202 model over the bin. 203 factor : number, optional 204 Factor of oversampling. Default factor = 10. 205 206 See Also 207 -------- 208 Gaussian1DKernel, Trapezoid1DKernel, RickerWavelet1DKernel 209 210 211 Examples 212 -------- 213 Kernel response function: 214 215 .. plot:: 216 :include-source: 217 218 import matplotlib.pyplot as plt 219 from astropy.convolution import Box1DKernel 220 box_1D_kernel = Box1DKernel(9) 221 plt.plot(box_1D_kernel, drawstyle='steps') 222 plt.xlim(-1, 9) 223 plt.xlabel('x [pixels]') 224 plt.ylabel('value') 225 plt.show() 226 227 """ 228 _separable = True 229 _is_bool = True 230 231 def __init__(self, width, **kwargs): 232 self._model = models.Box1D(1. / width, 0, width) 233 self._default_size = _round_up_to_odd_integer(width) 234 kwargs['mode'] = 'linear_interp' 235 super().__init__(**kwargs) 236 self._truncation = 0 237 self.normalize() 238 239 240class Box2DKernel(Kernel2D): 241 """ 242 2D Box filter kernel. 243 244 The Box filter or running mean is a smoothing filter. It is not isotropic 245 and can produce artifacts when applied repeatedly to the same data. 246 247 The generated kernel is normalized so that it integrates to 1. 248 249 By default the Box kernel uses the ``linear_interp`` discretization mode, 250 which allows non-shifting, even-sized kernels. This is achieved by 251 weighting the edge pixels with 1/2. 252 253 254 Parameters 255 ---------- 256 width : number 257 Width of the filter kernel. 258 mode : str, optional 259 One of the following discretization modes: 260 * 'center' 261 Discretize model by taking the value 262 at the center of the bin. 263 * 'linear_interp' (default) 264 Discretize model by performing a bilinear interpolation 265 between the values at the corners of the bin. 266 * 'oversample' 267 Discretize model by taking the average 268 on an oversampled grid. 269 * 'integrate' 270 Discretize model by integrating the 271 model over the bin. 272 factor : number, optional 273 Factor of oversampling. Default factor = 10. 274 275 276 See Also 277 -------- 278 Gaussian2DKernel, Tophat2DKernel, RickerWavelet2DKernel, Ring2DKernel, 279 TrapezoidDisk2DKernel, AiryDisk2DKernel, Moffat2DKernel 280 281 Examples 282 -------- 283 Kernel response: 284 285 .. plot:: 286 :include-source: 287 288 import matplotlib.pyplot as plt 289 from astropy.convolution import Box2DKernel 290 box_2D_kernel = Box2DKernel(9) 291 plt.imshow(box_2D_kernel, interpolation='none', origin='lower', 292 vmin=0.0, vmax=0.015) 293 plt.xlim(-1, 9) 294 plt.ylim(-1, 9) 295 plt.xlabel('x [pixels]') 296 plt.ylabel('y [pixels]') 297 plt.colorbar() 298 plt.show() 299 """ 300 _separable = True 301 _is_bool = True 302 303 def __init__(self, width, **kwargs): 304 self._model = models.Box2D(1. / width ** 2, 0, 0, width, width) 305 self._default_size = _round_up_to_odd_integer(width) 306 kwargs['mode'] = 'linear_interp' 307 super().__init__(**kwargs) 308 self._truncation = 0 309 self.normalize() 310 311 312class Tophat2DKernel(Kernel2D): 313 """ 314 2D Tophat filter kernel. 315 316 The Tophat filter is an isotropic smoothing filter. It can produce 317 artifacts when applied repeatedly on the same data. 318 319 The generated kernel is normalized so that it integrates to 1. 320 321 Parameters 322 ---------- 323 radius : int 324 Radius of the filter kernel. 325 mode : str, optional 326 One of the following discretization modes: 327 * 'center' (default) 328 Discretize model by taking the value 329 at the center of the bin. 330 * 'linear_interp' 331 Discretize model by performing a bilinear interpolation 332 between the values at the corners of the bin. 333 * 'oversample' 334 Discretize model by taking the average 335 on an oversampled grid. 336 * 'integrate' 337 Discretize model by integrating the 338 model over the bin. 339 factor : number, optional 340 Factor of oversampling. Default factor = 10. 341 342 343 See Also 344 -------- 345 Gaussian2DKernel, Box2DKernel, RickerWavelet2DKernel, Ring2DKernel, 346 TrapezoidDisk2DKernel, AiryDisk2DKernel, Moffat2DKernel 347 348 Examples 349 -------- 350 Kernel response: 351 352 .. plot:: 353 :include-source: 354 355 import matplotlib.pyplot as plt 356 from astropy.convolution import Tophat2DKernel 357 tophat_2D_kernel = Tophat2DKernel(40) 358 plt.imshow(tophat_2D_kernel, interpolation='none', origin='lower') 359 plt.xlabel('x [pixels]') 360 plt.ylabel('y [pixels]') 361 plt.colorbar() 362 plt.show() 363 364 """ 365 def __init__(self, radius, **kwargs): 366 self._model = models.Disk2D(1. / (np.pi * radius ** 2), 0, 0, radius) 367 self._default_size = _round_up_to_odd_integer(2 * radius) 368 super().__init__(**kwargs) 369 self._truncation = 0 370 371 372class Ring2DKernel(Kernel2D): 373 """ 374 2D Ring filter kernel. 375 376 The Ring filter kernel is the difference between two Tophat kernels of 377 different width. This kernel is useful for, e.g., background estimation. 378 379 The generated kernel is normalized so that it integrates to 1. 380 381 Parameters 382 ---------- 383 radius_in : number 384 Inner radius of the ring kernel. 385 width : number 386 Width of the ring kernel. 387 mode : str, optional 388 One of the following discretization modes: 389 * 'center' (default) 390 Discretize model by taking the value 391 at the center of the bin. 392 * 'linear_interp' 393 Discretize model by performing a bilinear interpolation 394 between the values at the corners of the bin. 395 * 'oversample' 396 Discretize model by taking the average 397 on an oversampled grid. 398 * 'integrate' 399 Discretize model by integrating the 400 model over the bin. 401 factor : number, optional 402 Factor of oversampling. Default factor = 10. 403 404 See Also 405 -------- 406 Gaussian2DKernel, Box2DKernel, Tophat2DKernel, RickerWavelet2DKernel, 407 TrapezoidDisk2DKernel, AiryDisk2DKernel, Moffat2DKernel 408 409 Examples 410 -------- 411 Kernel response: 412 413 .. plot:: 414 :include-source: 415 416 import matplotlib.pyplot as plt 417 from astropy.convolution import Ring2DKernel 418 ring_2D_kernel = Ring2DKernel(9, 8) 419 plt.imshow(ring_2D_kernel, interpolation='none', origin='lower') 420 plt.xlabel('x [pixels]') 421 plt.ylabel('y [pixels]') 422 plt.colorbar() 423 plt.show() 424 """ 425 def __init__(self, radius_in, width, **kwargs): 426 radius_out = radius_in + width 427 self._model = models.Ring2D(1. / (np.pi * (radius_out ** 2 - radius_in ** 2)), 428 0, 0, radius_in, width) 429 self._default_size = _round_up_to_odd_integer(2 * radius_out) 430 super().__init__(**kwargs) 431 self._truncation = 0 432 433 434class Trapezoid1DKernel(Kernel1D): 435 """ 436 1D trapezoid kernel. 437 438 The generated kernel is normalized so that it integrates to 1. 439 440 Parameters 441 ---------- 442 width : number 443 Width of the filter kernel, defined as the width of the constant part, 444 before it begins to slope down. 445 slope : number 446 Slope of the filter kernel's tails 447 mode : str, optional 448 One of the following discretization modes: 449 * 'center' (default) 450 Discretize model by taking the value 451 at the center of the bin. 452 * 'linear_interp' 453 Discretize model by linearly interpolating 454 between the values at the corners of the bin. 455 * 'oversample' 456 Discretize model by taking the average 457 on an oversampled grid. 458 * 'integrate' 459 Discretize model by integrating the 460 model over the bin. 461 factor : number, optional 462 Factor of oversampling. Default factor = 10. 463 464 See Also 465 -------- 466 Box1DKernel, Gaussian1DKernel, RickerWavelet1DKernel 467 468 Examples 469 -------- 470 Kernel response: 471 472 .. plot:: 473 :include-source: 474 475 import matplotlib.pyplot as plt 476 from astropy.convolution import Trapezoid1DKernel 477 trapezoid_1D_kernel = Trapezoid1DKernel(17, slope=0.2) 478 plt.plot(trapezoid_1D_kernel, drawstyle='steps') 479 plt.xlabel('x [pixels]') 480 plt.ylabel('amplitude') 481 plt.xlim(-1, 28) 482 plt.show() 483 """ 484 _is_bool = False 485 486 def __init__(self, width, slope=1., **kwargs): 487 self._model = models.Trapezoid1D(1, 0, width, slope) 488 self._default_size = _round_up_to_odd_integer(width + 2. / slope) 489 super().__init__(**kwargs) 490 self._truncation = 0 491 self.normalize() 492 493 494class TrapezoidDisk2DKernel(Kernel2D): 495 """ 496 2D trapezoid kernel. 497 498 The generated kernel is normalized so that it integrates to 1. 499 500 Parameters 501 ---------- 502 radius : number 503 Width of the filter kernel, defined as the width of the constant part, 504 before it begins to slope down. 505 slope : number 506 Slope of the filter kernel's tails 507 mode : str, optional 508 One of the following discretization modes: 509 * 'center' (default) 510 Discretize model by taking the value 511 at the center of the bin. 512 * 'linear_interp' 513 Discretize model by performing a bilinear interpolation 514 between the values at the corners of the bin. 515 * 'oversample' 516 Discretize model by taking the average 517 on an oversampled grid. 518 * 'integrate' 519 Discretize model by integrating the 520 model over the bin. 521 factor : number, optional 522 Factor of oversampling. Default factor = 10. 523 524 See Also 525 -------- 526 Gaussian2DKernel, Box2DKernel, Tophat2DKernel, RickerWavelet2DKernel, 527 Ring2DKernel, AiryDisk2DKernel, Moffat2DKernel 528 529 Examples 530 -------- 531 Kernel response: 532 533 .. plot:: 534 :include-source: 535 536 import matplotlib.pyplot as plt 537 from astropy.convolution import TrapezoidDisk2DKernel 538 trapezoid_2D_kernel = TrapezoidDisk2DKernel(20, slope=0.2) 539 plt.imshow(trapezoid_2D_kernel, interpolation='none', origin='lower') 540 plt.xlabel('x [pixels]') 541 plt.ylabel('y [pixels]') 542 plt.colorbar() 543 plt.show() 544 545 """ 546 _is_bool = False 547 548 def __init__(self, radius, slope=1., **kwargs): 549 self._model = models.TrapezoidDisk2D(1, 0, 0, radius, slope) 550 self._default_size = _round_up_to_odd_integer(2 * radius + 2. / slope) 551 super().__init__(**kwargs) 552 self._truncation = 0 553 self.normalize() 554 555 556class RickerWavelet1DKernel(Kernel1D): 557 """ 558 1D Ricker wavelet filter kernel (sometimes known as a "Mexican Hat" 559 kernel). 560 561 The Ricker wavelet, or inverted Gaussian-Laplace filter, is a 562 bandpass filter. It smooths the data and removes slowly varying 563 or constant structures (e.g. Background). It is useful for peak or 564 multi-scale detection. 565 566 This kernel is derived from a normalized Gaussian function, by 567 computing the second derivative. This results in an amplitude 568 at the kernels center of 1. / (sqrt(2 * pi) * width ** 3). The 569 normalization is the same as for `scipy.ndimage.gaussian_laplace`, 570 except for a minus sign. 571 572 .. note:: 573 574 See https://github.com/astropy/astropy/pull/9445 for discussions 575 related to renaming of this kernel. 576 577 Parameters 578 ---------- 579 width : number 580 Width of the filter kernel, defined as the standard deviation 581 of the Gaussian function from which it is derived. 582 x_size : int, optional 583 Size in x direction of the kernel array. Default = ⌊8*width +1⌋. 584 mode : str, optional 585 One of the following discretization modes: 586 * 'center' (default) 587 Discretize model by taking the value 588 at the center of the bin. 589 * 'linear_interp' 590 Discretize model by linearly interpolating 591 between the values at the corners of the bin. 592 * 'oversample' 593 Discretize model by taking the average 594 on an oversampled grid. 595 * 'integrate' 596 Discretize model by integrating the 597 model over the bin. 598 factor : number, optional 599 Factor of oversampling. Default factor = 10. 600 601 602 See Also 603 -------- 604 Box1DKernel, Gaussian1DKernel, Trapezoid1DKernel 605 606 Examples 607 -------- 608 Kernel response: 609 610 .. plot:: 611 :include-source: 612 613 import matplotlib.pyplot as plt 614 from astropy.convolution import RickerWavelet1DKernel 615 ricker_1d_kernel = RickerWavelet1DKernel(10) 616 plt.plot(ricker_1d_kernel, drawstyle='steps') 617 plt.xlabel('x [pixels]') 618 plt.ylabel('value') 619 plt.show() 620 621 """ 622 _is_bool = True 623 624 def __init__(self, width, **kwargs): 625 amplitude = 1.0 / (np.sqrt(2 * np.pi) * width ** 3) 626 self._model = models.RickerWavelet1D(amplitude, 0, width) 627 self._default_size = _round_up_to_odd_integer(8 * width) 628 super().__init__(**kwargs) 629 self._truncation = np.abs(self._array.sum() / self._array.size) 630 631 632class RickerWavelet2DKernel(Kernel2D): 633 """ 634 2D Ricker wavelet filter kernel (sometimes known as a "Mexican Hat" 635 kernel). 636 637 The Ricker wavelet, or inverted Gaussian-Laplace filter, is a 638 bandpass filter. It smooths the data and removes slowly varying 639 or constant structures (e.g. Background). It is useful for peak or 640 multi-scale detection. 641 642 This kernel is derived from a normalized Gaussian function, by 643 computing the second derivative. This results in an amplitude 644 at the kernels center of 1. / (pi * width ** 4). The normalization 645 is the same as for `scipy.ndimage.gaussian_laplace`, except 646 for a minus sign. 647 648 .. note:: 649 650 See https://github.com/astropy/astropy/pull/9445 for discussions 651 related to renaming of this kernel. 652 653 Parameters 654 ---------- 655 width : number 656 Width of the filter kernel, defined as the standard deviation 657 of the Gaussian function from which it is derived. 658 x_size : int, optional 659 Size in x direction of the kernel array. Default = ⌊8*width +1⌋. 660 y_size : int, optional 661 Size in y direction of the kernel array. Default = ⌊8*width +1⌋. 662 mode : str, optional 663 One of the following discretization modes: 664 * 'center' (default) 665 Discretize model by taking the value 666 at the center of the bin. 667 * 'linear_interp' 668 Discretize model by performing a bilinear interpolation 669 between the values at the corners of the bin. 670 * 'oversample' 671 Discretize model by taking the average 672 on an oversampled grid. 673 * 'integrate' 674 Discretize model by integrating the 675 model over the bin. 676 factor : number, optional 677 Factor of oversampling. Default factor = 10. 678 679 680 See Also 681 -------- 682 Gaussian2DKernel, Box2DKernel, Tophat2DKernel, Ring2DKernel, 683 TrapezoidDisk2DKernel, AiryDisk2DKernel, Moffat2DKernel 684 685 Examples 686 -------- 687 Kernel response: 688 689 .. plot:: 690 :include-source: 691 692 import matplotlib.pyplot as plt 693 from astropy.convolution import RickerWavelet2DKernel 694 ricker_2d_kernel = RickerWavelet2DKernel(10) 695 plt.imshow(ricker_2d_kernel, interpolation='none', origin='lower') 696 plt.xlabel('x [pixels]') 697 plt.ylabel('y [pixels]') 698 plt.colorbar() 699 plt.show() 700 """ 701 _is_bool = False 702 703 def __init__(self, width, **kwargs): 704 amplitude = 1.0 / (np.pi * width ** 4) 705 self._model = models.RickerWavelet2D(amplitude, 0, 0, width) 706 self._default_size = _round_up_to_odd_integer(8 * width) 707 super().__init__(**kwargs) 708 self._truncation = np.abs(self._array.sum() / self._array.size) 709 710 711class AiryDisk2DKernel(Kernel2D): 712 """ 713 2D Airy disk kernel. 714 715 This kernel models the diffraction pattern of a circular aperture. 716 717 The generated kernel is normalized so that it integrates to 1. 718 719 Parameters 720 ---------- 721 radius : float 722 The radius of the Airy disk kernel (radius of the first zero). 723 x_size : int, optional 724 Size in x direction of the kernel array. Default = ⌊8*radius + 1⌋. 725 y_size : int, optional 726 Size in y direction of the kernel array. Default = ⌊8*radius + 1⌋. 727 mode : str, optional 728 One of the following discretization modes: 729 * 'center' (default) 730 Discretize model by taking the value 731 at the center of the bin. 732 * 'linear_interp' 733 Discretize model by performing a bilinear interpolation 734 between the values at the corners of the bin. 735 * 'oversample' 736 Discretize model by taking the average 737 on an oversampled grid. 738 * 'integrate' 739 Discretize model by integrating the 740 model over the bin. 741 factor : number, optional 742 Factor of oversampling. Default factor = 10. 743 744 See Also 745 -------- 746 Gaussian2DKernel, Box2DKernel, Tophat2DKernel, RickerWavelet2DKernel, 747 Ring2DKernel, TrapezoidDisk2DKernel, Moffat2DKernel 748 749 Examples 750 -------- 751 Kernel response: 752 753 .. plot:: 754 :include-source: 755 756 import matplotlib.pyplot as plt 757 from astropy.convolution import AiryDisk2DKernel 758 airydisk_2D_kernel = AiryDisk2DKernel(10) 759 plt.imshow(airydisk_2D_kernel, interpolation='none', origin='lower') 760 plt.xlabel('x [pixels]') 761 plt.ylabel('y [pixels]') 762 plt.colorbar() 763 plt.show() 764 """ 765 _is_bool = False 766 767 def __init__(self, radius, **kwargs): 768 self._model = models.AiryDisk2D(1, 0, 0, radius) 769 self._default_size = _round_up_to_odd_integer(8 * radius) 770 super().__init__(**kwargs) 771 self.normalize() 772 self._truncation = None 773 774 775class Moffat2DKernel(Kernel2D): 776 """ 777 2D Moffat kernel. 778 779 This kernel is a typical model for a seeing limited PSF. 780 781 The generated kernel is normalized so that it integrates to 1. 782 783 Parameters 784 ---------- 785 gamma : float 786 Core width of the Moffat model. 787 alpha : float 788 Power index of the Moffat model. 789 x_size : int, optional 790 Size in x direction of the kernel array. Default = ⌊8*radius + 1⌋. 791 y_size : int, optional 792 Size in y direction of the kernel array. Default = ⌊8*radius + 1⌋. 793 mode : str, optional 794 One of the following discretization modes: 795 * 'center' (default) 796 Discretize model by taking the value 797 at the center of the bin. 798 * 'linear_interp' 799 Discretize model by performing a bilinear interpolation 800 between the values at the corners of the bin. 801 * 'oversample' 802 Discretize model by taking the average 803 on an oversampled grid. 804 * 'integrate' 805 Discretize model by integrating the 806 model over the bin. 807 factor : number, optional 808 Factor of oversampling. Default factor = 10. 809 810 See Also 811 -------- 812 Gaussian2DKernel, Box2DKernel, Tophat2DKernel, RickerWavelet2DKernel, 813 Ring2DKernel, TrapezoidDisk2DKernel, AiryDisk2DKernel 814 815 Examples 816 -------- 817 Kernel response: 818 819 .. plot:: 820 :include-source: 821 822 import matplotlib.pyplot as plt 823 from astropy.convolution import Moffat2DKernel 824 moffat_2D_kernel = Moffat2DKernel(3, 2) 825 plt.imshow(moffat_2D_kernel, interpolation='none', origin='lower') 826 plt.xlabel('x [pixels]') 827 plt.ylabel('y [pixels]') 828 plt.colorbar() 829 plt.show() 830 """ 831 _is_bool = False 832 833 def __init__(self, gamma, alpha, **kwargs): 834 # Compute amplitude, from 835 # https://en.wikipedia.org/wiki/Moffat_distribution 836 amplitude = (alpha - 1.0) / (np.pi * gamma * gamma) 837 self._model = models.Moffat2D(amplitude, 0, 0, gamma, alpha) 838 self._default_size = _round_up_to_odd_integer(4.0 * self._model.fwhm) 839 super().__init__(**kwargs) 840 self.normalize() 841 self._truncation = None 842 843 844class Model1DKernel(Kernel1D): 845 """ 846 Create kernel from 1D model. 847 848 The model has to be centered on x = 0. 849 850 Parameters 851 ---------- 852 model : `~astropy.modeling.Fittable1DModel` 853 Kernel response function model 854 x_size : int, optional 855 Size in x direction of the kernel array. Default = ⌊8*width +1⌋. 856 Must be odd. 857 mode : str, optional 858 One of the following discretization modes: 859 * 'center' (default) 860 Discretize model by taking the value 861 at the center of the bin. 862 * 'linear_interp' 863 Discretize model by linearly interpolating 864 between the values at the corners of the bin. 865 * 'oversample' 866 Discretize model by taking the average 867 on an oversampled grid. 868 * 'integrate' 869 Discretize model by integrating the 870 model over the bin. 871 factor : number, optional 872 Factor of oversampling. Default factor = 10. 873 874 Raises 875 ------ 876 TypeError 877 If model is not an instance of `~astropy.modeling.Fittable1DModel` 878 879 See also 880 -------- 881 Model2DKernel : Create kernel from `~astropy.modeling.Fittable2DModel` 882 CustomKernel : Create kernel from list or array 883 884 Examples 885 -------- 886 Define a Gaussian1D model: 887 888 >>> from astropy.modeling.models import Gaussian1D 889 >>> from astropy.convolution.kernels import Model1DKernel 890 >>> gauss = Gaussian1D(1, 0, 2) 891 892 And create a custom one dimensional kernel from it: 893 894 >>> gauss_kernel = Model1DKernel(gauss, x_size=9) 895 896 This kernel can now be used like a usual Astropy kernel. 897 """ 898 _separable = False 899 _is_bool = False 900 901 def __init__(self, model, **kwargs): 902 if isinstance(model, Fittable1DModel): 903 self._model = model 904 else: 905 raise TypeError("Must be Fittable1DModel") 906 super().__init__(**kwargs) 907 908 909class Model2DKernel(Kernel2D): 910 """ 911 Create kernel from 2D model. 912 913 The model has to be centered on x = 0 and y = 0. 914 915 Parameters 916 ---------- 917 model : `~astropy.modeling.Fittable2DModel` 918 Kernel response function model 919 x_size : int, optional 920 Size in x direction of the kernel array. Default = ⌊8*width +1⌋. 921 Must be odd. 922 y_size : int, optional 923 Size in y direction of the kernel array. Default = ⌊8*width +1⌋. 924 mode : str, optional 925 One of the following discretization modes: 926 * 'center' (default) 927 Discretize model by taking the value 928 at the center of the bin. 929 * 'linear_interp' 930 Discretize model by performing a bilinear interpolation 931 between the values at the corners of the bin. 932 * 'oversample' 933 Discretize model by taking the average 934 on an oversampled grid. 935 * 'integrate' 936 Discretize model by integrating the 937 model over the bin. 938 factor : number, optional 939 Factor of oversampling. Default factor = 10. 940 941 Raises 942 ------ 943 TypeError 944 If model is not an instance of `~astropy.modeling.Fittable2DModel` 945 946 See also 947 -------- 948 Model1DKernel : Create kernel from `~astropy.modeling.Fittable1DModel` 949 CustomKernel : Create kernel from list or array 950 951 Examples 952 -------- 953 Define a Gaussian2D model: 954 955 >>> from astropy.modeling.models import Gaussian2D 956 >>> from astropy.convolution.kernels import Model2DKernel 957 >>> gauss = Gaussian2D(1, 0, 0, 2, 2) 958 959 And create a custom two dimensional kernel from it: 960 961 >>> gauss_kernel = Model2DKernel(gauss, x_size=9) 962 963 This kernel can now be used like a usual astropy kernel. 964 965 """ 966 _is_bool = False 967 _separable = False 968 969 def __init__(self, model, **kwargs): 970 self._separable = False 971 if isinstance(model, Fittable2DModel): 972 self._model = model 973 else: 974 raise TypeError("Must be Fittable2DModel") 975 super().__init__(**kwargs) 976 977 978class PSFKernel(Kernel2D): 979 """ 980 Initialize filter kernel from astropy PSF instance. 981 """ 982 _separable = False 983 984 def __init__(self): 985 raise NotImplementedError('Not yet implemented') 986 987 988class CustomKernel(Kernel): 989 """ 990 Create filter kernel from list or array. 991 992 Parameters 993 ---------- 994 array : list or array 995 Filter kernel array. Size must be odd. 996 997 Raises 998 ------ 999 TypeError 1000 If array is not a list or array. 1001 `~astropy.convolution.KernelSizeError` 1002 If array size is even. 1003 1004 See also 1005 -------- 1006 Model2DKernel, Model1DKernel 1007 1008 Examples 1009 -------- 1010 Define one dimensional array: 1011 1012 >>> from astropy.convolution.kernels import CustomKernel 1013 >>> import numpy as np 1014 >>> array = np.array([1, 2, 3, 2, 1]) 1015 >>> kernel = CustomKernel(array) 1016 >>> kernel.dimension 1017 1 1018 1019 Define two dimensional array: 1020 1021 >>> array = np.array([[1, 1, 1], [1, 2, 1], [1, 1, 1]]) 1022 >>> kernel = CustomKernel(array) 1023 >>> kernel.dimension 1024 2 1025 """ 1026 def __init__(self, array): 1027 self.array = array 1028 super().__init__(self._array) 1029 1030 @property 1031 def array(self): 1032 """ 1033 Filter kernel array. 1034 """ 1035 return self._array 1036 1037 @array.setter 1038 def array(self, array): 1039 """ 1040 Filter kernel array setter 1041 """ 1042 if isinstance(array, np.ndarray): 1043 self._array = array.astype(np.float64) 1044 elif isinstance(array, list): 1045 self._array = np.array(array, dtype=np.float64) 1046 else: 1047 raise TypeError("Must be list or array.") 1048 1049 # Check if array is odd in all axes 1050 if has_even_axis(self): 1051 raise_even_kernel_exception() 1052 1053 # Check if array is bool 1054 ones = self._array == 1. 1055 zeros = self._array == 0 1056 self._is_bool = bool(np.all(np.logical_or(ones, zeros))) 1057 1058 self._truncation = 0.0 1059 1060 1061@deprecated('4.0', alternative='RickerWavelet1DKernel') 1062class MexicanHat1DKernel(RickerWavelet1DKernel): 1063 pass 1064 1065 1066@deprecated('4.0', alternative='RickerWavelet2DKernel') 1067class MexicanHat2DKernel(RickerWavelet2DKernel): 1068 pass 1069