diff --git a/src/Trajectory.jl b/src/Trajectory.jl index ae2c2aae..ee9e70d8 100644 --- a/src/Trajectory.jl +++ b/src/Trajectory.jl @@ -3,8 +3,7 @@ Trajectory constructor data type. -Defaults to reading with the Chemfiles infrastructure, except for DCD and PDB trajectory -files, if the "PDBTraj" option is provided. +Defaults to reading with the Chemfiles infrastructure, except for DCD and PDB trajectory files, if the "PDBTraj" option is provided. See memory issue (https://github.com/chemfiles/Chemfiles.jl/issues/44) @@ -21,6 +20,9 @@ end stream(traj::Trajectory) = traj.stream.st set_stream!(traj::Trajectory, st) = traj.stream.st = st +# Trajectory formats +const trajectory_formats = String[] + # Include specific predefined trajectory types include("./trajectory_formats/ChemFiles.jl") include("./trajectory_formats/NamdDCD.jl") @@ -33,20 +35,28 @@ function Trajectory( format::String="", chemfiles=false, ) + if !isempty(format) && !(format in trajectory_formats) + throw(ArgumentError("""\n + Trajectory format not properly set: $format + Available trajectory formats: $(join(trajectory_formats, ", ")) + + """)) + end filename = expanduser(filename) # expand tilde on Unix systems, to username - if !chemfiles && (format == "dcd" || split(filename, '.')[end] == "dcd") - trajectory = NamdDCD(filename, solute, solvent) - elseif !chemfiles && format == "PDBTraj" - trajectory = PDBTraj(filename, solute, solvent) + file_extension = filename[findlast(==('.'), filename)+1:end] + trajectory = if !chemfiles && (format == "dcd" || file_extension == "dcd") + NamdDCD(filename, solute, solvent) + elseif !chemfiles && (format == "PDBTraj" || file_extension =="pdb") + PDBTraj(filename, solute, solvent) else - trajectory = ChemFile(filename, solute, solvent, format=format) + ChemFile(filename, solute, solvent, format=format) end return trajectory end # If only one selection is provided, assume that the solute and the solvent are the same Trajectory(filename::String, solvent::AtomSelection; format::String="", chemfiles=false) = - Trajectory(filename, solvent, solvent, format=format, chemfiles=chemfiles) + Trajectory(filename, solvent, solvent; format, chemfiles) #= convert_unitcell(unitcell::Union{SVector{3}, SMatrix{3,3}}) @@ -108,6 +118,9 @@ end protein = AtomSelection(select(atoms, "protein"), nmols=1) tmao = AtomSelection(select(atoms, "resname TMAO"), natomspermol=14) + # Throw if trajectory type is not available + @test_throws ArgumentError Trajectory("$(Testing.data_dir)/PDB/trajectory.pdb", protein, tmao; format="XXX") + # NAMD DCD file traj = Trajectory("$(Testing.data_dir)/NAMD/trajectory.dcd", protein, tmao) @test traj.nframes == 20 @@ -124,6 +137,9 @@ end @test ComplexMixtures.natoms(traj.solute) == 1463 @test ComplexMixtures.natoms(traj.solvent) == 2534 @test ComplexMixtures.convert_unitcell(ComplexMixtures.getunitcell(traj)) ≈ SVector(84.47962951660156, 84.47962951660156, 84.47962951660156) + # Test determining format from file extension + traj = Trajectory("$(Testing.data_dir)/PDB/trajectory.pdb", protein, tmao) + @test ComplexMixtures.natoms(traj.solute) == 1463 # Chemfiles with NAMD traj = Trajectory( diff --git a/src/mddf.jl b/src/mddf.jl index bef36c38..3d4d5c1f 100644 --- a/src/mddf.jl +++ b/src/mddf.jl @@ -183,7 +183,7 @@ function mddf( trajectory_file::String, solute::AtomSelection, solvent::AtomSelection, - options::Options; + options::Options=Options(); trajectory_format::String="", chemfiles::Bool=false, kargs...) @@ -194,7 +194,7 @@ end function mddf( trajectory_file::String, solute_and_solvent::AtomSelection, - options::Options; + options::Options=Options(); trajectory_format::String="", chemfiles::Bool=false, kargs... @@ -528,7 +528,7 @@ function coordination_number( trajectory_file::String, solute::AtomSelection, solvent::AtomSelection, - options::Options; + options::Options=Options(); trajectory_format::String="", chemfiles::Bool=false, kargs...) @@ -540,7 +540,7 @@ end function coordination_number( trajectory_file::String, solute_and_solvent::AtomSelection, - options::Options; + options::Options=Options(); trajectory_format::String="", chemfiles::Bool=false, kargs... @@ -608,9 +608,6 @@ end # Test wrong frame_weights input @test_throws ArgumentError mddf(trajectory_file, protein, water, Options(); frame_weights=[1.0], trajectory_format) - # The coordination_number(::String, args...; kargs...) is a placeholder for the docs only - @test_throws ArgumentError coordination_number("string.txt", 1.0) - # Self correlation atoms = readPDB("$data_dir/toy/self_monoatomic.pdb") atom = AtomSelection(select(atoms, "resname WAT and model 1"), natomspermol=1) @@ -709,6 +706,38 @@ end @test R.density.solvent ≈ 2 / R.volume.total @test R.density.solvent_bulk ≈ 0.5 / R.volume.bulk end +end + +@testitem "mddf - available methods" begin + using ComplexMixtures + using PDBTools: readPDB, select + using ComplexMixtures.Testing: data_dir + atoms = readPDB("$data_dir/toy/cross.pdb") + protein = AtomSelection(select(atoms, "protein and model 1"), nmols=1) + water = AtomSelection(select(atoms, "resname WAT and model 1"), natomspermol=3) + trajectory_file = "$data_dir/toy/cross.pdb" + + atom = AtomSelection(select(atoms, "resname WAT and model 1"), natomspermol=1) + t = Trajectory(trajectory_file, atom; format="PDBTraj") + + options = Options() + trajectory_file = "$data_dir/toy/self_monoatomic.pdb" + atom = AtomSelection(select(atoms, "resname WAT and model 1"), natomspermol=1) + r1 = mddf(trajectory_file, atom, atom) + r2 = mddf(trajectory_file, atom, atom, Options()) + @test coordination_number(r1) ≈ coordination_number(r2) + r1 = mddf(trajectory_file, atom) + r2 = mddf(trajectory_file, atom, Options()) + @test coordination_number(r1) ≈ coordination_number(r2) + r1 = coordination_number(trajectory_file, atom, atom) + r2 = coordination_number(trajectory_file, atom, atom, Options()) + @test coordination_number(r1) ≈ coordination_number(r2) + r1 = coordination_number(trajectory_file, atom) + r2 = coordination_number(trajectory_file, atom, Options()) + @test coordination_number(r1) ≈ coordination_number(r2) + + # The coordination_number(::String, args...; kargs...) is a placeholder for the docs only + @test_throws ArgumentError coordination_number("string.txt", 1.0) end diff --git a/src/trajectory_formats/ChemFiles.jl b/src/trajectory_formats/ChemFiles.jl index 8e9097fb..482ae46c 100644 --- a/src/trajectory_formats/ChemFiles.jl +++ b/src/trajectory_formats/ChemFiles.jl @@ -5,6 +5,9 @@ # import Chemfiles +# Add this format to list of available formats +push!(trajectory_formats, "Chemfiles") + #""" # #$(TYPEDEF) diff --git a/src/trajectory_formats/NamdDCD.jl b/src/trajectory_formats/NamdDCD.jl index 726f648b..555bc243 100644 --- a/src/trajectory_formats/NamdDCD.jl +++ b/src/trajectory_formats/NamdDCD.jl @@ -1,3 +1,6 @@ +# Add this format to list of available formats +push!(trajectory_formats, "dcd") + # # Structure to contain DCD trajectories produces with Namd. # diff --git a/src/trajectory_formats/PDBTraj.jl b/src/trajectory_formats/PDBTraj.jl index d0cfaed9..e44aa7a1 100644 --- a/src/trajectory_formats/PDBTraj.jl +++ b/src/trajectory_formats/PDBTraj.jl @@ -1,3 +1,6 @@ +# Add this format to list of available files +push!(trajectory_formats, "PDBTraj") + #""" # #$(TYPEDEF)