На примере Example1, который моделирует реальный физический эксперимент, мы изучим как задается геометрия в Geant4. Эксперимент схематично нарисован ниже.
Начнем с того, что изучим файл, содержащий функцию main, который называется Polyethylene.cc. В Geant4 подразумевается определенная структура программы, в которой происходит инициализация объектов в определенной последовательности. Если мы забежим вперед и, скажем, запустим расчет до того, как мы создали геометрию то, само собой, программа вылетит с ошибкой. Если говорить точнее, она даже не соберется.
В начале файла содержится большое количество необходимых инклюдов классов, которые нам будут нужны в дальнейшем. Пока стоит обратить внимание на следующий кусок:
#ifdef G4MULTITHREADED
#include "G4MTRunManager.hh"
#else
#include "G4RunManager.hh"
#endif
И ниже:
#ifdef G4MULTITHREADED
G4MTRunManager* runManager = new G4MTRunManager;
runManager->SetNumberOfThreads(1);
#else
G4RunManager* runManager = new G4RunManager;
#endif
В первом куске содержится директива препроцессору, сообщающая о том, что если мы включили параметр G4MULTITHREADED на стадии cmake, то будет включен заголовок G4MTRunManager.hh (многопоточная версия класса G4RunManager), в другом случае будет включен G4RunManager. G4RunManager - самый главный класс в Geant4, который управляет всем: инициализацией геометрии, подключением физики, продвижением частицы через геометрию и т.д. Его функция Initialize() вызывается в main последней, уже после того как все условия эксперимента заданы.
Во втором куске инициализируется класс G4RunManager или, в многопоточной версии, класс G4MTRunManager с последующим заданием количества потоков.
Теперь мы дошли до этапа инициализации класса, содержащего геометрию.
PolyethyleneDetectorConstruction* massWorld = new PolyethyleneDetectorConstruction;
massWorld->RegisterParallelWorld(new PolyethyleneParallelWorld("PolyethyleneParallelWorld"));
runManager->SetUserInitialization(massWorld);
В первой строчке объект класса с геометрией, собственно, инициализируется. Вторую строчку мы пока пропускаем - это предмет отдельного занятия. В третьей строчке указатель на объект класса с геометрией передается в G4RunManager. Теперь он знает все о нашей геометрии.
Теперь поговорим о том, что же это за класс PolyethyleneDetectorConstruction. Вы, наверно, обратили внимание, что его название явно связано с названием эксперимента. Это не случайно. Класс геометрии - это класс, который мы
обязательно должны написать сами, наследуя класс G4VUserDetectorConstruction, содержащий пустую функцию Construct().
Давайте убедимся в этом сами. Наводим курсор на подсвеченное название класса PolyethyleneDetectorConstruction и нажимаем F2. Мы видим файл заголовка этого класса и вот такую строчку:
class PolyethyleneDetectorConstruction : public G4VUserDetectorConstruction
Эта конструкция и означает, что мы унаследовали класс G4VUserDetectorConstruction. Теперь если мы проведем ту же операцию с названием класса G4VUserDetectorConstruction, то окажемся в его файле заголовка и увидим такую строчку:
virtual G4VPhysicalVolume* Construct() = 0;
Она означает, что материнский класс G4VUserDetectorConstruction содержит пустую функцию Construct(), "= 0" всегда говорит о том, что функцию надо определить в дочернем классе
обязательно.
Настало время задать материалы, которые вместе с геометрией описываются в функции Construct() нашего класса PolyethyleneDetectorConstruction. Для того, чтобы выделить подзадачу задания материалов отдельно (что является хорошим стилем в объектно ориентированном программировании) я создала свою собственную фунцию InitializeMaterials() в нашем классе в разделе private (так как никакому другому классу, кроме конкретно этого, эта функция не понадобится). Проверьте это сами. Обратите внимание, что функция ничего не возвращает, void.
Поговорим о том, как задаются материалы в Geant4. Теоретически, можно создать любые элементы, изотопы и материалы самим, вручную. Но обычно это не нужно, так как Geant4 содержит большую базу данных элементов и материалов. Класс, который обращается к этой базе данных (которая называется NIST) и передает указатели на нужные нам элементы и материалы называется G4NistManager. Это класс static и singleton, что означает что он доступен отовсюду из нашей программы и может быть инициализирован только
один раз. Что означает, что может быть только один объект этого класса в нашей программе, и понятно почему - материал не может быть определен в двух разных местах, неоднозначно.
Итак, если класс G4NistManager static, то это означает, что мы можем прямо сейчас, находясь в функции InitializeMaterials() получить на него указатель. Смотрим в код файла PolyethyleneDetectorConstruction.cc (Напоминаю, что список всех файлов, относящихся к нашей программе, содержится слева в окне QtCreator).
void PolyethyleneDetectorConstruction::InitializeMaterials()
{
G4NistManager* nistManager = G4NistManager::Instance();
Всегда, когда мы будем иметь дело со static файлом, мы будем получать на него указатель такоей же конструкцией через Instance() вместо new. new создает новый объект класса, а, повторюсь, объект этого класса должен быть только один.
Теперь посомотрим на следующие строчки.
G4Element* H = nistManager->FindOrBuildElement(1);
Здесь мы получаем от G4NistManager указатель на необходимый элемент (в данном случае водород) путем вызова его функции FindOrBuildElement, которая принимает как аргумент номер элемента. Дальше идут точно такие же строчки, которые возвращают нам указатели на другие элементы, которые понадобятся нам при создании материалов. Обратите внимание, что мы
не создаем элементы - мы всего лишь получаем на них указатели. Сами элементы уже созданы внутри G4NistManager, второй раз создавать их не имеет смысла.
После того, как мы получили указатели на все необходимые нам элементы, мы можем начать конструировать материалы. Опять обращу внимание, если мы используем материал из базы данных, то нам не надо его создавать через new, в противном случае мы будем использовать new.
Итак, первый пример получения указателя на материал из базы данных полностью аналогичен получению указателя на элемент из базы данных:
G4Material* Air = nistManager->FindOrBuildMaterial("G4_AIR");
Дальше я добавляю его в список материалов, который является private членом нашего класса PolyethyleneDetectorConstruction для того, чтобы мы могли воспользоваться этими материалами из другой функции нашего класса, а именно, из функции Construct().
Посмотрим еще ниже по коду и увидим такой кусок:
G4Material* Brass = new G4Material("Brass", 8.50*g/cm3, 4);
Brass->AddElement(Zn, 0.354);
Brass->AddElement(Cu, 0.6175);
Brass->AddElement(Pb, 0.025);
Brass->AddElement(Fe, 0.0035);
MaterialMap["Brass"] = Brass;
Здесь мы создаем материал латунь своими руками. Сделать так мне пришлось, потому что мне нужен был определенный сплав с определенной массовой долей каждого элемента и плотностью. Первой строчкой мы создаем новый объект класса G4Material через new. Перейдите по F2 в файл заголовка этого класса и просмотрите его функции в разделе public. Повторюсь, что к функциям из раздела private мы не имеем доступа извне. Можно увидеть что функция AddElement, которую мы вызываем впоследствии, принимает как аргумент указатель на элемент и его массовую долю в материале. В файле заголовка содержаться и другие варианты этой функции. Компилятор сам определяет какая функция имеется в виду по аргументам, которые вы в нее передаете.
Теперь перейдем к созданию геометрии и возвращаемся в функцию Construct(). Мы видим, что первой строчкой вызывается функция InitializeMaterials(). Пропускаем следующие две строчки, они касаются визуализации, и о них мы поговорим отдельно.
Теперь мы находимся на участке кода, который создает мировой объем. Это первое и обязательное, что мы дожны сделать при создании геометрии. Важно, чтобы все другие объемы, которые мы поместим внутрь мирового не выходили за его пределы, иначе программа выдаст ошибку.
Первой строчкой мы создаем solid. Это просто форма, которая содержит информацию только о собственном имени "World" и своих линейных размерах. Обратите внимание, что линейные размеры в Geant4 - это полуширины и радиусы. То есть, размер объекта отсчитывается от его центра. Итак, этой строчкой мы создали объект world класса G4Box с именем "World" и размерами 6х6х6 метров. С единицами измерения в Geant4 все просто, они пишутся через звездочку после цифры. Если в начале файла не стоит using namespace CLHEP, то вместо *m мы будем писать *CLHEP::m. CLHEP - файл, а если быть точным namespace, который содаржит определения и взаимосвязи между всевозможными единицами измерения.
Вторая строчка создает логический объем, который содержит информацию о материале, которым он заполнен, своем названии, и главное, содержит указатель на solid. Обратите внимание, что в конструктор класса G4LogicalVolume мы передаем указатель на solid - world. Получается, как будто, вложенная структура. И сейчас мы добавим еще один, третий слой.
G4VPhysicalVolume - класс, который сожержит информацию о расположении логического объема в пространстве, его повороте, материнском логическом объеме (в случае мирового объема - материнского объема нет), своем названии, и еще нескольких других параметрах. Каких - узнайте сами, нажав F2 и перейдя в файл заголовка. И главное, в конструкторе передается указатель на соответствующий ему логический объем.
Обратите внимание, что конструктор G4PVPlacement отличасется от названия класса G4VPhysicalVolume. Это связано с тем, что класс G4PVPlacement - наследник класса G4VPhysicalVolume. Такова особенность языка C++. Для понимания этого нюанса вам будет необходимо продвинуться чуть дальше в понимании принципов объектно ориентированных языков и их качества - полиморфизма. Пока давайте просто запомним, что когда мы создаем физический объем, то используем конструктор G4PVPlacement.
Настало время обратиться к рисунку нашего эксперимента и обратить внимание, что весь объем состоит из повторяющихся структур, а именно тонкие листы латуни чередуются с толстыми слоями полиэтилена. Один лист латуни, заключенный между полуслоями полиэтилена составляет один канал. Поясню физический смысл канала. В данном эксперименте на слоистую структуру падал протонный пучок, и заряд, вследствие остановившихся протонов, накапливался в каналах. Лист латуни является проводником и был подключен к амперметру, который регистрировал ток. Такая конструкция наблюдалась в каждом канале. В итоге, к концу эксперимента у нас имелись данные о распределении заряда по глубине.
Начнем с того, что "нарисуем" один лист латуни с полуразмерами 7,5 х 7,5 х 0,00127 см. Вот как выглядел бы этот кусок кода, если бы нам нужно было нарисовать