关于粒子收集的一些问题 - 核能革新 ChinaNet
热图推荐
    查看: 7735|回复: 2
    打印 上一主题 下一主题

    关于粒子收集的一些问题

    [复制链接]

    39

    主题

    49

    帖子

    152

    积分

    QQ游客

    积分
    152
    跳转到指定楼层
    楼主
    发表于 2015-4-8 09:55:18 | 只看该作者 回帖奖励 |倒序浏览 |阅读模式
    大家好:
        初学编程对次级离子的收集不甚了解,还请高手指点。比如说,在许多的出射粒子中我只想要电子的信息,包括能谱,角分布。我该怎么实现呢,具体的该看哪些例子呢?
    Geant4/geant4.9.0/examples/Novice/ExN07
    用Particle filter这个method
    比较简单的可以在steppingAction里面这么做
    G4Track * theTrack = theStep->GetTrack();

      // 判断是不是gamma光子
      G4ParticleDefinition * particleType = theTrack->GetDefinition();
    if(particleType==G4Gamma::GammaDefinition())
    ……

    如果要判断电子改成 G4Electron::ElectronDefinition()就行了
    Re:关于粒子收集的一些问题


    给你一段我的代码作参考吧  

    void testSteppingAction::UserSteppingAction(const G4Step* theStep)  
    //这里是从G4自带的G4UserSteppingAction派生了一个类,  
    //这个函数就是替换G4UserSteppingAction中的UserSteppingAction(const G4Step* theStep)函数
      
    {  
    G4Track * theTrack = theStep->GetTrack();  

      // 判断是不是gamma光子  
      G4ParticleDefinition * particleType = theTrack->GetDefinition();  
      if(particleType==G4Gamma::GammaDefinition())  
      {  
      G4VPhysicalVolume * thePrePV = theTrack->GetVolume();  
      G4VPhysicalVolume * theNextPV = theTrack->GetNextVolume();  
      G4String thePrePVname = thePrePV->GetName();  
       
      //判断光子是否进入统计位置  
      if(theNextPV!=NULL)  
      {         
        G4String theNextPVname = theNextPV->GetName();  
        //判断光子是否进入探测面  
        if((theNextPVname(0,10)=="Dec_Sphere")&&  
          (thePrePVname(0,10)!="Dec_Sphere"))  //Dec_Sphere是我设定的探测面名字  
        {  
          //G4cout<<"enter the Detector volume"<<G4endl; //for Debug only!  
        G4double KineEnergy=theTrack->GetKineticEnergy();  //获取当前step的光子动能  
        G4ThreeVector GammaPosition=theTrack->GetPosition();  //获取当前step位置  
        G4double KineTst=KineEnergy/MeV;  //将能量单位转换成MeV  
    //-============================================================================================//  

    //统计能谱  
        G4double gammaEnergy = KineEnergy;  
        EnergyID=int(gammaEnergy/keV); //将当前光子归入某个能量道  
        //G4cout<<"Energy is"<< EnergyID <<G4endl; //for Debug only!  
        theRun->CountID[EnergyID]+=1;  //该能量道的光子数+1,   
                                                              //这里调用了RunAction中的数组  
                                                              //CountID[MaxEnergyID]进行数据传递  
        //G4cout<<"Counts of"<< EnergyID <<"is"<< theRun->CountID[EnergyID] <<G4endl;  
                                                                                                              //for Debug only!  
    //-============================================================================================//  
    //统计角分布  
        KineTst=KineTst*CalculateAirCoef(KineTst);  
        //上面的CalculateAirCoef()是自己定义的一个计算不同能量X射线空气质量能量吸收系数的函数  
        //上面的式子就是将能量转换成空气中X射线能量沉积  
        G4double zPosition=GammaPosition.z();  //获取位置的z方向(主束方向)坐标  
        zPosition=zPosition/mm; //坐标换成mm单位  
        G4double cosPhi=zPosition/1000; //计算极角余弦,  
                                                              //这里1000是我设定的球状探测面的内径,外径是1000.1,  
                                                              //单位mm  
        G4int Phi=int(acos(cosPhi)/deg); //计算出极角  
        //G4cout<<"Phi is "<< Phi <<G4endl; //for Debug only!  

        theRun->NumID[Phi]+=1;    //该极角内光子统计数+1  
        theRun->EDeposit[Phi]+=KineTst; //该极角内光子能量沉积+当前track到达探测面的光子能量沉积  
        G4double Area=cos(Phi*deg)-cos((Phi+1)*deg);  //计算每个极角占有的归一化面积  
        theRun->Dose[Phi]=theRun->EDeposit[Phi]/Area;  //将能量沉积按面积归一化,即转化成相对剂量  
        //G4cout <<"\t"<<Phi<<"\t"<<NumID[Phi]<<"\t\t"<<EDeposit[Phi]<<"\t"<<Dose[Phi]<<endl;  
                                                                                                                  //for Debug only!  
      }  
      }  
    }  
    }



    回复

    使用道具 举报

    39

    主题

    49

    帖子

    152

    积分

    QQ游客

    积分
    152
    沙发
     楼主| 发表于 2015-4-8 09:55:29 | 只看该作者
    你好,我看了N07的粒子了,可是我要的信息,找不到 Particle filte  
    上面是你的问题
    下面是N07中particle filter的段代码
    可以只收集光子,电子等等

    G4SDParticleFilter* gammaFilter = new G4SDParticleFilter

    (filterName="gammaFilter",particleName="gamma");
      G4SDParticleFilter* electronFilter = new G4SDParticleFilter

    (filterName="electronFilter",particleName="e-");
      G4SDParticleFilter* positronFilter = new G4SDParticleFilter

    (filterName="positronFilter",particleName="e+");
      G4SDParticleFilter* epFilter = new G4SDParticleFilter

    (filterName="epFilter");
      epFilter->add(particleName="e-");
      epFilter->add(particleName="e+");
    G4VPrimitiveScorer* primitive;
        primitive = new G4PSEnergyDeposit("eDep",j);
        det->RegisterPrimitive(primitive);
        primitive = new G4PSNofSecondary("nGamma",j);
        primitive->SetFilter(gammaFilter);
        det->RegisterPrimitive(primitive);
        primitive = new G4PSNofSecondary("nElectron",j);
        primitive->SetFilter(electronFilter);
        det->RegisterPrimitive(primitive);
        primitive = new G4PSNofSecondary("nPositron",j);
        primitive->SetFilter(positronFilter);
        det->RegisterPrimitive(primitive);
        primitive = new G4PSMinKinEAtGeneration("minEkinGamma",j);
        primitive->SetFilter(gammaFilter);
        det->RegisterPrimitive(primitive);
        primitive = new G4PSMinKinEAtGeneration("minEkinElectron",j);
        primitive->SetFilter(electronFilter);
        det->RegisterPrimitive(primitive);
        primitive = new G4PSMinKinEAtGeneration("minEkinPositron",j);
        primitive->SetFilter(positronFilter);
        det->RegisterPrimitive(primitive);
        primitive = new G4PSTrackLength("trackLength",j);
        primitive->SetFilter(epFilter);
        det->RegisterPrimitive(primitive);
        primitive = new G4PSNofStep("nStep",j);
        primitive->SetFilter(epFilter);
        det->RegisterPrimitive(primitive);
    回复 支持 反对

    使用道具 举报

    39

    主题

    49

    帖子

    152

    积分

    QQ游客

    积分
    152
    板凳
     楼主| 发表于 2015-4-8 09:55:39 | 只看该作者
    theRun->CountID[EnergyID]+=1;       //该能量道的光子数+1,   
                                                              //这里调用了RunAction中的数组
                                                              //CountID[MaxEnergyID]进行数据传递

    theRun没有定义阿,要怎么调用呢,能不能把调用讲讲~~
    具体就是说SteppingAction和RunAction的头文件要怎么写呢
    回复 支持 反对

    使用道具 举报

      关注我们
    • 微信公众号:
    • NuclearNet
    • 扫描二维码加关注

    Powered by Discuz! X3.2 © 2001-2013 Comsenz Inc.

    联系我们|网站声明|中国核网-核能领域第一垂直门户网站

    快速回复 返回顶部 返回列表