new HeliMax version

.gitignore

@@ -1,5 +1,9 @@
 __pycache__
 inout_files/
+inout_files_control/
+inout_files_heli_only/
+inout_files_heli_shield/
+inout_files_shield_only/
 *.txt
 *.h5
 *.log

README.md

@@ -1,38 +1,77 @@
 # HeliMax
-This repository aims to show how to generate gprMax input files using NumPy arrays, and how to perform simulation using these input files. The models generated and performed are inspired by an airborne ice Radar setup from ETH Zurich. This work takes place is a master thesis project. The ultimate goal of these simulations is to determine if a shielding between the GPR and the helicopter would be effective and if so, what geometry and properties it should have. The framework can be used to perform any type of layer-based models.
-
-Public description of the project : https://vaw.ethz.ch/en/research/glaciology/research-projects.html
-
-# The setup
-The setup is composed of multiple python files, which the most important one is `main.py`. The folder `inout_files` stores  input and output files, as well as the figures. This repository should be downloaded and placed into the gprMax folder. For example :
-
-`home/user/gprMax/HeliMax`
-
- You must have a gprMax environment correctly installed : https://github.com/gprMax/gprMax
- You must install additionnal requirments (´conda requirment.txt´)
-
-# How it works
-The script is executed using the terminal. `python main.py` generates the model, you can check your geometry and the files generated. `python main.py --run` regenerate all the files and launch the gprMax simulation.
-
-## Generate models
-The models can be generated using simple NumPy arrays. Any kind of geometry can be generated, on few conditions:
-1. The array is composed of integer (not physical properties)
-2. Every integer correspond to a specific material
-3. All materials are described in a `.txt` file (see the file for the structure)
-4. The order of the materials matters
-5. The matrix of integer is stored in a `.h5` file
-6. The `.in` file contains all the information, plus the `.txt` and `.h5` file paths
-
-## gprMax simulations
-The documentation of gprMax give some simple commands to run simulations. For example:
-
-`python -m gprMax user_models/cylinder_Ascan_2D.in`, 
-
-where `cylinder_Ascan.in` is a basic input file from gprMax. In certain cases, it's an advantage to use a GPU to run gprMax simulations. This can be done by using the `-gpu` argument.
-
-# Multi-GPU clusters
-For multi-GPU users, the installation requires some attention. In the same conda environment used for gprMax, some packages must be installed. The three main packages to install are `pycuda`, `openmpy` and `mpi4py`. This must be done using:
-
-`conda install -c conda-forge pycuda mpi4py openmpi`
-
-Then, in certain cases, it will be necessary to remove tha package `mpich` which is installed by conda when installing the package `mpi4py`. `pip` is not used as it will struggle to build wheels to install the packages. 
+
+This repository contains a Python program designed to run gprMax simulations specifically tailored for glaciology research. It utilizes ground-penetrating radar (GPR) simulation to understand and analyze sub-glacial features, including bedrock roughness and ice thickness. The program supports model creation, simulation execution, file merging, and result plotting, with the capability to specify GPU usage for computational efficiency.
+
+## Features
+
+- **Model Generation**: Dynamically generates simulation models based on user-defined parameters.
+- **Simulation Execution**: Runs simulations using specified GPUs to leverage computational resources effectively.
+- **File Merging**: Merges simulation output files for consolidated analysis.
+- **Result Plotting**: Plots the simulation results for visual analysis.
+
+## Requirements
+
+- Python 3.x
+- gprMax
+- tqdm
+- multiprocessing
+
+Please ensure you have gprMax installed and configured correctly on your system. You can find installation instructions for gprMax at [gprMax's GitHub page](https://github.com/gprMax/gprMax).
+
+## Installation
+
+Clone this repository to your local machine, make sure that it's in you folder `gprMax`:
+
+```bash
+git clone https://github.com/bart-inho/helimax.git
+cd helimax
+```
+
+Ensure you have the necessary Python dependencies installed:
+
+```bash
+pip install -r requirements.txt
+```
+
+## Usage
+
+The program is executed via the command line, offering several options to control its operation:
+
+- `--run`: Execute the simulation.
+- `--model`: Generate the simulation model.
+- `--merge`: Merge simulation output files.
+- `--plot`: Plot the simulation results.
+- `--name`: Specify the path for input and output files. This argument is required.
+- `--gpus`: Specify comma-separated GPU IDs to use for the simulation (e.g., `0,1,2,3`).
+
+### Examples
+
+Generate a model:
+
+```bash
+python main.py --model --name path/to/files
+```
+
+Run simulations using GPUs 0 and 1:
+
+```bash
+python main.py --run --name path/to/files --gpus 0,1
+```
+
+Merge simulation files and plot results:
+
+```bash
+python main.py --merge --plot --name path/to/files
+```
+
+## Customization
+
+You can customize the simulation parameters such as antenna spacing, measurement number, and dimensions directly within the script. For more advanced configurations, edit the model generation logic as needed.
+
+## Contributing
+
+Contributions to improve the program or address issues are welcome. Please feel free to submit a pull request or open an issue.
+
+## Acknowledgments
+
+This project utilizes [gprMax](https://github.com/gprMax/gprMax), an open-source project for simulating GPR. We acknowledge the developers and contributors of gprMax for their valuable work.

main.py

@@ -1,91 +1,87 @@
 # This program is used to run the gprMax simulation applied to glaciology.
-from services.file_service import FileService
 from services.folder_init import InitializeFolders
-from models.simulation_model import SimulationModel
-from models.material import Material
+from models.model_generation_logic import ModelGenerationLogic
 from simulations.simulation_runner import SimulationRunner
-from simulations.simulation_plot_profile import PlotProfile  
+from simulations.simulation_plot_profile import PlotProfile
 import argparse
+from tqdm import tqdm  # Import tqdm
+from concurrent.futures import ThreadPoolExecutor
+from multiprocessing import Pool
 
 def main():
     parser = argparse.ArgumentParser()
     parser.add_argument('--run', action='store_true', help='run the simulation')
+    parser.add_argument('--model', action='store_true', help='create the model')
+    parser.add_argument('--merge', action='store_true', help='merge the simulation files')
     parser.add_argument('--plot', action='store_true', help='plot the simulation')
     parser.add_argument('--rough', action='store_true', help='rough bedrock')
+    parser.add_argument('--name', type=str, help='Path to input and output files', required=True)
+    parser.add_argument('--gpus', type=str, help='Comma-separated list of GPU IDs to use', default='0')
+
     args = parser.parse_args()
 
     # Initialize folders
-    InitializeFolders.check_and_create_directories()
-
-    # Initialize Materials
-    # Change the material names in the "Material" class 
-    freespace = Material(1. , 0.   , 1., 0., 'freespace') # Free space
-    glacier   = Material(3.2, 5.e-8, 1., 0., 'glacier'  ) # Glacier
-    bedrock   = Material(5. , 1.e-2, 1., 0., 'bedrock'  ) # Bedrock
-    metal     = Material(1. , 'inf', 1., 0., 'metal'    ) # Helico
+    model_name    = 'model'
+    path_to_files = args.name
+    InitializeFolders.check_and_create_directories(path_to_files)
 
-    # Initialize SimulationModel
-    model_name    = 'test_roughness'
-    inout_files   = 'inout_files/'
-    path_to_files = inout_files + model_name
+    antenna_spacing = 4  # Change antenna spacing in [m] here
 
-    # Generate model
-    model = SimulationModel(model_name, 
-                            45, 10, 100, 
-                            [0.08, 0.08, 0.08], # Change discretisation if needed here
-                            [freespace, glacier, bedrock, metal], # Change name of materials here
-                            inout_files)
-
-    # Generate base model
-    model.generate_base_glacier()
-    model.generate_curved_bedrock_glacier([-10, 5, -200], # center of the curvature [m]
-                             100,            # radius of the curvature [m]
-                             args.rough)
+    dis = 0.06
+    measurement_number = 200 # number of traces
 
-    measurement_number = 24 # number of traces
-    antenna_spacing    = 4  # Change antenna spacing in [m] here
-    measurement_step   = model.calculate_measurment_step(measurement_number, 
-                                                         antenna_spacing) # Change antenna spacing in m here
-
-    # Add antenna positions
-    transceiver1 = [round(25 * model.discrete[0]), # 25 cells of buffer (20 minimum)    
-                    round(model.y_size/2        ),
-                    round(model.z_size/5-.5     )]
-
-    receiver1    = [round(25 * model.discrete[0] + antenna_spacing), # 25 cells of buffer (20 minimum)
-                    round(model.y_size/2                          ),
-                    round(model.z_size/5-.5                       )]
-
-    #Plot initial model
-    model.plot_initial_model(transceiver1, receiver1)
+    x_m = 60
+    y_m = 10
+    z_m = 100
 
-    # Call FileService to write files
-    FileService.write_materials_file(model.path + model.name + '_materials', 
-                                     model.materials)
-
-    FileService.write_h5_file(model.path + model.name + '_h5', 
-                              model)
+    measurement_step = (x_m-7.5) / measurement_number
+    print(measurement_step)
 
-    FileService.write_input_file(model, 
-                                path_to_files, 
-                                path_to_files + '_materials', 
-                                path_to_files + '_h5', 
-                                25e6,   # Change frequency in Hz here
-                                transceiver1, receiver1, 
-                                measurement_step, 
-                                1000e-9) # Change time window in s here
-
-    # Run simulation
+    antenna_x = round(30 * dis) # 30 cells of buffer (20 minimum
+    antenna_y = round(y_m/2)
+    antenna_z = round(20)
+
+    rope_length = 15 # length of the rope in [m]
+
+    if args.model:
+        with ThreadPoolExecutor(max_workers=18) as executor:
+            import concurrent.futures
+            futures = [executor.submit(ModelGenerationLogic.model_generation_logic, idx, 
+                                        model_name, 
+                                        x_m, y_m, z_m, 
+                                        dis, antenna_spacing, 
+                                        rope_length, 
+                                        path_to_files, 
+                                        antenna_x+(idx*measurement_step), antenna_y, antenna_z,
+                                        measurement_step,
+                                        args) for idx in range(measurement_number)]
+
+            # Process the futures as they complete (to maintain order, if necessary)
+            for future in tqdm(concurrent.futures.as_completed(futures), total=measurement_number, desc='Creating Models'):
+                result = future.result()  # result returned from create_model
+
     if args.run:
-        simulation_runner = SimulationRunner(model)
-        simulation_runner.run_simulation(measurement_number)
-        simulation_runner.merge_files(True)
-
-    # Plot profile
+        gpu_ids = [int(gpu) for gpu in args.gpus.split(',')]
+
+        # Prepare arguments for each simulation
+        total_gpus = len(gpu_ids)
+        simulation_args = [(args.name, 'model', idx, gpu_ids[idx % total_gpus]) for idx in range(200)]
+
+        # Use multiprocessing Pool to run simulations in parallel
+        with Pool(processes=total_gpus) as pool:
+            list(tqdm(pool.imap(SimulationRunner.run_simulation, simulation_args), total=len(simulation_args)))
+
+        print("All simulations completed.")
+
+    if args.merge:
+        # Merge the files
+        SimulationRunner.merge_files(False, path_to_files, model_name)
+
     if args.plot:
-        plot_profile = PlotProfile(model.path + model.name + '_merged.out', 'Ey')
+        # Plot the simulation
+        plot_profile = PlotProfile(path_to_files + model_name + '_merged.out', 'Ey')
         plot_profile.get_output_data()
         plot_profile.plot()
-
+        
 if __name__ == "__main__":
     main()

models/model_generation_logic.py

@@ -0,0 +1,83 @@
+from services.file_service import FileService
+from models.simulation_model import SimulationModel
+from models.material import Material
+
+class ModelGenerationLogic:
+    def model_generation_logic(idx, 
+                    model_name, 
+                    x_m, y_m, z_m, dis, 
+                    antenna_spacing, 
+                    rope_length, 
+                    path_to_files, 
+                    antenna_x, antenna_y, antenna_z,
+                    measurement_step,
+                    args):
+
+        # Initialize Materials
+        # Change the material names in the "Material" class 
+        freespace = Material(1. , 0.   , 1., 0., 'freespace') # Free space
+        glacier   = Material(3.2, 5.e-8, 1., 0., 'glacier'  ) # Glacier
+        bedrock   = Material(5. , 1.e-2, 1., 0., 'bedrock'  ) # Bedrock
+        metal     = Material(1. , 'inf', 1., 0., 'metal'    ) # Helico
+        shield    = Material(6. , 1.e-2, 1., 0, 'shield'    ) # Shield
+
+        # Initialize SimulationModel
+        if idx == 0:
+            model_name_it = model_name
+        elif idx > 0:
+            model_name_it = model_name + str(idx)
+
+        # Generate model
+        model = SimulationModel(model_name_it, 
+                                x_m, y_m, z_m, 
+                                [dis, dis, dis], # Change discretisation if needed here
+                                [freespace, glacier, bedrock, metal, shield], # Change name of materials here
+                                path_to_files)
+
+        # Generate base model
+        model.generate_base_glacier()
+        model.generate_curved_bedrock_glacier([-10, 5, -200], # center of the curvature [m]
+                                            100,            # radius of the curvature [m]
+                                            args.rough)
+
+        model.add_CI_shield(antenna_x, antenna_y, antenna_z, dis, antenna_spacing, dis)
+        # model.add_DMC_Tx(antenna_x, antenna_y, antenna_z, dis, antenna_spacing, dis)
+        # model.add_DMC_Rx(antenna_x, antenna_y, antenna_z, dis, antenna_spacing, dis)
+
+        # model.add_3D_oval_shape([antenna_x , antenna_y, antenna_z - rope_length], # center of the curvature [m]
+        #                         [10, .5, .5],            # radius of the curvature [m]])
+        #                         dis*20)
+        # model.add_3D_oval_shape([antenna_x + 3.5, antenna_y, antenna_z - rope_length - 1.5], # center of the curvature [m]
+        #                         [20, .1, .1],            # radius of the curvature [m]])
+        #                         dis*100)
+        # model.add_3D_oval_shape([antenna_x + antenna_spacing - 1, antenna_y, antenna_z - rope_length], # center of the curvature [m]
+        #                         [5, 2, 2],            # radius of the curvature [m]])
+        #                         dis*3)
+
+        transceiver1 = [antenna_x, # 25 cells of buffer (20 minimum)    
+                        antenna_y,
+                        antenna_z]
+
+        receiver1    = [antenna_x + antenna_spacing, # 25 cells of buffer (20 minimum)
+                        antenna_y,
+                        antenna_z]
+
+        #Plot initial model each 20 steps
+        # if idx%20 == 0:
+            # model.plot_initial_model(transceiver1, receiver1)
+
+        # Call FileService to write files
+        FileService.write_materials_file(model.path + model.name + '_materials', 
+                                        model.materials)
+
+        FileService.write_h5_file(model.path + model.name + '_h5', 
+                                model)
+
+        FileService.write_input_file(model, 
+                                    path_to_files+model_name_it, 
+                                    path_to_files+model_name_it + '_materials', 
+                                    path_to_files+model_name_it + '_h5', 
+                                    25e6,   # Change frequency in Hz here
+                                    transceiver1, receiver1, 
+                                    measurement_step, 
+                                    1000e-9) # Change time window in s here