SourcePopulation.h 7.07 KB
Newer Older
Andrew Strong's avatar
Andrew Strong committed
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161


 
#include"Distribution.h"
#include"Skymap.h"        //AWS20090107
#include<vector>          //AWS20100817
#include<valarray>        //AWS20100817

class SourcePopulation
{
 public:
  // parameters of population
  // luminosity in photons s-1 in given energy band

  double pi, rtd, dtr,  kpc_to_cm;

  char title[200]; //AWS20140102  was 100

  int luminosity_function_type;// luminosity-function: 1=delta-function, 2=power-law
  double L;  // luminosity if delta-function

  //  double density0,alpha_R,alpha_z;      // parameters of exponential disk
  double density0,alpha_R,beta_R,zscale;    // parameters of source distribution
  double alpha_z;     
 
  double oversample;//  source oversampling to avoid statistical fluctuations   AWS20060118


  double L_min,L_max,alpha_L;           // parameters of power-law luminosity function

  int    spectral_model;                // 1 = 2 breaks fixed,  2 = 2 breaks with dispersion AWS20100824

  double spectrum_g,spectrum_norm_E;    // spectral index and normalization  energy of sources
  double spectrum_g_0 ;//       spectrum with two breaks                           AWS20051219
  double spectrum_br0 ;
  double spectrum_g_1 ;
  double spectrum_br1 ;
  double spectrum_g_2 ;


  double spectrum_g_0_sigma ;//       spectrum with two breaks: Gaussian distribution AWS20100823
  double spectrum_br0_sigma ;
  double spectrum_g_1_sigma ;
  double spectrum_br1_sigma ;
  double spectrum_g_2_sigma ;

  int   n_E_bin;                        //     number of energy bins
  double *E_bin_low,*E_bin_high;        // boundaries of energy bins
  double  E_ref_low, E_ref_high;        // range for which luminosity is defined


  int   n_E;                           // number of energies for skymap spectra
  double *E;                           //           energies for skymap spectra

  double NL(double L);
  double density(double R, double z); // number per kpc^3

  // derived quantities
  long  n_sample;
  double *sample_R, *sample_x, *sample_y, *sample_z;
  double *sample_theta,*sample_d,*sample_l,*sample_b,*sample_flux;
  double *sample_flux100,*sample_flux1000,*sample_flux10000;               //AWS20120615 integral fluxes
  double *sample_flux_spectrum_E_ref_low,*sample_flux_spectrum_E_ref_high; //AWS20111001
  double sample_d_min,sample_d_max,sample_flux_min,sample_flux_max;
  double sample_total_flux,sample_total_L;
  double sample_selected_total_flux, sample_selected_total_flux_sublimit, sample_selected_total_flux_soplimit;
  double sample_selected_sourcecnt , sample_selected_sourcecnt_sublimit, sample_selected_sourcecnt_soplimit;
  double sample_selected_d_min,         sample_selected_d_max,         sample_selected_d_av;
  double sample_selected_d_min_sublimit,sample_selected_d_max_sublimit,sample_selected_d_av_sublimit;
  double sample_selected_d_min_soplimit,sample_selected_d_max_soplimit,sample_selected_d_av_soplimit;
  double  *sample_l_min,*sample_b_min, *sample_l_max,*sample_b_max, *sample_l_av,*sample_b_av;
  

  double **sample_flux_spectra, **sample_binned_flux_spectra; //AWS20100817 
  double *sample_spectrum_g_0,*sample_spectrum_br0,*sample_spectrum_g_1,*sample_spectrum_br1,*sample_spectrum_g_2; //AWS20100823 parameters per source

  double   latprof_l_min, latprof_l_max, latprof_dlat ; // for profiles
  double  longprof_b_min,longprof_b_max,longprof_dlong;

  double long_min1,long_max1; // for logN-logS
  double long_min2,long_max2;
  double  lat_min1, lat_max1;
  double  lat_min2, lat_max2;

  double *longprof_intensity,*longprof_intensity_sublimit ;
  double *longprof_sourcecnt,*longprof_sourcecnt_sublimit ;
  int     longprof_nlong;
  double * latprof_intensity,* latprof_intensity_sublimit ;
  double * latprof_sourcecnt,* latprof_sourcecnt_sublimit ;
  int      latprof_nlat; 
  
  double flux_detection_limit;
  long i_flux_min,i_flux_max;



  double *dlnN_dlnS;
  double *dlnN_dlnS_deg;

  double *dlnN_dlnS_int;
  double *dlnN_dlnS_int_deg;

  double *FS,*FS_int;

  double lnS_min,lnS_max,dlnS;
  double dlnN_dlnS_total_N, dlnN_dlnS_total_S;
  int   n_dlnN_dlnS;

 
  double *dlnN_dlnS_soplimit;  // N(S) above flux limit  AWS20170123
  double *dlnN_dlnS_sublimit;  // N(S) below flux limit  AWS20170123

  double *dlnN_dlnS_int_soplimit;                      //AWS20170123
  double *dlnN_dlnS_int_sublimit;                      //AWS20170123

  double *FS_soplimit,*FS_int_soplimit;                //AWS20170123
  double *FS_sublimit,*FS_int_sublimit;                //AWS20170123

  double R0 ;// Sun-Galactic centre in kpc

// skymaps of source counts N(S)
  Distribution skymap_sourcecnt;
  Distribution skymap_sourcecnt_sublimit;
  Distribution skymap_sourcecnt_soplimit;

// skymaps integrated over energy ranges
  Distribution skymap_intensity;
  Distribution skymap_intensity_sublimit;
  Distribution skymap_intensity_soplimit;

// skymaps for full energy spectrum
  Distribution skymap_intensity_spectrum;
  Distribution skymap_intensity_spectrum_sublimit;
  Distribution skymap_intensity_spectrum_soplimit;

// ---  healpix versions of the same skymaps ---
  int healpix_skymap_order; // order of healpix maps

// skymaps of source counts N(S)
  Skymap  <double> healpix_skymap_sourcecnt;
  Skymap  <double> healpix_skymap_sourcecnt_sublimit;
  Skymap  <double> healpix_skymap_sourcecnt_soplimit;

// skymaps integrated over energy ranges
  Skymap  <double> healpix_skymap_intensity;
  Skymap  <double> healpix_skymap_intensity_sublimit;
  Skymap  <double> healpix_skymap_intensity_soplimit;

// skymaps for full energy spectrum
  Skymap  <double> healpix_skymap_intensity_spectrum;
  Skymap  <double> healpix_skymap_intensity_spectrum_sublimit;
  Skymap  <double> healpix_skymap_intensity_spectrum_soplimit;

// Fermi-LAT source sensitivity maps                       AWS20170112
  Healpix_Map <double> healpix_skymap_Fermi_sensitivity; //AWS20170112
  char* Fermi_sensitivity_file;                          //AWS20170112

  int verbose;

  ofstream *txt_stream;

162
163
164
165
  // control idl plotting
  int IDL_plot_control;                                  //AWS20200623
  int IDL_psym;                                          //AWS20200623

Andrew Strong's avatar
Andrew Strong committed
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
  int init();
  int print(int options=1);
  int gen_sample(unsigned seed); 
  int analyse_sample();
  int gen_sample_EGRET_catalogue(int options);
  int gen_sample_SPI_catalogue  (char *directory, char *filename, int options);
  int gen_sample_Fermi_catalogue(char *directory, char *filename, int options);//AWS20081215
  
  int gen_Fermi_catalogue_from_sample(char *directory, char *filename, int options);//AWS20110829

  int read_sample_MSPs               (char *directory, char *filename, int options);//AWS20120521


  SourcePopulation operator=(SourcePopulation); // copy constructor

  double integratepowbr(double A0,double g0,double g1,double g2, double br0,double br1, double E1,double E2);
  double     shapepowbr(double A0,double g0,double g1,double g2, double br0,double br1, double E);

  double flux_cutoff_powerlaw(double E                , double g, double Ecutoff,double Emin_ref, double Emax_ref, double flux_ref);
  double flux_cutoff_powerlaw(double Emin, double Emax, double g, double Ecutoff,double Emin_ref, double Emax_ref, double flux_ref);

 };