Commit 7fc2b89c5fb24a716546237f60f19d4e259c087e
1 parent
924e3120
Transform tab to spaces.
Showing
1 changed file
with
49 additions
and
49 deletions
src/ctbutterfly/ctbutterfly.cpp
| ... | ... | @@ -19,7 +19,7 @@ |
| 19 | 19 | * * |
| 20 | 20 | ***************************************************************************/ |
| 21 | 21 | /** |
| 22 | - * @file ctbutterfly.hpp | |
| 22 | + * @file ctbutterfly.cpp | |
| 23 | 23 | * @brief Butterfly computation tool |
| 24 | 24 | * @author Michael Mayer |
| 25 | 25 | */ |
| ... | ... | @@ -232,7 +232,7 @@ void ctbutterfly::run(void) |
| 232 | 232 | GModels models = m_obs.models(); |
| 233 | 233 | |
| 234 | 234 | // Compute covariance matrix if necessary |
| 235 | - if (m_covariance.size() == 0) { | |
| 235 | + if (m_covariance.size() == 0) { | |
| 236 | 236 | |
| 237 | 237 | // Write header |
| 238 | 238 | if (logTerse()) { |
| ... | ... | @@ -241,14 +241,14 @@ void ctbutterfly::run(void) |
| 241 | 241 | } |
| 242 | 242 | |
| 243 | 243 | // Evaluate curvature matrix at the actual parameters |
| 244 | - GOptimizerPars pars = models.pars(); | |
| 245 | - GObservations::likelihood likelihood = m_obs.function(); | |
| 246 | - likelihood.eval(pars); | |
| 244 | + GOptimizerPars pars = models.pars(); | |
| 245 | + GObservations::likelihood likelihood = m_obs.function(); | |
| 246 | + likelihood.eval(pars); | |
| 247 | 247 | |
| 248 | 248 | // Get covariance matrix by inverting the curvature matrix |
| 249 | 249 | m_covariance = likelihood.curvature()->invert(); |
| 250 | 250 | |
| 251 | - } | |
| 251 | + } | |
| 252 | 252 | |
| 253 | 253 | // Write header |
| 254 | 254 | if (logTerse()) { |
| ... | ... | @@ -268,73 +268,73 @@ void ctbutterfly::run(void) |
| 268 | 268 | // Initialise number of gradients |
| 269 | 269 | int num_gradient = 0; |
| 270 | 270 | |
| 271 | - // Get the energy of current bin | |
| 271 | + // Get the energy of current bin | |
| 272 | 272 | GEnergy energy = m_ebounds.elogmean(i); |
| 273 | 273 | |
| 274 | - // Header for verbose logging | |
| 275 | - if (logExplicit()) { | |
| 274 | + // Header for verbose logging | |
| 275 | + if (logExplicit()) { | |
| 276 | 276 | std::string msg = "Computing butterfly for bin number "+ |
| 277 | 277 | gammalib::str(i)+" at "+energy.print(); |
| 278 | 278 | log << std::endl; |
| 279 | 279 | log.header2(msg); |
| 280 | - } | |
| 280 | + } | |
| 281 | 281 | |
| 282 | - // Initialise model flux value | |
| 283 | - double model_flux = 0.0; | |
| 282 | + // Initialise model flux value | |
| 283 | + double model_flux = 0.0; | |
| 284 | 284 | |
| 285 | - // Loop over models | |
| 286 | - for (int j = 0; j < models.size(); ++j) { | |
| 285 | + // Loop over models | |
| 286 | + for (int j = 0; j < models.size(); ++j) { | |
| 287 | 287 | |
| 288 | - // Check wether model is a skymodel | |
| 289 | - GModelSky* skymodel = dynamic_cast<GModelSky*>(models[j]); | |
| 288 | + // Check wether model is a skymodel | |
| 289 | + GModelSky* skymodel = dynamic_cast<GModelSky*>(models[j]); | |
| 290 | 290 | |
| 291 | - // Yes ... | |
| 292 | - if (skymodel != NULL) { | |
| 291 | + // Yes ... | |
| 292 | + if (skymodel != NULL) { | |
| 293 | 293 | |
| 294 | - // Skip spatial models | |
| 295 | - num_gradient += skymodel->spatial()->size(); | |
| 294 | + // Skip spatial models | |
| 295 | + num_gradient += skymodel->spatial()->size(); | |
| 296 | 296 | |
| 297 | - // Get pointer to spectral model | |
| 298 | - GModelSpectral* spectral = skymodel->spectral(); | |
| 297 | + // Get pointer to spectral model | |
| 298 | + GModelSpectral* spectral = skymodel->spectral(); | |
| 299 | 299 | |
| 300 | - // Set flux value of source of interest | |
| 301 | - if (skymodel->name() == m_srcname) { | |
| 300 | + // Set flux value of source of interest | |
| 301 | + if (skymodel->name() == m_srcname) { | |
| 302 | 302 | model_flux = spectral->eval_gradients(energy, time); |
| 303 | - } // endif: model was source of interest | |
| 304 | - | |
| 305 | - // Loop over model parameters, get gradients | |
| 306 | - // and assign them to the vector | |
| 307 | - for (int k = 0; k < spectral->size(); ++k) { | |
| 308 | - grad[num_gradient] = (*spectral)[k].gradient(); | |
| 309 | - num_gradient++; | |
| 310 | - } // endfor: looped over spectral parameters | |
| 303 | + } // endif: model was source of interest | |
| 304 | + | |
| 305 | + // Loop over model parameters, get gradients | |
| 306 | + // and assign them to the vector | |
| 307 | + for (int k = 0; k < spectral->size(); ++k) { | |
| 308 | + grad[num_gradient] = (*spectral)[k].gradient(); | |
| 309 | + num_gradient++; | |
| 310 | + } // endfor: looped over spectral parameters | |
| 311 | 311 | |
| 312 | - // Skip temporal models | |
| 313 | - num_gradient += skymodel->temporal()->size(); | |
| 312 | + // Skip temporal models | |
| 313 | + num_gradient += skymodel->temporal()->size(); | |
| 314 | 314 | |
| 315 | - } // endif: model was sky model | |
| 315 | + } // endif: model was sky model | |
| 316 | 316 | |
| 317 | - else { | |
| 317 | + else { | |
| 318 | 318 | |
| 319 | - // Skip other models (e.g. GModelData instances) | |
| 320 | - num_gradient += models[j]->size(); | |
| 319 | + // Skip other models (e.g. GModelData instances) | |
| 320 | + num_gradient += models[j]->size(); | |
| 321 | 321 | |
| 322 | - } // endelse: Model was not GModelSky | |
| 322 | + } // endelse: Model was not GModelSky | |
| 323 | 323 | |
| 324 | - } // endfor: Looped over models | |
| 324 | + } // endfor: Looped over models | |
| 325 | 325 | |
| 326 | - // Multiply covariance to the gradient vector | |
| 327 | - GVector vector = m_covariance * grad; | |
| 326 | + // Multiply covariance to the gradient vector | |
| 327 | + GVector vector = m_covariance * grad; | |
| 328 | 328 | |
| 329 | - // Get the error from the scalar product | |
| 330 | - double error = std::sqrt(grad * vector); | |
| 329 | + // Get the error from the scalar product | |
| 330 | + double error = std::sqrt(grad * vector); | |
| 331 | 331 | |
| 332 | - // Store flux, value and energy for saving | |
| 333 | - m_fluxes.push_back(model_flux); | |
| 334 | - m_energies.push_back(energy.MeV()); | |
| 335 | - m_errors.push_back(error); | |
| 332 | + // Store flux, value and energy for saving | |
| 333 | + m_fluxes.push_back(model_flux); | |
| 334 | + m_energies.push_back(energy.MeV()); | |
| 335 | + m_errors.push_back(error); | |
| 336 | 336 | |
| 337 | - // Log information | |
| 337 | + // Log information | |
| 338 | 338 | if (logExplicit()) { |
| 339 | 339 | log << " Flux at "; |
| 340 | 340 | log << energy.print(); | ... | ... |