Базовый, простой пример на Geant4

Рассмотрим простой пример на Geant4 и разберем его, что бы понять, как можно составлять проекты на Джанте.

Переменные среды Geant4

Сам проект является по сути самостоятельной программой, для которой Geant является библиотекой. Поэтому, при сборке проекта необходимо, чтобы ему были доступны переменные среды от Geant4.

В Ubuntu переменные среды задаются командой

source {папка установки Geant4}/bin/geant4.sh

Пример используемой команды в Ubuntu 16.04

source /home/black/geant4/geant4.10.3-install/bin/geant4.sh

При компилировании проекта вначале нужно установить переменные среды в консоли и только затем можно компилировать. Однако проще, что бы эти переменные были доступны сразу при загрузке компьютера. Для этого можно добавить эту команду в конец файла .bashrc в домашней директории пользователя

~.bashrc

У меня это

/home/black/.bashrc

Так как название файла начинается на точку, что означает, что он является скрытым, поэтому нужно указать в файловом менеджере показ скрытых файлов или использовать файловый менеджер midnight commander из консоли

mc

Сборка с cmake

Проект Geant4 собирается с использованием Cmake, для нашего проекта использовался следующий файл с комментариями на русском

# $Id: CMakeLists.txt 2018-01-08 14:39:15Z Виктор Гавриловец $
#CMakeList.txt, с его помощью собирается проект Geant4, файл основан на примере
#B1 входящим в исходники Geant4
#Как правило здесь настраивать ничего не надо

#----------------------------------------------------------------------------
# Запуск проекта, проверка необходимой версии Cmake
cmake_minimum_required(VERSION 2.6 FATAL_ERROR)
# Название проекта
project(example1)

#----------------------------------------------------------------------------
# Поиск пакетов Geant4, включение всех доступных интерфейсов пользователя
# (Qt, GAG, tcsh, csh) и драйверов визуализации по умолчанию.
# Можно использовать переменную окружения WITH_GEANT4_UIVIS, установить в OFF (выключено)
# через командную строку или с помощью ccmake/cmake-gui для сборки проекта только в
# пакетном режимом
#
option(WITH_GEANT4_UIVIS "Build example with Geant4 UI and Vis drivers" ON)
if(WITH_GEANT4_UIVIS)
  find_package(Geant4 REQUIRED ui_all vis_all)
else()
  find_package(Geant4 REQUIRED)
endif()

#----------------------------------------------------------------------------
# Подключение директорий  из Geant4
# Подключение директорий проекта
#
include(${Geant4_USE_FILE})
include_directories(${PROJECT_SOURCE_DIR}/include)


#----------------------------------------------------------------------------
# Расположение искодного кода и заголовочных файлов этого проекта
# Заголовочные файлы подключаются, значит что они будут показаны в
# IDE - среде разработки
file(GLOB sources ${PROJECT_SOURCE_DIR}/src/*.cpp)
file(GLOB headers ${PROJECT_SOURCE_DIR}/include/*.hh)

# Иногда возникает проблема, что не находит файлы проекта
# это потому, что у файлов исходников могут быть другие расширения, например

# file(GLOB sources ${PROJECT_SOURCE_DIR}/src/*.cc)
# file(GLOB headers ${PROJECT_SOURCE_DIR}/include/*.h)

#----------------------------------------------------------------------------
# Добавление исполняемого файла и линковка его с библиотеками Geant4
#
add_executable(example1 example1.cpp ${sources} ${headers})
target_link_libraries(example1 ${Geant4_LIBRARIES})

#----------------------------------------------------------------------------
# Копирование скриптов в директорию сборки - туда, где мы собираем проект.
# Таким образом мы сможем запустить проект в директории сборки - текущей рабочей
# директории, что нужно например, при отладке.
#
set(EXAMPLE1_SCRIPTS
  run1.mac
  vis.mac
  )

foreach(_script ${EXAMPLE1_SCRIPTS})
  configure_file(
    ${PROJECT_SOURCE_DIR}/${_script}
    ${PROJECT_BINARY_DIR}/${_script}
    COPYONLY
    )
endforeach()

#----------------------------------------------------------------------------
# Для внутеренего использования в Geant4, но не имеет влияния если проект
# компилируется для отдельного использования
#
add_custom_target(ex1 DEPENDS example1)

#----------------------------------------------------------------------------
# Установка исполняемого файла в 'bin' директорию по пути
# установки CMAKE_INSTALL_PREFIX
#
install(TARGETS example1 DESTINATION bin)

С помощью Cmake настраивается сам проект, определяется где что должно лежать и где системе сборки искать файлы.

Основной файл примера

Далее рассмотрим основной файл проекта, содержащий функцию main() и с которого начинается выполнение программы.

// Подключаем заголовочные файлы
#include "G4RunManager.hh" // RunManager, класс из ядра Geant4,
//должен быть включен обязательно
#include "G4UImanager.hh" // Менеджер взаимодействия с пользователем
#include "ExG4DetectorConstruction01.hh" // Конструкция и структура детектора,
//должна определяться пользователем
#include "QBBC.hh" // Подключение физического листа из Geant4, где определяется
// используемые физические процессы и частицы
#include "ExG4ActionInitialization01.hh" // Пользовательский класс
//для задания начального источника частиц
int main()
{
// Создание класса G4RunManager, он контролирует выполнение программы и
// управляет событиями (запуском начальных частиц) при запуске проекта
G4RunManager* runManager = new G4RunManager;
// Установка обязательных инициализирующих классов
// Создание геометрии, материала детектора и мишени
// т.е. определение моделируемой установки
runManager->SetUserInitialization(new ExG4DetectorConstruction01);
// Создание физического листа - набора моделируемых частиц и физических процессов
// которые используются в данном моделировании.
// Используется готовый физический лист из Geant4
runManager->SetUserInitialization(new QBBC);
// Объявление начальных частиц (параметры пучка) и
// подключение прочих классов, используемых
// для получения данных о частицах в процессе моделирования
runManager->SetUserInitialization(new ExG4ActionInitialization01);
// Инициализация ядра Geant4
runManager->Initialize();
// Получение указателя на менеджера взаимодействия с пользователем.
// Указатель нужен, для отправки команд пользователя
// во время выполнения программы
G4UImanager* UI = G4UImanager::GetUIpointer();
// Запуск моделирования, запускаем 3 частицы
runManager->BeamOn(3);
// Окончание работы, вызов деструктора (удаление) G4RunManager
delete runManager;
return 0;
}

Для запуска примера нужно определить как минимум еще три обязательных класса. Заголовочные файлы .hh разместим в папку include, файлы реализации .cpp, разместим в директории src.

Задание геометрии системы

В классе ExG4DetectorConstruction01 опишем геометрию системы

#ifndef ExG4DetectorConstruction01_h
#define ExG4DetectorConstruction01_h 1

#include "G4VUserDetectorConstruction.hh"
#include "globals.hh"

class G4VPhysicalVolume;
class G4LogicalVolume;


/// \brief The ExG4DetectorConstruction01 class Класс геометрии установки,
///  объявление материалов и детекторов
class ExG4DetectorConstruction01 : public G4VUserDetectorConstruction
{
  public:
    //Конструктор, вызывается при создании экземпляра класса
    //Обычно используется для задания начальных значений и значений по умолчанию
    //при создании геометрии и материалов
    ExG4DetectorConstruction01();
    //Деструктор, вызывается при удалении экземпляра класса
    //Обычно используется для освобождения памяти инициализированных массивов внутри класса
    virtual ~ExG4DetectorConstruction01();
    //Объявление и создание детекторов и среды
    virtual G4VPhysicalVolume* Construct();
};
#endif

И файл реализации класса

#include "ExG4DetectorConstruction01.hh"

#include "G4RunManager.hh"
#include "G4NistManager.hh"
#include "G4Box.hh"
#include "G4LogicalVolume.hh"
#include "G4PVPlacement.hh"
#include "G4SystemOfUnits.hh"

// Конструктор класса объявления материалов и геометрии всей моделируемой системы
ExG4DetectorConstruction01::ExG4DetectorConstruction01()
: G4VUserDetectorConstruction()
{ }

// Деструктор
ExG4DetectorConstruction01::~ExG4DetectorConstruction01()
{ }

// Функция определения материалов и геометрии всей системы,
// должна возвращать физический объем - ссылку на экземпляр класса G4VPhysicalVolume
// Геометрию проектировать будем следующую: пучок протонов попадает на мишень
// вольфрамовый лист толщиной 1 мм, а за мишенью поставим детектор
// таких же размеров, он будет регистрировать то, что в него попало.
G4VPhysicalVolume* ExG4DetectorConstruction01::Construct()
{
  // Для простоты используем предопределенные в Geant4 материалы
  // Так объявляется менеджер, из которого можно извлечь
  // ранее предопределенные материалы
  G4NistManager* nist = G4NistManager::Instance();

  // Определяем размеры детектора
  G4double det_sizeXY = 25*cm, det_sizeZ = 0.1*cm;

  // Материал детектора, выбираем из базы данных вольфрам
  G4Material* det_mat = nist->FindOrBuildMaterial("G4_W");

  // Опция для включения/выключения проверки перекрытия объемов
  G4bool checkOverlaps = true;

  //
  // World
  // Объем мира, самый большой объем, включающий остальные, аналог экспериментального
  // зала
  G4double world_sizeXY = 30*cm;//Размер по x и y здесь будут одинаковы - ширина и высота
  G4double world_sizeZ  = 20*cm;//Размер по z - толщина
  // Выбор материала для мира из предопределенных в Geant4, берем воздух
  G4Material* world_mat = nist->FindOrBuildMaterial("G4_AIR");

  // Создание объема для мира установки (экспериментального зала),
  // определяется сама форма объема,
  // берем параллелепипед, это просто геометрическая фигура
  G4Box* solidWorld =
    new G4Box("World",                       //its name, название объема
       0.5*world_sizeXY, 0.5*world_sizeXY, 0.5*world_sizeZ);     //its size, его размеры
  // указываются половины размеров высоты, ширины и глубины

  // Логический объем, здесь подключается материал, из которого сделан объем
  G4LogicalVolume* logicWorld =
    new G4LogicalVolume(solidWorld,          //its solid, геометрический объем, объявлен выше
                        world_mat,           //its material, материал объема
                        "World");            //its name, название логического объема
                                             //совпадает с названием объема, но
                                             //для Geant4 это разные объекты
                                             //геомертический объем и логический объем

  //Физический объем, а теперь наш логический объем помещаем в "ральный" мир
  G4VPhysicalVolume* physWorld =
    new G4PVPlacement(0,                     //no rotation, нет вращения
                      G4ThreeVector(),       //at (0,0,0), расположение в центре (0,0,0)
                      logicWorld,            //its logical volume, логический объем этого физического
                      "World",               //its name, название физического объема
                      0,                     //its mother  volume, материнский объем, этот самый первый, поэтому 0
                      false,                 //no boolean operation, без логических (булевых) операций
                      0,                     //copy number, номер копии
                      checkOverlaps);        //overlaps checking, флаг проверки перекрытия объемов

  // Детектор, для него также используем параллелепипед

  G4Box* solidDet =
    new G4Box("Detector",                    //its name, имя
        0.5*det_sizeXY, 0.5*det_sizeXY, 0.5*det_sizeZ); //its size, размеры

  //Логический объем
  G4LogicalVolume* logicDet =
    new G4LogicalVolume(solidDet,            //its solid, объем
                        det_mat,             //its material, указываем материал детектора
                        "Detector");         //its name, его имя

  //Физический объем детектора
  new G4PVPlacement(0,                       //no rotation, так же без вращения
                    G4ThreeVector(0,0,5*cm), //at (0,0,5 см) положение центра детектора, он смещен на 5 см от центра объема World
                    logicDet,                //its logical volume, подключаем логический объем
                    "Detector",              //its name, имя физического объема
                    logicWorld,              //its mother  volume, родительский логический объем, помещаем в world!
                    false,                   //no boolean operation, без булевых операций
                    0,                       //copy number, номер копии
                    checkOverlaps);          //overlaps checking, флаг проверки перекрытия объемов

  // Для мишени, на которую будет падать пучок, возьмем геометрические размеры как
  // у детектора, параллелепипед - лист вольфрама.
  //Логический объем
  G4LogicalVolume* logicTar =
    new G4LogicalVolume(solidDet,            //its solid, объем
                        det_mat,             //its material, указываем материал мишени
                        "Target");         //its name, его имя

  //Физический объем мишени
  new G4PVPlacement(0,                       //no rotation, так же без вращения
                    G4ThreeVector(0,0,-5*cm),//at (0,0,-5 см) положение центра мишени, смещена на 5 см от центра объема World
                    logicTar,                //its logical volume, подключаем логический объем
                    "Target",              //its name, имя физического объема
                    logicWorld,              //its mother  volume, родительский логический объем, world!
                    false,                   //no boolean operation, без булевых операций
                    0,                       //copy number, номер копии
                    checkOverlaps);
  //Всегда возвращает физический объем
  return physWorld;
}
}

Инициализация пользовательских классов

Класс, который обязательно должен быть реализованн

#ifndef ExG4ActionInitialization01_h
#define ExG4ActionInitialization01_h 1

#include "G4VUserActionInitialization.hh"

/// Обязательный класс, который должен быть объявлен в проекте Geant4
/// Имя класса может быть другим, и он должен наследоваться от
/// класса G4VUserActionInitialization

class ExG4ActionInitialization01 : public G4VUserActionInitialization
{
  public:
    ExG4ActionInitialization01();//Конструктор
    virtual ~ExG4ActionInitialization01();//Деструктор
    virtual void Build() const;//Создание источника первичных частиц

};

#endif

В данном файле регистрируются пользовательские файлы, связанные с заданием источника начальных частиц или пучка. Далее приводится его реализация.

#include "ExG4ActionInitialization01.hh"
#include "ExG4PrimaryGeneratorAction01.hh"//Подключаем обязательный класс
//в котором описываются источник начальных частиц

/// Обязательный класс, который должен быть объявлен в проекте Geant4
/// Имя класса может быть другим, и он должен наследоваться от
/// класса G4VUserActionInitialization
///
/// Конструктор
ExG4ActionInitialization01::ExG4ActionInitialization01()
 : G4VUserActionInitialization()
{}
//Деструктор, ничего не объявляли, поэтому оставим пустым
ExG4ActionInitialization01::~ExG4ActionInitialization01()
{}
//Создание источника первичных частиц
void ExG4ActionInitialization01::Build() const
{
SetUserAction(new ExG4PrimaryGeneratorAction01);//Задается источник первичных частиц
// через обязательный класс ExG4PrimaryGeneratorAction01
}

Задание первичных частиц

Класс описывающий источник начальных частиц

#ifndef ExG4PrimaryGeneratorAction_h
#define ExG4PrimaryGeneratorAction_h 1

#include "G4VUserPrimaryGeneratorAction.hh"
#include "G4ParticleGun.hh"
#include "globals.hh"

class G4ParticleGun;
class G4Event;

/// Класс определения источника первичных частиц
class ExG4PrimaryGeneratorAction01 : public G4VUserPrimaryGeneratorAction
{
  public:
    ///
    /// \brief ExG4PrimaryGeneratorAction01 Конструктор
    ///
    ExG4PrimaryGeneratorAction01();
    ///
    /// \brief ~ExG4PrimaryGeneratorAction01 Деструктор
    ///
    virtual ~ExG4PrimaryGeneratorAction01();
    ///
    /// \brief GeneratePrimaries Метод из базового класса, задает параметры источника начальных частиц
    ///
    virtual void GeneratePrimaries(G4Event*);
    ///
    /// \brief GetParticleGun Метод для доступа к источнику частиц (пушке частиц ;) )
    /// \return Указатель на источник частиц
    ///
    const G4ParticleGun* GetParticleGun() const { return fParticleGun; }

  private:
    G4ParticleGun*  fParticleGun; //указатель на источник частиц
};
#endif

И реализация этого класса

// Подключаем необходимые заголовочные файлы
#include "ExG4PrimaryGeneratorAction01.hh"
// Источник частиц
#include "G4ParticleGun.hh"
// Таблица указателей на данные частиц
#include "G4ParticleTable.hh"
// Данные частиц
#include "G4ParticleDefinition.hh"
// Подключаем систему едениц измерения, что бы пользоваться еденицами cm, MeV.
#include "G4SystemOfUnits.hh"

///
/// \brief ExG4PrimaryGeneratorAction01::ExG4PrimaryGeneratorAction01 Формирование начального пучка.
/// Класс, в котором описывается положение, тип, энергия, направление вылета
/// и распределение начальных частиц
///
ExG4PrimaryGeneratorAction01::ExG4PrimaryGeneratorAction01()
: G4VUserPrimaryGeneratorAction(),
  fParticleGun(0)
{
  // Данные, которые здесь задаются после могут быть изменены через команды в mac файле,
  // или в консольном режиме.
  // По умолчанию у источника частиц поставим 1 частицу
  G4int n_particle = 1;
  fParticleGun  = new G4ParticleGun(n_particle);

  // Получаем встроеную в Geant4 таблицу частиц
  G4ParticleTable* particleTable = G4ParticleTable::GetParticleTable();
  G4String particleName;
  // Ищем частицу, в нашем случае протон
  G4ParticleDefinition* particle
    = particleTable->FindParticle(particleName="proton");
  // Устанавливаем полученную частицу в качестве испускаемого типа начальных частиц в источнике
  fParticleGun->SetParticleDefinition(particle);
  // Устанавливаем направление движение частицы по (x,y,z)
  // Здесь установливается направление вдоль оси Z
  fParticleGun->SetParticleMomentumDirection(G4ThreeVector(0.,0.,1.));
  // Установка начальной энергии испускаемых частиц в 50 МэВ
  fParticleGun->SetParticleEnergy(50*MeV);
}

/// Деструктор, удаляем созданный в конструкторе экземпляр класса источника G4ParticleGun
/// для очистки памяти
ExG4PrimaryGeneratorAction01::~ExG4PrimaryGeneratorAction01()
{
  delete fParticleGun;
}

/// Эта функция вызывается в начале каждого события и генерирует начальные частицы.
/// Определенные здесь условия переписывают команды из файла mac и команды из
/// консольного режима.
void ExG4PrimaryGeneratorAction01::GeneratePrimaries(G4Event* anEvent)
{
  // Задаем координаты где будет источник начальных частиц
  G4double x0 = 0;
  G4double y0 = 0;
  // Устанавливаем в положение -10 см, 0 означает цент мирового объема
  G4double z0 = -0.5 * 20*cm;
  // Команда устанавливает позицию источника начальных частиц
  fParticleGun->SetParticlePosition(G4ThreeVector(x0,y0,z0));
  // Генерация события
  fParticleGun->GeneratePrimaryVertex(anEvent);
}

Когда соберем пример на Geant4, получим в консоли примерно такой вывод

*************************************************************
 Geant4 version Name: geant4-10-03-patch-01    (24-February-2017)
                      Copyright : Geant4 Collaboration
                      Reference : NIM A 506 (2003), 250-303
                            WWW : http://cern.ch/geant4
*************************************************************

<<< Reference Physics List QBBC 
Checking overlaps for volume Detector ... OK! 
Checking overlaps for volume Target ... OK! 
### Adding tracking cuts for neutron  TimeCut(ns)= 10000  KinEnergyCut(MeV)= 0

phot:   for  gamma    SubType= 12  BuildTable= 0
      LambdaPrime table from 200 keV to 100 TeV in 61 bins 
      ===== EM models for the G4Region  DefaultRegionForTheWorld ======
       PhotoElectric :  Emin=        0 eV    Emax=      100 TeV   AngularGenSauterGavrila  FluoActive

compt:   for  gamma    SubType= 13  BuildTable= 1
      Lambda table from 100 eV  to 1 MeV, 7 bins per decade, spline: 1
      LambdaPrime table from 1 MeV to 100 TeV in 56 bins 
      ===== EM models for the G4Region  DefaultRegionForTheWorld ======
       Klein-Nishina :  Emin=        0 eV    Emax=      100 TeV

...................
...................
Еще много строчек и


---------------------------------------------------
                           Hadronic Processes for triton

  Process: hadElastic
        Model:              hElasticLHEP: 0 eV /n ---> 100 TeV/n
     Cr_sctns:            GheishaElastic: 0 eV  ---> 100 TeV

  Process: tInelastic
        Model:  Binary Light Ion Cascade: 0 eV /n ---> 4 GeV/n
        Model:                      FTFP: 2 GeV/n ---> 100 TeV/n
     Cr_sctns: Glauber-Gribov nucleus nucleus: 0 eV  ---> 2.88022e+295 J  
     Cr_sctns:          GheishaInelastic: 0 eV  ---> 100 TeV

================================================================
=======================================================================
======       Pre-compound/De-excitation Physics Parameters     ========
=======================================================================
Type of pre-compound inverse x-section              3
Type of de-excitation inverse x-section             3
Min excitation energy (keV)                         0.01
Level density (1/MeV)                               0.1
Time limit for long lived isomeres (ns)             1000
Use new data files                                  1
Use complete data files                             0
Correlated gamma emission flag                      0
=======================================================================
Number of de-excitation channels                    8

Данные пример только выводит информацию о используемой версии Geant4 и о подключенных процессах, моделях и в принципе показывает, что он запускается. Как увидеть треки частиц, взаимодействия и вывести информацию о них в файл будет идти речь в следующем уроке с полным примером по Geant.

Весь представленный код примера находится здесь.