Source code for SimEx.Calculators.CrystFELPhotonDiffractor

""":module CrystFELPhotonDiffractor: Module that holds the CrystFELPhotonDiffractor class.  """
##########################################################################
#                                                                        #
# Copyright (C) 2017-2020 Carsten Fortmann-Grote, Juncheng E             #
# Contact: Carsten Fortmann-Grote <carsten.grote@xfel.eu>                #
#                                                                        #
# This file is part of simex_platform.                                   #
# simex_platform is free software: you can redistribute it and/or modify #
# it under the terms of the GNU General Public License as published by   #
# the Free Software Foundation, either version 3 of the License, or      #
# (at your option) any later version.                                    #
#                                                                        #
# simex_platform is distributed in the hope that it will be useful,      #
# but WITHOUT ANY WARRANTY; without even the implied warranty of         #
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the          #
# GNU General Public License for more details.                           #
#                                                                        #
# You should have received a copy of the GNU General Public License      #
# along with this program.  If not, see <http://www.gnu.org/licenses/>.  #
#                                                                        #
##########################################################################
import h5py
import numpy
import os, sys
import subprocess,shlex
import tempfile

from SimEx.Calculators.AbstractPhotonDiffractor import AbstractPhotonDiffractor
from SimEx.Parameters.CrystFELPhotonDiffractorParameters import CrystFELPhotonDiffractorParameters
from SimEx.Parameters.PhotonBeamParameters import propToBeamParameters
from SimEx.Parameters.DetectorGeometry import detectorGeometryFromFile
from SimEx.Utilities import ParallelUtilities
from SimEx.Utilities.EntityChecks import checkAndSetInstance
from SimEx.Utilities.Units import electronvolt, meter
from SimEx.Utilities import IOUtilities

[docs]class CrystFELPhotonDiffractor(AbstractPhotonDiffractor): """ :class CrystFELPhotonDiffractor: Represents simulation of photon diffraction by crystals using CrystFEL.oattern_sim. """ def __init__(self, parameters=None, input_path=None, output_path=None): """ :param parameters: Parameters of the calculation (not data). :type parameters: CrystFELPhotonDiffractorParameters :param input_path: Path to hdf5 file holding the input data. :type input_path: str :param output_path: Path to hdf5 file for output. :type output_path: str """ # Set default parameters if no parameters given. if parameters is None: raise ValueError( "Parameters must at least specify the sample structure, e.g. in pdb format.") else: self.__parameters = checkAndSetInstance( CrystFELPhotonDiffractorParameters, parameters, None ) # Dummy input path since none required (we don't read from a beam propgation yet.) if input_path is None: input_path = os.path.dirname(__file__) else: if parameters.beam_parameters is not None: raise RuntimeError("Beam parameters can only be passed through the parameters object or via the input_path argument.") parameters.beam_parameters = propToBeamParameters( input_path ) if output_path is None: output_path = os.path.join(os.path.dirname(__file__), "diffr") # Init base class. super(CrystFELPhotonDiffractor, self).__init__(parameters,input_path,output_path) ### FIXME self.__provided_data = [ '/data/data', '/data/diffr', '/data/angle', '/history/parent/detail', '/history/parent/parent', '/info/package_version', '/info/contact', '/info/data_description', '/info/method_description', '/params/geom/detectorDist', '/params/geom/pixelWidth', '/params/geom/pixelHeight', '/params/geom/mask', '/params/beam/photonEnergy', '/params/beam/photons', '/params/beam/focusArea', '/params/info', ]
[docs] def expectedData(self): """ Query for the data expected by the Diffractor. """ return None
[docs] def providedData(self): """ Query for the data provided by the Diffractor. """ return None
[docs] def computeNTasks(self): """ Calculate the number of MPI tasks as function of available resources and assigned cpus per task.""" resources=ParallelUtilities.getParallelResourceInfo() ncores=resources['NCores'] nnodes=resources['NNodes'] if nnodes > 1: raise RuntimeError( "Backengine does not support MPI parallelism. " ) if self.parameters.cpus_per_task=="MAX": np=nnodes else: np=max(int(ncores/int(self.parameters.cpus_per_task)),1) return np
# Useless now, just for compatibility def _run_geom(self): """ Perform the actual calls to pattern_sim with multi-panel .geom file. """ # Setup directory structure as needed. if not os.path.isdir( self.output_path ): os.makedirs( self.output_path ) output_file_base = os.path.join( self.output_path, "diffr_out") if self.parameters.number_of_diffraction_patterns == 1: output_file_base += "_0000001.h5" geom_filename = self.parameters.detector_geometry # Setup command, minimum set first. # Distribute patterns over available processes in round-robin. command_sequence = ['pattern_sim', '-p%s' % self.parameters.sample, '--geometry=%s' % geom_filename, '--output=%s' % (output_file_base), '--number=%d' % (self.parameters.number_of_diffraction_patterns) ] # Handle random rotation as requested. if self.parameters.uniform_rotation is True: command_sequence.append('--random-orientation') command_sequence.append('--really-random') if self.parameters.beam_parameters is not None: command_sequence.append('--photon-energy=%f' % (self.parameters.beam_parameters.photon_energy.m_as(electronvolt))) command_sequence.append('--beam-bandwidth=%f' % (self.parameters.beam_parameters.photon_energy_relative_bandwidth)) nphotons = self.parameters.beam_parameters.pulse_energy / self.parameters.beam_parameters.photon_energy command_sequence.append('--nphotons=%e' % (nphotons)) command_sequence.append('--beam-radius=%e' % (self.parameters.beam_parameters.beam_diameter_fwhm.m_as(meter)/2.)) command_sequence.append('--spectrum=%s' % (self.parameters.beam_parameters.photon_energy_spectrum_type.lower())) if self.parameters.beam_parameters.photon_energy_spectrum_type.lower() == "sase": command_sequence.append('--sample-spectrum=512') # Handle intensities list if present. if self.parameters.intensities_file is not None: command_sequence.append('--intensities=%s' % (self.parameters.intensities_file)) # Handle powder if present. if self.parameters.powder is True: command_sequence.append('--powder=%s' % (os.path.join(self.output_path, "powder.h5"))) # Handle size range if present. if self.parameters.crystal_size_min is not None: command_sequence.append('--min-size=%f' % (self.parameters.crystal_size_min.m_as(1e-9*meter) )) if self.parameters.crystal_size_max is not None: command_sequence.append('--max-size=%f' % (self.parameters.crystal_size_max.m_as(1e-9*meter) )) # Handle gpu acceleration. if self.parameters.gpus_per_task > 0: # Check if crystfel was built with opencv support. # Get pattern_sim's path. pattern_sim_path = subprocess.check_output(shlex.split("which pattern_sim"))[:-1].decode('utf-8') # Get list of dynamic dependencies. ldd = subprocess.check_output(shlex.split("ldd "+pattern_sim_path)).decode('utf-8') if "libOpenCL.so.1" in ldd: command_sequence.append('--gpu') if 'SIMEX_VERBOSE' in os.environ: print("Pattern_sim call: "+ " ".join(command_sequence)) # Run the backengine command. proc = subprocess.Popen(command_sequence) proc.wait() return proc.returncode
[docs] def backengine(self): """ This method drives the backengine CrystFEL.pattern_sim.""" # pattern_sim backengine does not support MPI. return self._run() # collect MPI arguments if self.parameters.forced_mpi_command=="": np=self.computeNTasks() mpicommand=ParallelUtilities.prepareMPICommandArguments(np) else: mpicommand=self.parameters.forced_mpi_command # Dump to a temporary file. fname = IOUtilities.getTmpFileName() self.dumpToFile(fname) mpicommand += " ".join([" ",sys.executable, __file__, fname]) args = shlex.split(mpicommand) if 'SIMEX_VERBOSE' in os.environ and os.environ['SIMEX_VERBOSE'] == 'MPI': print("MPI command: "+mpicommand) try: proc = subprocess.Popen(args,universal_newlines=True) proc.wait() except: raise finally: os.remove(fname) return proc.returncode
def _run(self): """ Perform the actual calls to pattern_sim. """ # Setup directory structure as needed. if not os.path.isdir( self.output_path ): os.makedirs( self.output_path ) output_file_base = os.path.join( self.output_path, "diffr_out") if self.parameters.number_of_diffraction_patterns == 1: output_file_base += "_0000001.h5" # Serialize geometry if necessary. if isinstance(self.parameters.detector_geometry, str) and os.path.isfile(self.parameters.detector_geometry): # Convert input .geom file into detector_geometry class input_file = self.parameters.detector_geometry self.parameters.detector_geometry = detectorGeometryFromFile(input_file) geom_file = tempfile.NamedTemporaryFile(suffix=".geom", delete=True) geom_filename = geom_file.name self.parameters.detector_geometry.serialize(stream=geom_filename, caller=self.parameters) # Setup command, minimum set first. # Distribute patterns over available processes in round-robin. command_sequence = ['pattern_sim', '-p%s' % self.parameters.sample, '--geometry=%s' % geom_filename, '--output=%s' % (output_file_base), '--number=%d' % (self.parameters.number_of_diffraction_patterns) ] # Handle random rotation as requested. if self.parameters.uniform_rotation is True: command_sequence.append('--random-orientation') command_sequence.append('--really-random') if self.parameters.beam_parameters is not None: command_sequence.append('--photon-energy=%f' % (self.parameters.beam_parameters.photon_energy.m_as(electronvolt))) command_sequence.append('--beam-bandwidth=%f' % (self.parameters.beam_parameters.photon_energy_relative_bandwidth)) nphotons = self.parameters.beam_parameters.pulse_energy / self.parameters.beam_parameters.photon_energy command_sequence.append('--nphotons=%e' % (nphotons)) command_sequence.append('--beam-radius=%e' % (self.parameters.beam_parameters.beam_diameter_fwhm.m_as(meter)/2.)) command_sequence.append('--spectrum=%s' % (self.parameters.beam_parameters.photon_energy_spectrum_type.lower())) if self.parameters.beam_parameters.photon_energy_spectrum_type.lower() == "sase": command_sequence.append('--sample-spectrum=512') # Handle intensities list if present. if self.parameters.intensities_file is not None: command_sequence.append('--intensities=%s' % (self.parameters.intensities_file)) # Handle powder if present. if self.parameters.powder is True: command_sequence.append('--powder=%s' % (os.path.join(self.output_path, "powder.h5"))) # Handle size range if present. if self.parameters.crystal_size_min is not None: command_sequence.append('--min-size=%f' % (self.parameters.crystal_size_min.m_as(1e-9*meter) )) if self.parameters.crystal_size_max is not None: command_sequence.append('--max-size=%f' % (self.parameters.crystal_size_max.m_as(1e-9*meter) )) # Handle gpu acceleration. if self.parameters.gpus_per_task > 0: # Check if crystfel was built with opencv support. # Get pattern_sim's path. pattern_sim_path = subprocess.check_output(shlex.split("which pattern_sim"))[:-1].decode('utf-8') # Get list of dynamic dependencies. ldd = subprocess.check_output(shlex.split("ldd "+pattern_sim_path)).decode('utf-8') if "libOpenCL.so.1" in ldd: command_sequence.append('--gpu') if 'SIMEX_VERBOSE' in os.environ: print("Pattern_sim call: "+ " ".join(command_sequence)) # Run the backengine command. proc = subprocess.Popen(command_sequence) proc.wait() return proc.returncode @property def data(self): """ Query for the field data. """ return self.__data def _readH5(self): """ """ """ Private method for reading the hdf5 input and extracting the parameters and data relevant to initialize the object. """ pass # Deprecated
[docs] def saveH5_geom(self): """ Method to save the output to a file. Creates links to h5 files that all contain only one pattern. """ # Path where individual h5 files are located. path_to_files = self.output_path # Rename files. _rename_files(path_to_files) print("Linking all patterns into %s.h5." % (self.output_path)) # Setup new file. with h5py.File( self.output_path + ".h5" , "w") as h5_outfile: data_group = h5_outfile.create_group("data") params_group = h5_outfile.create_group("params") beam_params_group = params_group.create_group("beam") geom_params_group = params_group.create_group("geom") beam_params_group["photonEnergy"] = self.parameters.beam_parameters.photon_energy.m_as(electronvolt) beam_params_group["photonEnergy"].attrs["unit_symbol"] = "eV" beam_params_group["photonEnergy"].attrs["unit_longname"] = "electronvolt" beam_params_group["focusArea"] = self.parameters.beam_parameters.beam_diameter_fwhm.m_as(meter)**2 beam_params_group["focusArea"].attrs["unit_symbol"] = "m^2" beam_params_group["focusArea"].attrs["unit_longname"] = "square_metre" # Files to read from. individual_files = [os.path.join( path_to_files, f ) for f in os.listdir( path_to_files ) ] individual_files.sort() # Loop over all individual files and link in the top level groups. for ind_file in individual_files: # Open file. with h5py.File( ind_file, 'r') as h5_infile: # Get file ID. file_ID = os.path.split(ind_file)[-1].split(".h5")[0].split("_")[-1] # Create group data_group.create_group(file_ID) # Links must be relative. relative_link_target = os.path.relpath(path=ind_file, start=os.path.dirname(os.path.dirname(ind_file))) # Link in the data. path_in_target = "/data/data" path_in_origin = "data/%s/data" % (file_ID) h5_outfile[path_in_origin] = h5py.ExternalLink(relative_link_target, path_in_target) # Close input file. h5_infile.close() # Close file. h5_outfile.close() # Reset output_path self.output_path = self.output_path+".h5"
[docs] def saveH5(self): """ Method to save the output to a file. Creates links to h5 files that all contain only one pattern. """ # Path where individual h5 files are located. path_to_files = self.output_path # Rename files. _rename_files(path_to_files) print("Linking all patterns into %s.h5." % (self.output_path)) # Setup new file. with h5py.File( self.output_path + ".h5" , "w") as h5_outfile: data_group = h5_outfile.create_group("data") params_group = h5_outfile.create_group("params") beam_params_group = params_group.create_group("beam") geom_params_group = params_group.create_group("geom") beam_params_group["photonEnergy"] = self.parameters.beam_parameters.photon_energy.m_as(electronvolt) beam_params_group["photonEnergy"].attrs["unit_symbol"] = "eV" beam_params_group["photonEnergy"].attrs["unit_longname"] = "electronvolt" beam_params_group["focusArea"] = self.parameters.beam_parameters.beam_diameter_fwhm.m_as(meter)**2 beam_params_group["focusArea"].attrs["unit_symbol"] = "m^2" beam_params_group["focusArea"].attrs["unit_longname"] = "square_metre" geom_params_group["detectorDist"] = self.parameters.detector_geometry.panels[0].distance_from_interaction_plane.m_as(meter) mask = self.parameters.detector_geometry.panels[0].mask if mask is None: mask = numpy.ones((self.parameters.detector_geometry.panels[0].number_of_pixels_slow, self.parameters.detector_geometry.panels[0].number_of_pixels_fast)) geom_params_group["mask"] = mask geom_params_group["pixelHeight"] = self.parameters.detector_geometry.panels[0].pixel_size.m_as(meter) geom_params_group["pixelWidth"] = self.parameters.detector_geometry.panels[0].pixel_size.m_as(meter) # Files to read from. individual_files = [os.path.join( path_to_files, f ) for f in os.listdir( path_to_files ) ] individual_files.sort() # Loop over all individual files and link in the top level groups. for ind_file in individual_files: # Open file. with h5py.File( ind_file, 'r') as h5_infile: # Get file ID. file_ID = os.path.split(ind_file)[-1].split(".h5")[0].split("_")[-1] # Create group data_group.create_group(file_ID) # Links must be relative. relative_link_target = os.path.relpath(path=ind_file, start=os.path.dirname(os.path.dirname(ind_file))) # Link in the data. path_in_target = "/data/data" path_in_origin = "data/%s/data" % (file_ID) h5_outfile[path_in_origin] = h5py.ExternalLink(relative_link_target, path_in_target) # Close input file. h5_infile.close() # Close file. h5_outfile.close() # Reset output_path self.output_path = self.output_path+".h5"
def _rename_files(path): """ """ """ Renames all files generated by pattern_sim to simex conform filenames.""" old_wd = os.getcwd() os.chdir( path ) original_files = os.listdir(".") original_files.sort() for i,f in enumerate(original_files): if not f.split(".")[-1] == "h5": continue new_filename = "%s_%07d.h5" % ("".join(f.split("-")[:-1]),i+1) print("Renaming %s to %s." % (f, new_filename)) os.rename(f, new_filename) os.chdir( old_wd ) if __name__ == '__main__': CrystFELPhotonDiffractor.runFromCLI()