<SYSTEM_TASK:>
Given the following text description, write Python code to implement the functionality described below step by step
<END_TASK>
<USER_TASK:>
Description:
 Фотометрия
Ноутбук с функциями для работы с фотометрией.
Step1:  Калибровки Bell et al. 2003
Калибровки Bell et al. (2003) https
Step2:  $$\log_{10}(M/L)=a_{\lambda} + b_{\lambda}\times Color$$
Step3:  Самосогласованные калибровки из McGaugh 2014
http
Step4:  Насколько я понял, последняя таблица самая точная - возьмем ее
Step5:  Компоненты галактики
Диск
Step6:  Абсолютная величина Солнца в нескольких полосах отсюда http
Step7:  S4G калибровки и фотометрия
В работе с $\rm{S^4G}$ обзором есть две особенности.
Первая из них - это более точная оценка M/L по калибровке из https
Step8:  Хотя как хорошее приближение можно использовать 0.6, отличие будет небольшим.
Вторая, более важная особенность заключается в том, что поверхностные яркости в декомпозициях даны в т.н. AB-magnitudes, которые надо доисправлять. Из письма
<ASSISTANT_TASK:>
Python Code:
from IPython.display import Image
import numpy as np
import math
%pylab
%matplotlib inline
Explanation: Фотометрия
Ноутбук с функциями для работы с фотометрией.
End of explanation
Image('../Bell_2003.png')
Explanation: Калибровки Bell et al. 2003
Калибровки Bell et al. (2003) https://ui.adsabs.harvard.edu/#abs/2003ApJS..149..289B/abstract
End of explanation
def bell_mass_to_light(color, band, color_str):
    '''Отношение масса светимость вычисляется по калибровке из статьи Bell E. 2003 Table7.'''
    coeffs = {'B-V' : {'B' : [-0.942, 1.737], 'V' : [-0.628, 1.305], 'R' : [-0.520, 1.094], 
                       'I' : [-0.399, 0.824], 'J' : [-0.261, 0.433], 'H' : [-0.209, 0.210], 'K' : [-0.206, 0.135]},
             'B-R' : {'B' : [-0.976, 1.111], 'V' : [-0.633, 0.816], 'R' : [-0.523, 0.683], 
                      'I' : [-0.405, 0.518], 'J' : [-0.289, 0.297], 'H' : [-0.262, 0.180], 'K' : [-0.264, 0.138]},
             'u-g' : { 'g': [-0.221, 0.485], 'r': [-0.099, 0.345], 'i': [-0.053, 0.268], 'z': [-0.105, 0.226]},
             'u-r' : { 'g': [-0.390, 0.417], 'r': [-0.223, 0.299], 'i': [-0.151, 0.233], 'z': [-0.178, 0.192]},
             'u-i' : { 'g': [-0.375, 0.359], 'r': [-0.212, 0.257], 'i': [-0.144, 0.201], 'z': [-0.171, 0.165]},
             'u-z' : { 'g': [-0.400, 0.332], 'r': [-0.232, 0.239], 'i': [-0.161, 0.187], 'z': [-0.179, 0.151]},
             'g-r' : { 'g': [-0.499, 1.519], 'r': [-0.306, 1.097], 'i': [-0.222, 0.864], 'z': [-0.223, 0.689]},
             'g-i' : { 'g': [-0.379, 0.914], 'r': [-0.220, 0.661], 'i': [-0.152, 0.518], 'z': [-0.175, 0.421]},
             'g-z' : { 'g': [-0.367, 0.698], 'r': [-0.215, 0.508], 'i': [-0.153, 0.402], 'z': [-0.171, 0.322]},
             'r-i' : { 'g': [-0.106, 1.982], 'r': [-0.022, 1.431], 'i': [0.006, 1.114], 'z': [-0.052, 0.923]},
             'r-z' : { 'g': [-0.124, 1.067], 'r': [-0.041, 0.780], 'i': [-0.018, 0.623], 'z': [-0.041, 0.463]}}
    a = coeffs[color_str][band][0]
    b = coeffs[color_str][band][1]
    return np.power(10, a + b * color)
Explanation: $$\log_{10}(M/L)=a_{\lambda} + b_{\lambda}\times Color$$
End of explanation
Image('../McGaugh_2014_table7.png')
Explanation: Самосогласованные калибровки из McGaugh 2014
http://iopscience.iop.org/article/10.1088/0004-6256/148/5/77/pdf
End of explanation
def mcgaugh_mass_to_light(color, band) :
    '''Отношение масса светимость вычисляется по калибровке из статьи McGaugh 2014 Table 7; возвращется 4 возможных значения. Цвет из B-V'''
    Bell_2003  = {'V' : [-0.628, 1.305], 'I' : [-0.275, 0.612], 'mu36' : [-0.322, -0.007]}
    Portinari_2004 = {'V' : [-0.654, 1.290], 'I' : [-0.321, 0.701], 'mu36' : [-0.594, 0.467]}
    Zibetti_2009 = {'V' : [-1.075, 1.837], 'I' : [-0.477, 1.004], 'mu36' : [-1.147, 1.289]}
    Into_2013 = {'V' : [-0.900, 1.627], 'I' : [-0.421, 0.898], 'mu36' : [-0.861, 0.849]}
    return np.round([np.power(10, Bell_2003[band][0] + Bell_2003[band][1] * color),
            np.power(10, Portinari_2004[band][0] + Portinari_2004[band][1] * color),
            np.power(10, Zibetti_2009[band][0] + Zibetti_2009[band][1] * color),
            np.power(10, Into_2013[band][0] + Into_2013[band][1] * color)], 2)
Explanation: Насколько я понял, последняя таблица самая точная - возьмем ее:
End of explanation
#TODO: добавить ссылку
def mu_bulge(r, mu_eff=None, r_eff=None, n=None):
    '''поверхностная яркость балджа с квадр. угл. секунды
    параметры - эфф. яркость, эфф. радиус, параметр Серсика, расстояние'''
    bn = 1.9992*n - 0.3271
    return mu_eff + 1.0857*bn*((r/r_eff)**(1./n) - 1)
def mu_bulge2(r, mu_eff=None, r_eff=None, n=None):
    '''второй вариант формулы для поверхностной яркости балджа с квадр. угл. секунды
    параметры те же'''
    nu_n = 2*n - 1./3  + 4./(405*n) + 46./(25515*n**2)
    return mu_eff + (2.5*nu_n/np.log(10))*((r/r_eff)**(1./n) - 1)
def mu_disc(r, mu0=None, h=None):
    '''поверхностная яркость для диска'''
    return mu0 + 1.0857*r/h
def total_mu_profile(*args):
    '''Пересчет суммарного профиля яркости в звездных величинах для нескольких компонент, например disc+bulge
    или disc1+ disc2. Важно, что они должны быть из одной полосы.'''
    points = len(args[0])
    sum_brightness = np.zeros(points)
    for arg in args:
        sum_brightness += np.power(10, -np.array(arg)/2.5)
    return -2.5*np.log10(sum_brightness)
Explanation: Компоненты галактики
Диск:
$$I_d(r) = I_{0,d}\exp{-\frac{r}{h}}$$ или $$\mu_d(r) = \mu_{0,d} + 1.0857\frac{r}{h}$$
Балдж:
$$I_b(r) = I_e\exp\large[-b_n\large[{(\frac{r}{r_e})}^{\frac{1}{n}}-1 \large] \large]$$ $I_e = I_{0,b}\exp{b_n}$ $$\mu_b(r) = \mu_{e} + 1.0857b_n\large[{(\frac{r}{r_e})}^{\frac{1}{n}}-1 \large]$$
Или в другом виде:
$$\mu_b(r) = \mu_{0} + \frac{2.5\nu_n}{\log 10}(\frac{r}{r_e})^{\frac{1}{n}}$$ где приближенно $\nu_n = 2n - \frac{1}{3} + \frac{4}{405n} = \frac{46}{25515n^2} + O(n^{-3})$ (из Решетникова http://www.astro.spbu.ru/staff/resh/Books/SurfPhot/node14.html). Отметим, что в формуле выше $\mu_e$, a ниже $\mu_0$ - надо аккуратно за этим следить ($\mu_e = \mu_0 + 2.5\nu_n/\log 10$).
End of explanation
sun_abs_mags = {'U' : 5.61,
                'B' : 5.48,
                'V' : 4.83,
                'R' : 4.42,
                'I' : 4.08,
                'J' : 3.64,
                'H' : 3.32,
                'K' : 3.28,
                '3.6' : 3.24, # Oh et al. 2008
                'u' : 6.77, #SDSS bands from http://mips.as.arizona.edu/~cnaw/sun.html
                'g' : 5.36,
                'r' : 4.67,
                'i' : 4.48,
                'z' : 4.42
               }
def surface_brightness(mu=None, band=None):
    '''Пересчет поверхностной яркости из mag/arcsec^2 в Lsun/pc^2 по формуле
    http://www.astro.spbu.ru/staff/resh/Books/SurfPhot/node3.html#sec1.1'''
    sun_mag = sun_abs_mags[band]
    return 4.255 * math.pow(10, 8 + 0.4 * (sun_mag - mu))
def surf_density(mu=None, M_to_L=None, band=None):
    '''Плотность для любого компонента галактики - пересчитываем mu в яркость и домножаем на M/L, размерность M_sun/pc^2'''
    return M_to_L*surface_brightness(mu=mu, band=band)
def disc_totmag(mu0d, h, scale):
    '''Полная звездная светимость диска, формула из Решетникова, масштаб нужен для перевода в кпк.'''
    return mu0d - 5*log10(h*scale) - 38.57
Explanation: Абсолютная величина Солнца в нескольких полосах отсюда http://www.ucolick.org/~cnaw/sun.html (на самом деле из Binney and Merrifield 1998, Galactic Astronomy, Table 2.1 (page 53)):
End of explanation
def s4g_mass_to_light(color36, color45):
    return np.power(10., -0.339*(color36 - color45) -0.336)
Explanation: S4G калибровки и фотометрия
В работе с $\rm{S^4G}$ обзором есть две особенности.
Первая из них - это более точная оценка M/L по калибровке из https://arxiv.org/abs/1410.0009 для цветов в $3.6\mu m$ и $4.5\mu m$, которые можно брать из https://irsa.ipac.caltech.edu/cgi-bin/Gator/nph-dd?catalog=s4gcat $$log(M/L) = −0.339\times([3.6]−[4.5])−0.336$$
End of explanation
def s4g_surf_density(mu=None, M_to_L=None):
    return M_to_L*np.power(10., 16.76 - 0.4*mu)/1000000.
Explanation: Хотя как хорошее приближение можно использовать 0.6, отличие будет небольшим.
Вторая, более важная особенность заключается в том, что поверхностные яркости в декомпозициях даны в т.н. AB-magnitudes, которые надо доисправлять. Из письма:
I think I know what can be the problem. Sasha, SB and magnitudes in S4G are in AB-values (AB-magnitudes should be converted to Vega magnitudes and then to luminosities, see for WISE here http://wise2.ipac.caltech.edu/docs/release/allsky/expsup/sec4_4h.html#example, Section 3). E.g. in https://arxiv.org/pdf/1304.6083v1.pdf, there is the formula A2, with which you can convert AB-magnitudes (absolute) to solar masses. Check your conversion formula (from SB to surface density) for S4G.  
Похоже, можно использовать формулу (А5) из https://arxiv.org/pdf/1304.6083v1.pdf сразу для исправления поверхностной плотности 
$$\log_{10}\Sigma_∗(M_{\circ}{kpc}^{−2}) = 16.76 − 0.4\mu_{3.6}(AB_{mag} {arcsec}^{−2})$$
End of explanation
<END_TASK>