Page Menu
Home
c4science
Search
Configure Global Search
Log In
Files
F81432797
molfile_interface.cpp
No One
Temporary
Actions
Download File
Edit File
Delete File
View Transforms
Subscribe
Mute Notifications
Award Token
Subscribers
None
File Metadata
Details
File Info
Storage
Attached
Created
Fri, Sep 6, 01:22
Size
31 KB
Mime Type
text/x-c
Expires
Sun, Sep 8, 01:22 (2 d)
Engine
blob
Format
Raw Data
Handle
20536912
Attached To
rLAMMPS lammps
molfile_interface.cpp
View Options
This document is not UTF8. It was detected as Shift JIS and converted to UTF8 for display.
/* -*- c++ -*- ----------------------------------------------------------
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
http://lammps.sandia.gov, Sandia National Laboratories
Steve Plimpton, sjplimp@sandia.gov
Copyright (2003) Sandia Corporation. Under the terms of Contract
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
certain rights in this software. This software is distributed under
the GNU General Public License.
See the README file in the top-level LAMMPS directory.
------------------------------------------------------------------------- */
/* ----------------------------------------------------------------------
Contributing author: Axel Kohlmeyer (Temple)
------------------------------------------------------------------------- */
#include "molfile_interface.h"
#include <sys/types.h>
#include <stdio.h>
#include <dirent.h>
#include <string.h>
#include <stdlib.h>
#include <ctype.h>
#if defined(_WIN32)
#include <windows.h>
#else
#include <dirent.h>
#include <dlfcn.h>
#endif
#include "molfile_plugin.h"
#if vmdplugin_ABIVERSION < 16
#error "unsupported VMD molfile plugin ABI version"
#endif
#define DEBUG 0
extern
"C"
{
typedef
int
(
*
initfunc
)(
void
);
typedef
int
(
*
regfunc
)(
void
*
,
vmdplugin_register_cb
);
typedef
int
(
*
finifunc
)(
void
);
typedef
struct
{
void
*
p
;
const
char
*
name
;
}
plugin_reginfo_t
;
// callback function for plugin registration.
static
int
plugin_register_cb
(
void
*
v
,
vmdplugin_t
*
p
)
{
plugin_reginfo_t
*
r
=
static_cast
<
plugin_reginfo_t
*>
(
v
);
// make sure we have the proper plugin type (native reader)
// for the desired file type (called "name" at this level)
if
((
strcmp
(
MOLFILE_PLUGIN_TYPE
,
p
->
type
)
==
0
)
&&
(
strcmp
(
r
->
name
,
p
->
name
)
==
0
)
)
{
r
->
p
=
static_cast
<
void
*>
(
p
);
}
return
0
;
}
/* periodic table of elements for translation of ordinal to atom type */
static
const
char
*
pte_label
[]
=
{
"X"
,
"H"
,
"He"
,
"Li"
,
"Be"
,
"B"
,
"C"
,
"N"
,
"O"
,
"F"
,
"Ne"
,
"Na"
,
"Mg"
,
"Al"
,
"Si"
,
"P"
,
"S"
,
"Cl"
,
"Ar"
,
"K"
,
"Ca"
,
"Sc"
,
"Ti"
,
"V"
,
"Cr"
,
"Mn"
,
"Fe"
,
"Co"
,
"Ni"
,
"Cu"
,
"Zn"
,
"Ga"
,
"Ge"
,
"As"
,
"Se"
,
"Br"
,
"Kr"
,
"Rb"
,
"Sr"
,
"Y"
,
"Zr"
,
"Nb"
,
"Mo"
,
"Tc"
,
"Ru"
,
"Rh"
,
"Pd"
,
"Ag"
,
"Cd"
,
"In"
,
"Sn"
,
"Sb"
,
"Te"
,
"I"
,
"Xe"
,
"Cs"
,
"Ba"
,
"La"
,
"Ce"
,
"Pr"
,
"Nd"
,
"Pm"
,
"Sm"
,
"Eu"
,
"Gd"
,
"Tb"
,
"Dy"
,
"Ho"
,
"Er"
,
"Tm"
,
"Yb"
,
"Lu"
,
"Hf"
,
"Ta"
,
"W"
,
"Re"
,
"Os"
,
"Ir"
,
"Pt"
,
"Au"
,
"Hg"
,
"Tl"
,
"Pb"
,
"Bi"
,
"Po"
,
"At"
,
"Rn"
,
"Fr"
,
"Ra"
,
"Ac"
,
"Th"
,
"Pa"
,
"U"
,
"Np"
,
"Pu"
,
"Am"
,
"Cm"
,
"Bk"
,
"Cf"
,
"Es"
,
"Fm"
,
"Md"
,
"No"
,
"Lr"
,
"Rf"
,
"Db"
,
"Sg"
,
"Bh"
,
"Hs"
,
"Mt"
,
"Ds"
,
"Rg"
};
static
const
int
nr_pte_entries
=
sizeof
(
pte_label
)
/
sizeof
(
char
*
);
/* corresponding table of masses. */
static
const
float
pte_mass
[]
=
{
/* X */
0.00000
,
1.00794
,
4.00260
,
6.941
,
9.012182
,
10.811
,
/* C */
12.0107
,
14.0067
,
15.9994
,
18.9984032
,
20.1797
,
/* Na */
22.989770
,
24.3050
,
26.981538
,
28.0855
,
30.973761
,
/* S */
32.065
,
35.453
,
39.948
,
39.0983
,
40.078
,
44.955910
,
/* Ti */
47.867
,
50.9415
,
51.9961
,
54.938049
,
55.845
,
58.9332
,
/* Ni */
58.6934
,
63.546
,
65.409
,
69.723
,
72.64
,
74.92160
,
/* Se */
78.96
,
79.904
,
83.798
,
85.4678
,
87.62
,
88.90585
,
/* Zr */
91.224
,
92.90638
,
95.94
,
98.0
,
101.07
,
102.90550
,
/* Pd */
106.42
,
107.8682
,
112.411
,
114.818
,
118.710
,
121.760
,
/* Te */
127.60
,
126.90447
,
131.293
,
132.90545
,
137.327
,
/* La */
138.9055
,
140.116
,
140.90765
,
144.24
,
145.0
,
150.36
,
/* Eu */
151.964
,
157.25
,
158.92534
,
162.500
,
164.93032
,
/* Er */
167.259
,
168.93421
,
173.04
,
174.967
,
178.49
,
180.9479
,
/* W */
183.84
,
186.207
,
190.23
,
192.217
,
195.078
,
196.96655
,
/* Hg */
200.59
,
204.3833
,
207.2
,
208.98038
,
209.0
,
210.0
,
222.0
,
/* Fr */
223.0
,
226.0
,
227.0
,
232.0381
,
231.03588
,
238.02891
,
/* Np */
237.0
,
244.0
,
243.0
,
247.0
,
247.0
,
251.0
,
252.0
,
257.0
,
/* Md */
258.0
,
259.0
,
262.0
,
261.0
,
262.0
,
266.0
,
264.0
,
269.0
,
/* Mt */
268.0
,
271.0
,
272.0
};
/*
* corresponding table of VDW radii.
* van der Waals radii are taken from A. Bondi,
* J. Phys. Chem., 68, 441 - 452, 1964,
* except the value for H, which is taken from R.S. Rowland & R. Taylor,
* J.Phys.Chem., 100, 7384 - 7391, 1996. Radii that are not available in
* either of these publications have RvdW = 2.00 ナ.
* The radii for Ions (Na, K, Cl, Ca, Mg, and Cs are based on the CHARMM27
* Rmin/2 parameters for (SOD, POT, CLA, CAL, MG, CES) by default.
*/
static
const
float
pte_vdw_radius
[]
=
{
/* X */
1.5
,
1.2
,
1.4
,
1.82
,
2.0
,
2.0
,
/* C */
1.7
,
1.55
,
1.52
,
1.47
,
1.54
,
/* Na */
1.36
,
1.18
,
2.0
,
2.1
,
1.8
,
/* S */
1.8
,
2.27
,
1.88
,
1.76
,
1.37
,
2.0
,
/* Ti */
2.0
,
2.0
,
2.0
,
2.0
,
2.0
,
2.0
,
/* Ni */
1.63
,
1.4
,
1.39
,
1.07
,
2.0
,
1.85
,
/* Se */
1.9
,
1.85
,
2.02
,
2.0
,
2.0
,
2.0
,
/* Zr */
2.0
,
2.0
,
2.0
,
2.0
,
2.0
,
2.0
,
/* Pd */
1.63
,
1.72
,
1.58
,
1.93
,
2.17
,
2.0
,
/* Te */
2.06
,
1.98
,
2.16
,
2.1
,
2.0
,
/* La */
2.0
,
2.0
,
2.0
,
2.0
,
2.0
,
2.0
,
/* Eu */
2.0
,
2.0
,
2.0
,
2.0
,
2.0
,
/* Er */
2.0
,
2.0
,
2.0
,
2.0
,
2.0
,
2.0
,
/* W */
2.0
,
2.0
,
2.0
,
2.0
,
1.72
,
1.66
,
/* Hg */
1.55
,
1.96
,
2.02
,
2.0
,
2.0
,
2.0
,
2.0
,
/* Fr */
2.0
,
2.0
,
2.0
,
2.0
,
2.0
,
1.86
,
/* Np */
2.0
,
2.0
,
2.0
,
2.0
,
2.0
,
2.0
,
2.0
,
2.0
,
/* Md */
2.0
,
2.0
,
2.0
,
2.0
,
2.0
,
2.0
,
2.0
,
2.0
,
/* Mt */
2.0
,
2.0
,
2.0
};
/* lookup functions */
static
const
char
*
get_pte_label
(
const
int
idx
)
{
if
((
idx
<
1
)
||
(
idx
>=
nr_pte_entries
))
return
pte_label
[
0
];
return
pte_label
[
idx
];
}
static
float
get_pte_mass
(
const
int
idx
)
{
if
((
idx
<
1
)
||
(
idx
>=
nr_pte_entries
))
return
pte_mass
[
0
];
return
pte_mass
[
idx
];
}
static
float
get_pte_vdw_radius
(
const
int
idx
)
{
if
((
idx
<
1
)
||
(
idx
>=
nr_pte_entries
))
return
pte_vdw_radius
[
0
];
#if 1
/* Replace with Hydrogen radius with an "all-atom" radius */
if
(
idx
==
1
)
return
1.0
;
/* H */
#else
/* Replace with old VMD atom radii values */
switch
(
idx
)
{
case
1
:
return
1.0
;
/* H */
case
6
:
return
1.5
;
/* C */
case
7
:
return
1.4
;
/* N */
case
8
:
return
1.3
;
/* O */
case
9
:
return
1.2
;
/* F */
case
15
:
return
1.5
;
/* P */
case
16
:
return
1.9
;
/* S */
}
#endif
return
pte_vdw_radius
[
idx
];
}
static
int
get_pte_idx_from_string
(
const
char
*
label
)
{
int
i
,
ind
;
char
atom
[
3
];
if
(
label
!=
NULL
)
{
/* zap string */
atom
[
0
]
=
atom
[
1
]
=
atom
[
2
]
=
'\0'
;
for
(
ind
=
0
,
i
=
0
;
(
ind
<
2
)
&&
(
label
[
i
]
!=
'\0'
);
i
++
)
{
if
(
label
[
i
]
!=
' '
)
{
atom
[
ind
]
=
toupper
(
label
[
i
]);
ind
++
;
}
}
if
(
ind
<
1
)
return
0
;
/* no non-whitespace characters */
for
(
i
=
0
;
i
<
nr_pte_entries
;
++
i
)
{
if
((
toupper
(
pte_label
[
i
][
0
])
==
atom
[
0
])
&&
(
toupper
(
pte_label
[
i
][
1
])
==
atom
[
1
]))
return
i
;
}
}
return
0
;
}
// directory traversal helper functions
#if defined(_WIN32)
// Win32 directory traversal handle
typedef
struct
{
HANDLE
h
;
WIN32_FIND_DATA
fd
;
char
*
name
;
char
*
searchname
;
int
dlen
;
}
dirhandle_t
;
// open a directory handle
static
dirhandle_t
*
my_opendir
(
const
char
*
dirname
)
{
dirhandle_t
*
d
;
int
len
;
if
(
dirname
==
NULL
)
return
NULL
;
d
=
new
dirhandle_t
;
len
=
2
+
strlen
(
dirname
);
d
->
name
=
new
char
[
len
];
strcpy
(
d
->
name
,
dirname
);
strcat
(
d
->
name
,
"
\\
"
);
d
->
dlen
=
len
;
len
+=
1
;
d
->
searchname
=
new
char
[
len
];
strcpy
(
d
->
searchname
,
dirname
);
strcat
(
d
->
searchname
,
"
\\
*"
);
d
->
h
=
FindFirstFile
(
d
->
searchname
,
&
(
d
->
fd
));
if
(
d
->
h
==
((
HANDLE
)(
-
1
)))
{
delete
[]
d
->
searchname
;
delete
[]
d
->
name
;
delete
d
;
return
NULL
;
}
return
d
;
}
// get next file name from directory handle
static
char
*
my_readdir
(
dirhandle_t
*
d
)
{
if
(
FindNextFile
(
d
->
h
,
&
(
d
->
fd
)))
{
return
d
->
fd
.
cFileName
;
}
return
NULL
;
}
// close directory handle
static
void
my_closedir
(
dirhandle_t
*
d
)
{
if
(
d
->
h
!=
NULL
)
{
FindClose
(
d
->
h
);
}
delete
[]
d
->
searchname
;
delete
[]
d
->
name
;
delete
d
;
}
// open a shared object file
static
void
*
my_dlopen
(
const
char
*
fname
)
{
return
(
void
*
)
LoadLibrary
(
fname
);
}
// resolve a symbol in shared object
static
void
*
my_dlsym
(
void
*
h
,
const
char
*
sym
)
{
return
(
void
*
)
GetProcAddress
((
HINSTANCE
)
h
,
sym
);
}
// close a shared object
static
int
my_dlclose
(
void
*
h
)
{
/* FreeLibrary returns nonzero on success */
return
!
FreeLibrary
((
HINSTANCE
)
h
);
}
#else
// Unix directory traversal handle
typedef
struct
{
DIR
*
d
;
char
*
name
;
int
dlen
;
}
dirhandle_t
;
// open a directory handle
static
dirhandle_t
*
my_opendir
(
const
char
*
dirname
)
{
dirhandle_t
*
d
;
int
len
;
if
(
dirname
==
NULL
)
return
NULL
;
d
=
new
dirhandle_t
;
len
=
2
+
strlen
(
dirname
);
d
->
name
=
new
char
[
len
];
strcpy
(
d
->
name
,
dirname
);
strcat
(
d
->
name
,
"/"
);
d
->
dlen
=
len
;
d
->
d
=
opendir
(
d
->
name
);
if
(
d
->
d
==
NULL
)
{
delete
[]
d
->
name
;
delete
d
;
return
NULL
;
}
return
d
;
}
// get next file name from directory handle
static
char
*
my_readdir
(
dirhandle_t
*
d
)
{
struct
dirent
*
p
;
if
((
p
=
readdir
(
d
->
d
))
!=
NULL
)
{
return
p
->
d_name
;
}
return
NULL
;
}
// close directory handle
static
void
my_closedir
(
dirhandle_t
*
d
)
{
if
(
d
->
d
!=
NULL
)
{
closedir
(
d
->
d
);
}
delete
[]
d
->
name
;
delete
d
;
return
;
}
// open a shared object file
static
void
*
my_dlopen
(
const
char
*
fname
)
{
return
dlopen
(
fname
,
RTLD_NOW
);
}
// resolve a symbol in shared object
static
void
*
my_dlsym
(
void
*
h
,
const
char
*
sym
)
{
return
dlsym
(
h
,
sym
);
}
// close a shared object
static
int
my_dlclose
(
void
*
h
)
{
return
dlclose
(
h
);
}
#endif
}
// end of extern "C" region
using
namespace
LAMMPS_NS
;
// constructor.
MolfileInterface
::
MolfileInterface
(
const
char
*
type
,
const
int
mode
)
:
_plugin
(
0
),
_dso
(
0
),
_ptr
(
0
),
_info
(
0
),
_natoms
(
0
),
_mode
(
mode
),
_caps
(
M_NONE
)
{
_name
=
new
char
[
5
];
strcpy
(
_name
,
"none"
);
_type
=
new
char
[
1
+
strlen
(
type
)];
strcpy
(
_type
,
type
);
}
// destructor.
MolfileInterface
::~
MolfileInterface
()
{
forget_plugin
();
if
(
_info
)
{
molfile_atom_t
*
a
=
static_cast
<
molfile_atom_t
*>
(
_info
);
delete
[]
a
;
_info
=
NULL
;
}
delete
[]
_name
;
delete
[]
_type
;
}
// register the best matching plugin in a given directory
int
MolfileInterface
::
find_plugin
(
const
char
*
pluginpath
)
{
dirhandle_t
*
dir
;
char
*
filename
,
*
ext
,
*
next
,
*
path
,
*
plugindir
;
int
retval
=
E_NONE
;
#if defined(_WIN32)
#define MY_PATHSEP ';'
#else
#define MY_PATHSEP ':'
#endif
if
(
pluginpath
==
NULL
)
return
E_DIR
;
plugindir
=
path
=
strdup
(
pluginpath
);
while
(
plugindir
)
{
// check if this a single directory or path.
next
=
strchr
(
plugindir
,
MY_PATHSEP
);
if
(
next
)
{
*
next
=
'\0'
;
++
next
;
}
dir
=
my_opendir
(
plugindir
);
if
(
!
dir
)
retval
=
(
retval
>
E_DIR
)
?
retval
:
E_DIR
;
// search for suitable file names and try to inspect them
while
(
dir
)
{
char
*
fullname
;
int
len
;
filename
=
my_readdir
(
dir
);
if
(
filename
==
NULL
)
break
;
// only look at .so files
ext
=
strrchr
(
filename
,
'.'
);
if
(
ext
==
NULL
)
continue
;
if
(
strcasecmp
(
ext
,
".so"
)
!=
0
)
continue
;
// construct full pathname of potential DSO
len
=
dir
->
dlen
;
len
+=
strlen
(
filename
);
fullname
=
new
char
[
len
];
strcpy
(
fullname
,
dir
->
name
);
strcat
(
fullname
,
filename
);
// try to register plugin at file name.
int
rv
=
load_plugin
(
fullname
);
if
(
rv
>
retval
)
retval
=
rv
;
delete
[]
fullname
;
}
if
(
dir
)
my_closedir
(
dir
);
plugindir
=
next
;
}
free
(
path
);
return
retval
;
}
// register the best matching plugin in a given directory
int
MolfileInterface
::
load_plugin
(
const
char
*
filename
)
{
void
*
dso
;
int
len
,
retval
=
E_NONE
;
// access shared object
dso
=
my_dlopen
(
filename
);
if
(
dso
==
NULL
)
return
E_FILE
;
// check for required plugin symbols
void
*
ifunc
=
my_dlsym
(
dso
,
"vmdplugin_init"
);
void
*
rfunc
=
my_dlsym
(
dso
,
"vmdplugin_register"
);
void
*
ffunc
=
my_dlsym
(
dso
,
"vmdplugin_fini"
);
if
(
ifunc
==
NULL
||
rfunc
==
NULL
||
ffunc
==
NULL
)
{
my_dlclose
(
dso
);
return
E_SYMBOL
;
}
// initialize plugin. skip plugin if it fails.
if
(((
initfunc
)(
ifunc
))())
{
my_dlclose
(
dso
);
return
E_SYMBOL
;
}
// pre-register plugin.
// the callback will be called for each plugin in the DSO and
// check the file type. plugin->name will change if successful.
plugin_reginfo_t
reginfo
;
reginfo
.
p
=
NULL
;
reginfo
.
name
=
_type
;
((
regfunc
)
rfunc
)(
&
reginfo
,
plugin_register_cb
);
// make some checks to see if the plugin is suitable or not.
molfile_plugin_t
*
plugin
=
static_cast
<
molfile_plugin_t
*>
(
reginfo
.
p
);
// if the callback found a matching plugin and copied the struct,
// its name element will point to a different location now.
if
(
plugin
==
NULL
)
{
retval
=
E_TYPE
;
// check if the ABI matches the one used to compile this code
}
else
if
(
plugin
->
abiversion
!=
vmdplugin_ABIVERSION
)
{
retval
=
E_ABI
;
// check if (basic) reading is supported
}
else
if
((
_mode
&
M_READ
)
&&
(
(
plugin
->
open_file_read
==
NULL
)
||
(
plugin
->
read_next_timestep
==
NULL
)
||
(
plugin
->
close_file_read
==
NULL
)
))
{
retval
=
E_MODE
;
// check if (basic) writing is supported
}
else
if
(
(
_mode
&
M_WRITE
)
&&
(
(
plugin
->
open_file_write
==
NULL
)
||
(
plugin
->
write_timestep
==
NULL
)
||
(
plugin
->
close_file_write
==
NULL
)
))
{
retval
=
E_MODE
;
// make some additional check, if we
// already have a plugin registered.
// NOTE: this has to come last.
}
else
if
(
_dso
&&
_plugin
)
{
molfile_plugin_t
*
p
;
p
=
static_cast
<
molfile_plugin_t
*>
(
_plugin
);
// check if the new plugin is of a newer major version
if
(
p
->
majorv
>
plugin
->
majorv
)
{
retval
=
E_VERSION
;
// check if the new plugin is of a newer minor version
}
else
if
(
(
p
->
majorv
==
plugin
->
majorv
)
&&
(
p
->
minorv
>=
plugin
->
minorv
)
)
{
retval
=
E_VERSION
;
}
}
// bingo! this one is a keeper.
if
(
retval
==
E_NONE
)
{
// make sure any existing plugin is wiped out
forget_plugin
();
delete
[]
_name
;
len
=
16
;
len
+=
strlen
(
plugin
->
prettyname
);
len
+=
strlen
(
plugin
->
author
);
_name
=
new
char
[
len
];
sprintf
(
_name
,
"%s v%d.%d by %s"
,
plugin
->
prettyname
,
plugin
->
majorv
,
plugin
->
minorv
,
plugin
->
author
);
// determine plugin capabilities
_caps
=
M_NONE
;
if
(
plugin
->
read_next_timestep
)
_caps
|=
M_READ
;
if
(
plugin
->
write_timestep
)
_caps
|=
M_WRITE
;
#if vmdplugin_ABIVERSION > 10
// required to tell if velocities are present
if
(
plugin
->
read_timestep_metadata
)
_caps
|=
M_RVELS
;
// we can always offer velocities. we may not know if
// they will be written by the plugin though.
if
(
plugin
->
write_timestep
)
_caps
|=
M_WVELS
;
#endif
if
(
plugin
->
read_structure
)
_caps
|=
M_RSTRUCT
;
if
(
plugin
->
write_structure
)
_caps
|=
M_WSTRUCT
;
if
(
plugin
->
read_bonds
)
_caps
|=
M_RBONDS
;
if
(
plugin
->
write_bonds
)
_caps
|=
M_WBONDS
;
if
(
plugin
->
read_angles
)
_caps
|=
M_RANGLES
;
if
(
plugin
->
write_angles
)
_caps
|=
M_WANGLES
;
if
(
plugin
->
read_volumetric_data
)
_caps
|=
M_RVOL
;
if
(
plugin
->
write_volumetric_data
)
_caps
|=
M_WVOL
;
if
(
_mode
&
M_WRITE
)
{
_mode
|=
(
_caps
&
M_WSTRUCT
);
_mode
|=
(
_caps
&
M_WVELS
);
}
else
if
(
_mode
&
M_READ
)
{
_mode
|=
(
_caps
&
M_RSTRUCT
);
_mode
|=
(
_caps
&
M_RVELS
);
}
_plugin
=
plugin
;
_dso
=
dso
;
return
E_MATCH
;
}
// better luck next time. clean up and return.
my_dlclose
(
dso
);
return
retval
;
}
// deregister a plugin and close or reset all associated objects.
void
MolfileInterface
::
forget_plugin
()
{
if
(
_ptr
)
close
();
if
(
_plugin
)
_plugin
=
NULL
;
if
(
_dso
)
{
void
*
ffunc
=
my_dlsym
(
_dso
,
"vmdplugin_fini"
);
if
(
ffunc
)
((
finifunc
)
ffunc
)();
my_dlclose
(
_dso
);
}
_dso
=
NULL
;
delete
[]
_name
;
_name
=
new
char
[
5
];
strcpy
(
_name
,
"none"
);
_caps
=
M_NONE
;
}
// open file for reading or writing
int
MolfileInterface
::
open
(
const
char
*
name
,
int
*
natoms
)
{
if
(
!
_plugin
||
!
_dso
||
!
natoms
)
return
E_FILE
;
molfile_plugin_t
*
p
=
static_cast
<
molfile_plugin_t
*>
(
_plugin
);
if
(
_mode
&
M_WRITE
)
_ptr
=
p
->
open_file_write
(
name
,
_type
,
*
natoms
);
else
if
(
_mode
&
M_READ
)
_ptr
=
p
->
open_file_read
(
name
,
_type
,
natoms
);
if
(
_ptr
==
NULL
)
return
E_FILE
;
_natoms
=
*
natoms
;
// we need to deal with structure information,
// so we allocate and initialize storage for it.
if
(
_mode
&
(
M_RSTRUCT
|
M_WSTRUCT
))
{
molfile_atom_t
*
a
=
new
molfile_atom_t
[
_natoms
];
_info
=
a
;
memset
(
_info
,
0
,
_natoms
*
sizeof
(
molfile_atom_t
));
for
(
int
i
=
0
;
i
<
_natoms
;
++
i
)
{
a
[
i
].
name
[
0
]
=
'X'
;
a
[
i
].
type
[
0
]
=
a
[
i
].
resname
[
0
]
=
a
[
i
].
segid
[
0
]
=
'U'
;
a
[
i
].
type
[
1
]
=
a
[
i
].
resname
[
1
]
=
a
[
i
].
segid
[
1
]
=
'N'
;
a
[
i
].
type
[
2
]
=
a
[
i
].
resname
[
2
]
=
a
[
i
].
segid
[
2
]
=
'K'
;
a
[
i
].
chain
[
0
]
=
'X'
;
}
}
return
E_NONE
;
}
// get of set atom structure information
int
MolfileInterface
::
structure
()
{
if
(
!
_plugin
||
!
_dso
)
return
E_FILE
;
molfile_plugin_t
*
p
=
static_cast
<
molfile_plugin_t
*>
(
_plugin
);
int
optflags
=
MOLFILE_NOOPTIONS
;
if
(
_mode
&
M_WSTRUCT
)
{
optflags
|=
(
_props
&
P_BFAC
)
?
MOLFILE_BFACTOR
:
0
;
optflags
|=
(
_props
&
P_OCCP
)
?
MOLFILE_OCCUPANCY
:
0
;
optflags
|=
(
_props
&
P_MASS
)
?
MOLFILE_MASS
:
0
;
optflags
|=
(
_props
&
P_CHRG
)
?
MOLFILE_CHARGE
:
0
;
optflags
|=
(
_props
&
P_RADS
)
?
MOLFILE_RADIUS
:
0
;
optflags
|=
(
_props
&
P_ATMN
)
?
MOLFILE_ATOMICNUMBER
:
0
;
molfile_atom_t
*
a
=
static_cast
<
molfile_atom_t
*>
(
_info
);
p
->
write_structure
(
_ptr
,
optflags
,
a
);
}
else
if
(
_mode
&
M_RSTRUCT
)
{
molfile_atom_t
*
a
=
static_cast
<
molfile_atom_t
*>
(
_info
);
p
->
read_structure
(
_ptr
,
&
optflags
,
a
);
// mandatory properties
_props
=
P_NAME
|
P_TYPE
|
P_RESN
|
P_RESI
|
P_SEGN
|
P_CHAI
;
// optional properties
_props
|=
(
optflags
&
MOLFILE_BFACTOR
)
?
P_BFAC
:
0
;
_props
|=
(
optflags
&
MOLFILE_OCCUPANCY
)
?
P_OCCP
:
0
;
_props
|=
(
optflags
&
MOLFILE_MASS
)
?
P_MASS
:
0
;
_props
|=
(
optflags
&
MOLFILE_CHARGE
)
?
P_CHRG
:
0
;
_props
|=
(
optflags
&
MOLFILE_RADIUS
)
?
P_RADS
:
0
;
_props
|=
(
optflags
&
MOLFILE_ATOMICNUMBER
)
?
P_ATMN
:
0
;
}
return
0
;
}
// safely close file
int
MolfileInterface
::
close
()
{
if
(
!
_plugin
||
!
_dso
||
!
_ptr
)
return
E_FILE
;
molfile_plugin_t
*
p
=
static_cast
<
molfile_plugin_t
*>
(
_plugin
);
if
(
_mode
&
M_WRITE
)
{
p
->
close_file_write
(
_ptr
);
}
else
if
(
_mode
&
M_READ
)
{
p
->
close_file_read
(
_ptr
);
}
if
(
_info
)
{
molfile_atom_t
*
a
=
static_cast
<
molfile_atom_t
*>
(
_info
);
delete
[]
a
;
_info
=
NULL
;
}
_ptr
=
NULL
;
_natoms
=
0
;
return
E_NONE
;
}
// read or write timestep
int
MolfileInterface
::
timestep
(
float
*
coords
,
float
*
vels
,
float
*
cell
,
double
*
simtime
)
{
if
(
!
_plugin
||
!
_dso
||
!
_ptr
)
return
1
;
molfile_plugin_t
*
p
=
static_cast
<
molfile_plugin_t
*>
(
_plugin
);
molfile_timestep_t
*
t
=
new
molfile_timestep_t
;
int
rv
;
if
(
_mode
&
M_WRITE
)
{
t
->
coords
=
coords
;
t
->
velocities
=
vels
;
if
(
cell
!=
NULL
)
{
t
->
A
=
cell
[
0
];
t
->
B
=
cell
[
1
];
t
->
C
=
cell
[
2
];
t
->
alpha
=
cell
[
3
];
t
->
beta
=
cell
[
4
];
t
->
gamma
=
cell
[
5
];
}
else
{
t
->
A
=
0.0f
;
t
->
B
=
0.0f
;
t
->
C
=
0.0f
;
t
->
alpha
=
90.0f
;
t
->
beta
=
90.0f
;
t
->
gamma
=
90.0f
;
}
if
(
simtime
)
t
->
physical_time
=
*
simtime
;
else
t
->
physical_time
=
0.0
;
rv
=
p
->
write_timestep
(
_ptr
,
t
);
}
else
{
// no coordinate storage => skip step
if
(
coords
==
NULL
)
{
rv
=
p
->
read_next_timestep
(
_ptr
,
_natoms
,
NULL
);
}
else
{
t
->
coords
=
coords
;
t
->
velocities
=
vels
;
t
->
A
=
0.0f
;
t
->
B
=
0.0f
;
t
->
C
=
0.0f
;
t
->
alpha
=
90.0f
;
t
->
beta
=
90.0f
;
t
->
gamma
=
90.0f
;
t
->
physical_time
=
0.0
;
rv
=
p
->
read_next_timestep
(
_ptr
,
_natoms
,
t
);
if
(
cell
!=
NULL
)
{
cell
[
0
]
=
t
->
A
;
cell
[
1
]
=
t
->
B
;
cell
[
2
]
=
t
->
C
;
cell
[
3
]
=
t
->
alpha
;
cell
[
4
]
=
t
->
beta
;
cell
[
5
]
=
t
->
gamma
;
}
if
(
simtime
)
*
simtime
=
t
->
physical_time
;
}
if
(
rv
==
MOLFILE_EOF
)
{
delete
t
;
return
1
;
}
}
delete
t
;
return
0
;
}
// functions to read properties from molfile structure
#define PROPUPDATE(PROP,ENTRY,VAL) \
if (propid == PROP) { VAL = a.ENTRY; }
#define PROPSTRCPY(PROP,ENTRY,VAL) \
if (propid == PROP) { strcpy(VAL,a.ENTRY); }
// single precision floating point props
static
float
read_float_property
(
molfile_atom_t
&
a
,
const
int
propid
)
{
float
prop
=
0.0f
;
int
iprop
=
0
;
PROPUPDATE
(
MolfileInterface
::
P_OCCP
,
occupancy
,
prop
);
PROPUPDATE
(
MolfileInterface
::
P_BFAC
,
bfactor
,
prop
);
PROPUPDATE
(
MolfileInterface
::
P_MASS
,
mass
,
prop
);
PROPUPDATE
(
MolfileInterface
::
P_CHRG
,
charge
,
prop
);
PROPUPDATE
(
MolfileInterface
::
P_RADS
,
radius
,
prop
);
PROPUPDATE
((
MolfileInterface
::
P_ATMN
|
MolfileInterface
::
P_MASS
),
atomicnumber
,
iprop
);
PROPUPDATE
((
MolfileInterface
::
P_ATMN
|
MolfileInterface
::
P_RADS
),
atomicnumber
,
iprop
);
if
(
propid
&
MolfileInterface
::
P_ATMN
)
{
if
(
propid
&
MolfileInterface
::
P_MASS
)
prop
=
get_pte_mass
(
iprop
);
if
(
propid
&
MolfileInterface
::
P_RADS
)
prop
=
get_pte_vdw_radius
(
iprop
);
}
return
prop
;
}
// integer and derived props
static
int
read_int_property
(
molfile_atom_t
&
a
,
const
int
propid
)
{
int
prop
=
0
;
const
char
*
sprop
;
PROPUPDATE
(
MolfileInterface
::
P_RESI
,
resid
,
prop
);
PROPUPDATE
(
MolfileInterface
::
P_ATMN
,
atomicnumber
,
prop
);
PROPUPDATE
((
MolfileInterface
::
P_ATMN
|
MolfileInterface
::
P_NAME
),
name
,
sprop
);
PROPUPDATE
((
MolfileInterface
::
P_ATMN
|
MolfileInterface
::
P_TYPE
),
type
,
sprop
);
if
(
propid
&
MolfileInterface
::
P_ATMN
)
{
if
(
propid
&
(
MolfileInterface
::
P_NAME
|
MolfileInterface
::
P_TYPE
))
prop
=
get_pte_idx_from_string
(
sprop
);
}
return
prop
;
}
// string and derived props
static
const
char
*
read_string_property
(
molfile_atom_t
&
a
,
const
int
propid
)
{
const
char
*
prop
=
NULL
;
int
iprop
=
0
;
PROPUPDATE
(
MolfileInterface
::
P_NAME
,
name
,
prop
);
PROPUPDATE
(
MolfileInterface
::
P_TYPE
,
type
,
prop
);
PROPUPDATE
(
MolfileInterface
::
P_RESN
,
resname
,
prop
);
PROPUPDATE
(
MolfileInterface
::
P_SEGN
,
segid
,
prop
);
PROPUPDATE
((
MolfileInterface
::
P_ATMN
|
MolfileInterface
::
P_NAME
),
atomicnumber
,
iprop
);
PROPUPDATE
((
MolfileInterface
::
P_ATMN
|
MolfileInterface
::
P_TYPE
),
atomicnumber
,
iprop
);
if
(
propid
&
MolfileInterface
::
P_ATMN
)
{
if
(
propid
&
(
MolfileInterface
::
P_NAME
|
MolfileInterface
::
P_TYPE
))
prop
=
get_pte_label
(
iprop
);
}
return
prop
;
}
#undef PROPUPDATE
#undef PROPSTRCPY
// functions to store properties into molfile structure
#define PROPUPDATE(PROP,ENTRY,VAL) \
if ((propid & PROP) == PROP) { a.ENTRY = VAL; plist |= PROP; }
#define PROPSTRCPY(PROP,ENTRY,VAL) \
if ((propid & PROP) == PROP) { strcpy(a.ENTRY,VAL); plist |= PROP; }
// floating point props
static
int
write_atom_property
(
molfile_atom_t
&
a
,
const
int
propid
,
const
float
prop
)
{
int
plist
=
MolfileInterface
::
P_NONE
;
PROPUPDATE
(
MolfileInterface
::
P_OCCP
,
occupancy
,
prop
);
PROPUPDATE
(
MolfileInterface
::
P_BFAC
,
bfactor
,
prop
);
PROPUPDATE
(
MolfileInterface
::
P_MASS
,
mass
,
prop
);
PROPUPDATE
(
MolfileInterface
::
P_CHRG
,
charge
,
prop
);
PROPUPDATE
(
MolfileInterface
::
P_RADS
,
radius
,
prop
);
return
plist
;
}
// double precision floating point props
static
int
write_atom_property
(
molfile_atom_t
&
a
,
const
int
propid
,
const
double
prop
)
{
return
write_atom_property
(
a
,
propid
,
static_cast
<
float
>
(
prop
));
}
// integer and derived props
static
int
write_atom_property
(
molfile_atom_t
&
a
,
const
int
propid
,
const
int
prop
)
{
int
plist
=
MolfileInterface
::
P_NONE
;
PROPUPDATE
(
MolfileInterface
::
P_RESI
,
resid
,
prop
);
PROPUPDATE
(
MolfileInterface
::
P_ATMN
,
atomicnumber
,
prop
);
PROPUPDATE
((
MolfileInterface
::
P_ATMN
|
MolfileInterface
::
P_MASS
),
mass
,
get_pte_mass
(
prop
));
PROPSTRCPY
((
MolfileInterface
::
P_ATMN
|
MolfileInterface
::
P_NAME
),
name
,
get_pte_label
(
prop
));
PROPSTRCPY
((
MolfileInterface
::
P_ATMN
|
MolfileInterface
::
P_TYPE
),
type
,
get_pte_label
(
prop
));
return
plist
;
}
// integer and derived props
static
int
write_atom_property
(
molfile_atom_t
&
a
,
const
int
propid
,
const
char
*
prop
)
{
int
plist
=
MolfileInterface
::
P_NONE
;
PROPSTRCPY
(
MolfileInterface
::
P_NAME
,
name
,
prop
);
PROPSTRCPY
(
MolfileInterface
::
P_TYPE
,
type
,
prop
);
PROPSTRCPY
(
MolfileInterface
::
P_RESN
,
resname
,
prop
);
PROPSTRCPY
(
MolfileInterface
::
P_SEGN
,
segid
,
prop
);
return
plist
;
}
#undef PROPUPDATE
#undef PROPSTRCPY
// set/get atom floating point property
int
MolfileInterface
::
property
(
int
propid
,
int
idx
,
float
*
prop
)
{
if
((
_info
==
NULL
)
||
(
prop
==
NULL
)
||
(
idx
<
0
)
||
(
idx
>=
_natoms
))
return
P_NONE
;
molfile_atom_t
*
a
=
static_cast
<
molfile_atom_t
*>
(
_info
);
if
(
_mode
&
M_WSTRUCT
)
_props
|=
write_atom_property
(
a
[
idx
],
propid
,
*
prop
);
if
(
_mode
&
M_RSTRUCT
)
*
prop
=
read_float_property
(
a
[
idx
],
propid
);
return
_props
;
}
// set/get per type floating point property
int
MolfileInterface
::
property
(
int
propid
,
int
*
types
,
float
*
prop
)
{
if
((
_info
==
NULL
)
||
(
types
==
NULL
)
||
(
prop
==
NULL
))
return
P_NONE
;
molfile_atom_t
*
a
=
static_cast
<
molfile_atom_t
*>
(
_info
);
if
(
_mode
&
M_WSTRUCT
)
{
for
(
int
i
=
0
;
i
<
_natoms
;
++
i
)
_props
|=
write_atom_property
(
a
[
i
],
propid
,
prop
[
types
[
i
]]);
}
// useless for reading.
if
(
_mode
&
M_RSTRUCT
)
return
P_NONE
;
return
_props
;
}
// set/get per atom floating point property
int
MolfileInterface
::
property
(
int
propid
,
float
*
prop
)
{
if
((
_info
==
NULL
)
||
(
prop
==
NULL
))
return
P_NONE
;
molfile_atom_t
*
a
=
static_cast
<
molfile_atom_t
*>
(
_info
);
if
(
_mode
&
M_WSTRUCT
)
{
for
(
int
i
=
0
;
i
<
_natoms
;
++
i
)
_props
|=
write_atom_property
(
a
[
i
],
propid
,
prop
[
i
]);
}
if
(
_mode
&
M_RSTRUCT
)
{
for
(
int
i
=
0
;
i
<
_natoms
;
++
i
)
prop
[
i
]
=
read_float_property
(
a
[
i
],
propid
);
}
return
_props
;
}
// set/get atom floating point property
int
MolfileInterface
::
property
(
int
propid
,
int
idx
,
double
*
prop
)
{
if
((
_info
==
NULL
)
||
(
prop
==
NULL
)
||
(
idx
<
0
)
||
(
idx
>=
_natoms
))
return
P_NONE
;
molfile_atom_t
*
a
=
static_cast
<
molfile_atom_t
*>
(
_info
);
if
(
_mode
&
M_WSTRUCT
)
return
write_atom_property
(
a
[
idx
],
propid
,
*
prop
);
if
(
_mode
&
M_RSTRUCT
)
*
prop
=
static_cast
<
double
>
(
read_float_property
(
a
[
idx
],
propid
));
return
_props
;
}
// set/get per type floating point property
int
MolfileInterface
::
property
(
int
propid
,
int
*
types
,
double
*
prop
)
{
if
((
_info
==
NULL
)
||
(
types
==
NULL
)
||
(
prop
==
NULL
))
return
P_NONE
;
molfile_atom_t
*
a
=
static_cast
<
molfile_atom_t
*>
(
_info
);
if
(
_mode
&
M_WSTRUCT
)
{
for
(
int
i
=
0
;
i
<
_natoms
;
++
i
)
_props
|=
write_atom_property
(
a
[
i
],
propid
,
prop
[
types
[
i
]]);
}
// useless for reading
if
(
_mode
&
M_RSTRUCT
)
return
P_NONE
;
return
_props
;
}
// set/get per atom floating point property
int
MolfileInterface
::
property
(
int
propid
,
double
*
prop
)
{
if
((
_info
==
NULL
)
||
(
prop
==
NULL
))
return
P_NONE
;
molfile_atom_t
*
a
=
static_cast
<
molfile_atom_t
*>
(
_info
);
if
(
_mode
&
M_WSTRUCT
)
{
for
(
int
i
=
0
;
i
<
_natoms
;
++
i
)
_props
|=
write_atom_property
(
a
[
i
],
propid
,
prop
[
i
]);
}
if
(
_mode
&
M_RSTRUCT
)
{
for
(
int
i
=
0
;
i
<
_natoms
;
++
i
)
prop
[
i
]
=
static_cast
<
double
>
(
read_float_property
(
a
[
i
],
propid
));
}
return
_props
;
}
#define INT_TO_STRING_BODY(IDX) \
buf[15] = 0; \
if (propid & P_NAME) \
_props |= write_atom_property(a[IDX],P_NAME,buf); \
if (propid & P_TYPE) \
_props |= write_atom_property(a[IDX],P_TYPE,buf); \
buf[7] = 0; \
if (propid & P_RESN) \
_props |= write_atom_property(a[IDX],P_RESN,buf); \
if (propid & P_SEGN) \
_props |= write_atom_property(a[IDX],P_SEGN,buf); \
buf[1] = 0; \
if (propid & P_CHAI) \
_props |= write_atom_property(a[IDX],P_CHAI,buf)
// set/get atom integer property
int
MolfileInterface
::
property
(
int
propid
,
int
idx
,
int
*
prop
)
{
if
((
_info
==
NULL
)
||
(
prop
==
NULL
)
||
(
idx
<
0
)
||
(
idx
>=
_natoms
))
return
P_NONE
;
molfile_atom_t
*
a
=
static_cast
<
molfile_atom_t
*>
(
_info
);
if
(
_mode
&
M_WSTRUCT
)
{
char
buf
[
64
];
_props
|=
write_atom_property
(
a
[
idx
],
propid
,
*
prop
);
if
(
propid
&
(
P_NAME
|
P_TYPE
|
P_RESN
|
P_SEGN
|
P_CHAI
))
{
sprintf
(
buf
,
"%d"
,
*
prop
);
INT_TO_STRING_BODY
(
idx
);
}
}
if
(
_mode
&
M_RSTRUCT
)
*
prop
=
read_int_property
(
a
[
idx
],
propid
);
return
_props
;
}
// set/get per type integer property
int
MolfileInterface
::
property
(
int
propid
,
int
*
types
,
int
*
prop
)
{
if
((
_info
==
NULL
)
||
(
types
==
NULL
)
||
(
prop
==
NULL
))
return
P_NONE
;
molfile_atom_t
*
a
=
static_cast
<
molfile_atom_t
*>
(
_info
);
if
(
_mode
&
M_WSTRUCT
)
{
char
buf
[
64
];
for
(
int
i
=
0
;
i
<
_natoms
;
++
i
)
_props
|=
write_atom_property
(
a
[
i
],
propid
,
prop
[
types
[
i
]]);
if
(
propid
&
(
P_NAME
|
P_TYPE
|
P_RESN
|
P_SEGN
|
P_CHAI
))
{
for
(
int
i
=
0
;
i
<
_natoms
;
++
i
)
{
sprintf
(
buf
,
"%d"
,
prop
[
types
[
i
]]);
INT_TO_STRING_BODY
(
i
);
}
}
}
// useless when reading
if
(
_mode
&
M_RSTRUCT
)
return
P_NONE
;
return
_props
;
}
// set/get per atom integer property
int
MolfileInterface
::
property
(
int
propid
,
int
*
prop
)
{
if
((
_info
==
NULL
)
||
(
prop
==
NULL
))
return
P_NONE
;
molfile_atom_t
*
a
=
static_cast
<
molfile_atom_t
*>
(
_info
);
if
(
_mode
&
M_WSTRUCT
)
{
char
buf
[
64
];
for
(
int
i
=
0
;
i
<
_natoms
;
++
i
)
_props
|=
write_atom_property
(
a
[
i
],
propid
,
prop
[
i
]);
if
(
propid
&
(
P_NAME
|
P_TYPE
|
P_RESN
|
P_SEGN
|
P_CHAI
))
{
for
(
int
i
=
0
;
i
<
_natoms
;
++
i
)
{
sprintf
(
buf
,
"%d"
,
prop
[
i
]);
INT_TO_STRING_BODY
(
i
);
}
}
}
if
(
_mode
&
M_RSTRUCT
)
{
for
(
int
i
=
0
;
i
<
_natoms
;
++
i
)
prop
[
i
]
=
read_int_property
(
a
[
i
],
propid
);
}
return
_props
;
}
#undef INT_TO_STRING_BODY
// set/get atom string property
int
MolfileInterface
::
property
(
int
propid
,
int
idx
,
char
*
prop
)
{
if
((
_info
==
NULL
)
||
(
prop
==
NULL
)
||
(
idx
<
0
)
||
(
idx
>=
_natoms
))
return
P_NONE
;
molfile_atom_t
*
a
=
static_cast
<
molfile_atom_t
*>
(
_info
);
if
(
_mode
&
M_WSTRUCT
)
{
_props
|=
write_atom_property
(
a
[
idx
],
propid
,
prop
);
}
if
(
_mode
&
M_RSTRUCT
)
strcpy
(
prop
,
read_string_property
(
a
[
idx
],
propid
));
return
_props
;
}
// set/get per type string property
int
MolfileInterface
::
property
(
int
propid
,
int
*
types
,
char
**
prop
)
{
if
((
_info
==
NULL
)
||
(
types
==
NULL
)
||
(
prop
==
NULL
))
return
P_NONE
;
molfile_atom_t
*
a
=
static_cast
<
molfile_atom_t
*>
(
_info
);
if
(
_mode
&
M_WSTRUCT
)
{
for
(
int
i
=
0
;
i
<
_natoms
;
++
i
)
{
_props
|=
write_atom_property
(
a
[
i
],
propid
,
prop
[
types
[
i
]]);
}
}
// useless when reading
if
(
_mode
&
M_RSTRUCT
)
return
P_NONE
;
return
_props
;
}
// set/get per atom string property
int
MolfileInterface
::
property
(
int
propid
,
char
**
prop
)
{
if
((
_info
==
NULL
)
||
(
prop
==
NULL
))
return
P_NONE
;
molfile_atom_t
*
a
=
static_cast
<
molfile_atom_t
*>
(
_info
);
if
(
_mode
&
M_WSTRUCT
)
{
for
(
int
i
=
0
;
i
<
_natoms
;
++
i
)
{
_props
|=
write_atom_property
(
a
[
i
],
propid
,
prop
[
i
]);
}
}
// not supported right now. XXX: should we use strdup() here?
if
(
_mode
&
M_RSTRUCT
)
return
P_NONE
;
return
_props
;
}
Event Timeline
Log In to Comment