Universita' di GenovaINFN Sezione di Genova  
AIRWATCH / EUSO Genova

O1_ClearSkyPropagator - source file

// ESAF : Euso Simulation and Analysis Framework
// $Id: O1_ClearSkyPropagator.cc,v 1.19 2005/05/02 16:51:55 moreggia Exp $
// Sylvain Moreggia created Jun,  2 2004

#include "O1_ClearSkyPropagator.hh"
#include "EarthVector.hh"
#include "ListPhotonsInAtmosphere.hh"
#include "Atmosphere.hh"
#include "BunchOfPhotons.hh"
#include "RadiativeProcessesCalculator.hh"

ClassImp(O1_ClearSkyPropagator)

//_____________________________________________________________________________
 O1_ClearSkyPropagator::O1_ClearSkyPropagator() : ClearSkyPropagator() {
    //
    // ctor (should not be used)
    //
    fOrder = 1;
}

//_____________________________________________________________________________
 O1_ClearSkyPropagator::O1_ClearSkyPropagator(const Ground* g) : ClearSkyPropagator(g) {
    //
    // ctor, copy RadiatvieTransfer ground description
    //
    fOrder = 1;
}

//_____________________________________________________________________________
 O1_ClearSkyPropagator::~O1_ClearSkyPropagator() {
    //
    // dtor
    //
}

//_____________________________________________________________________________
 Medium O1_ClearSkyPropagator::Go(BunchOfPhotons& bunch,ListPhotonsInAtmosphere& dummy) const {
    //
    // Transport a BunchOfPhotons in clear sky conditions
    // One step propagation, bunch spectrum is convoluted with transmission coefficients
    // No SinglePhoton created in this method
    //
    
    #ifdef DEBUG
    Msg(EsafMsg::Debug) << "O1_CS_Propagator::GO()" << MsgDispatch;
    #endif
    
    // initializations
    static UInt_t bunch_id = 0;
    const BunchOfPhotons& b = bunch;
    EarthVector impact;
    Medium rtn = NONE;
        
    // if fluo, not propagated
    if(b.GetType() == Fluo) {
	bunch.SetFate(1);
        return rtn;
    }

    // if first treatment for this bunch,
    // get next impact and next medium encountered
    // get final impact and final medium encountered (for whole propagation) -> allow (or not) ground reflexion
    //TOFIX : valid algorithm if a SINGLE intermediate medium exists
    if(bunch_id != b.GetId()) {
        bunch_id = b.GetId();
	fFinal_medium = GetFinalImpact(b,impact);
	bunch.SetFinalImpact(impact);
	if(b.GetFinalImpact().Z() == HUGE) Msg(EsafMsg::Warning) << "Bunches should always have an "impact"" << MsgDispatch;
        rtn = GetNextImpact(b,impact);
        bunch.SetNextImpact(impact);
	if(impact.Z() == HUGE) Msg(EsafMsg::Warning) << "Bunches should always have an "impact"" << MsgDispatch;
	if((b.GetNextImpact() - b.GetFinalImpact()).Mag() == 0) rtn = fFinal_medium;  // if no intermediate medium, final bunch propagation step
    }
    // otherwise go until final impact
    else {
        impact = b.GetFinalImpact();
	rtn = fFinal_medium;
    }
    
    // if bunch is underground, its simulation is stopped now
    if(b.GetFinalImpact().Z() == -HUGE) {
        bunch.SetFate(2);
	return NONE;
    }

    // if bunch showerpos is within clouds, switch off the present clear sky propagation
    if(rtn == CLOUDY && (impact - b.GetPos()).Mag() < TOLERANCE) {
        return rtn;
    }

    // initializations
    Double_t conv = 1.;
    Int_t nb = b.GetWlSpectrum().GetAxis().GetNbins();
    vector<Double_t*> coeff;
    Double_t* TotalTrans = new Double_t[nb];
    Double_t* RayleighTrans = new Double_t[nb];
    Double_t* OzoneTrans = new Double_t[nb];
    Double_t* AerosolTrans = new Double_t[nb];
    if(!TotalTrans) Msg(EsafMsg::Panic) << "Pb of memory allocation for TotalTrans" << MsgDispatch;
    if(!RayleighTrans) Msg(EsafMsg::Panic) << "Pb of memory allocation for RayleighTrans" << MsgDispatch;
    if(!OzoneTrans) Msg(EsafMsg::Panic) << "Pb of memory allocation for OzoneTrans" << MsgDispatch;
    if(!AerosolTrans) Msg(EsafMsg::Panic) << "Pb of memory allocation for AerosolTrans" << MsgDispatch;
    coeff.push_back(TotalTrans);
    coeff.push_back(RayleighTrans);
    coeff.push_back(OzoneTrans);
    coeff.push_back(AerosolTrans);

    // calculations
    fCalcul->Trans(b,impact,coeff);
    //bunch.ConvoluteTransSpectra(nb,coeff);
    conv = bunch.ConvoluteWl(coeff[0]);
    bunch.Convolute(conv);
    bunch.AddToPosTof(impact - b.GetPos());

    // cleaning of coeff
    for(UInt_t m=0; m<coeff.size(); m++) {
	if(coeff[m]) delete[] coeff[m];
	coeff[m] = 0;
    }
    coeff.clear();

    if(rtn == GROUND) bunch.SetFate(4);
    // in the case impact is out of FoV, rtn set to NONE to stop the propagation
    if(rtn == CLEARSKY) {
        rtn = NONE;
	bunch.SetFate(5);
    }
    
    return rtn;
} 

//_____________________________________________________________________________
 Double_t O1_ClearSkyPropagator::GoToDetector(EsafSpectrum& spec,const EarthVector& pos) const {
    //
    // method used in Reconstruction
    // convoluates the given spectrum by transmission coefficients between pos and EUSO detector
    // spectrum is normalized and the norm. coeff is returned
    //
    
    // initializations
    BunchOfPhotons b(1.,0.,pos,pos,0.,0.,spec,pos,Cerenkov);
    Double_t rtn = 1.;
    Int_t nb = b.GetWlSpectrum().GetAxis().GetNbins();
    vector<Double_t*> coeff;
    Double_t* TotalTrans = new Double_t[nb];
    Double_t* RayleighTrans = new Double_t[nb];
    Double_t* OzoneTrans = new Double_t[nb];
    Double_t* AerosolTrans = new Double_t[nb];
    if(!TotalTrans) Msg(EsafMsg::Panic) << "Pb of memory allocation for TotalTrans" << MsgDispatch;
    if(!RayleighTrans) Msg(EsafMsg::Panic) << "Pb of memory allocation for RayleighTrans" << MsgDispatch;
    if(!OzoneTrans) Msg(EsafMsg::Panic) << "Pb of memory allocation for OzoneTrans" << MsgDispatch;
    if(!AerosolTrans) Msg(EsafMsg::Panic) << "Pb of memory allocation for AerosolTrans" << MsgDispatch;
    coeff.push_back(TotalTrans);
    coeff.push_back(RayleighTrans);
    coeff.push_back(OzoneTrans);
    coeff.push_back(AerosolTrans);
    
    // set the next impact to EUSO position
    b.SetNextImpact(EUSO());

    // calculations
    fCalcul->Trans(b,EUSO(),coeff);
    rtn = b.ConvoluteWl(coeff[0]);      // normalization coeff.
    b.Convolute(rtn);
    b.GetWlSpectrum().Copy(spec);       // spec is filled here

    // cleaning of coeff
    for(UInt_t m=0; m<coeff.size(); m++) {
	if(coeff[m]) delete[] coeff[m];
	coeff[m] = 0;
    }
    coeff.clear();
    
    return rtn;
}

//_____________________________________________________________________________
 void O1_ClearSkyPropagator::Reset() {
    //
    // reset imbricated objects
    //
    fCalcul->Reset();
}

About Us | EUSO Official Website | Web pages created by Roberto Pesce and Alessandro Thea - Last Update Wed Nov 16 16:57:39 2005 Wed Nov 16 16:29:22 2005