Previous Thread
Next Thread
Print Thread
Page 1 of 3 1 2 3
Joined: May 2022
Posts: 15
J
Forum Member
OP Offline
Forum Member
J
Joined: May 2022
Posts: 15
Hello, I try to calculate rmds between a pdb file and NAMD coor file.

Code
bomlev -5
wrnlev  5

stream toppar.str
stream job.str

open unit 10 read card name setup.psf
read psf xplor card unit 10 append

open unit 10 read card name setup.pdb
read coord pdb unit 10

open read unform unit 22 name dynami.1.coor

define heavy select (segid HETE .and. .not. type h*) end
 
coor copy comp select heavy end
 
coor orient rms select heavy end

set rmsd ?rms

    
open write unit 30 card name rmsd_first-last-frames_coor.txt
write title unit 30
* @rmsd

stop


I obtain value like 1.310232404E-15

Does anyone know where the error is?

Joined: Sep 2003
Posts: 8,621
Likes: 24
rmv Online Content
Forum Member
Online Content
Forum Member
Joined: Sep 2003
Posts: 8,621
Likes: 24
From the CHARMM documentation, it appears that the command to read a NAMD file may be

READ NAMD FILE "filename"


Rick Venable
computational chemist

1 member likes this: judkil21
Joined: Sep 2003
Posts: 4,860
Likes: 10
Forum Member
Online Content
Forum Member
Joined: Sep 2003
Posts: 4,860
Likes: 10
As far as I can tell you did not actually read more than one coordinate set. So getting the rmsd between two identical coordinate sets is indeed expected to return zero, within numerical (double) precision.
And I hope you understand that BOMLEV -5 is an absolute No-No if you want to get help ....


Lennart Nilsson
Karolinska Institutet
Stockholm, Sweden
1 member likes this: judkil21
Joined: May 2022
Posts: 15
J
Forum Member
OP Offline
Forum Member
J
Joined: May 2022
Posts: 15
In setup.pdb and dynami.1.coor there are : a protein, a ligand in a water box.

I made a simulation setup.pdb is the first frame and dynami.1.coor is the last frame.

The value of rmsd should be between 0 and 3. Here is the value I get 2.581522112E-15

From my readings, Charmm works with main set and comp set but it still doesn't work

Joined: May 2022
Posts: 15
J
Forum Member
OP Offline
Forum Member
J
Joined: May 2022
Posts: 15
CHARMM> coor copy comp select heavy end
SELRPN> 20 atoms have been selected out of 32599
SELECTED COORDINATES COPIED TO THE COMPARISON SET.


CHARMM> coor orient rms select heavy end
SELRPN> 20 atoms have been selected out of 32599
CENTER OF ATOMS BEFORE TRANSLATION -0.73425 1.86415 0.97435
CENTER OF REFERENCE COORDINATE SET -0.73425 1.86415 0.97435
NET TRANSLATION OF ROTATED ATOMS 0.00000 0.00000 0.00000
ROTATION MATRIX
1.000000 -0.000000 -0.000000
-0.000000 1.000000 0.000000
0.000000 -0.000000 1.000000
VERY LITTLE ROTATION: NO AXIS FOUND
TOTAL SQUARE DIFF IS 0.0000 DENOMINATOR IS 20.0000
THUS RMS DIFF IS 0.000000
ALL COORDINATES ORIENTED IN THE MAIN SET BASED ON SELECTED ATOMS.


CHARMM> set rmsd ?rms
RDCMND substituted energy or value "?RMS" to "2.581522112E-15"
Parameter: RMSD <- "2.581522112E-15"

CHARMM>

CHARMM> open write unit 30 card name rmsd_first-last-frames_coor.txt
VOPEN> Attempting to open::rmsd_first-last-frames_coor.txt::
OPNLGU> Unit 30 opened for WRITE access to rmsd_first-last-frames_coor.txt

CHARMM> write title unit 30
RDTITL> * 2.581522112E-15

Joined: Sep 2003
Posts: 4,860
Likes: 10
Forum Member
Online Content
Forum Member
Joined: Sep 2003
Posts: 4,860
Likes: 10
Your statement "open read unform ....." only open the file, no reading is done.
Is dynami.1.coor a binary file? It is better to work with formatted coordinate files, either in pdb-format or the charmm format.
To avoid the confusion the obsolete construct with an open statement followe by a read or write should be avoided (open is still necessary in some cases):
read psf card append name setup.psf (I am not sure that xplor is needed)
read coor pdb name setup.pdb
read namd file dynami.1.coor

Note also that
coor copy comp select heavy end
copies coordinates from the main set (which came from setup.pdb) to the comparison set, and that the following RMS calculation is therefore performed on two identical coordinate sets.


Lennart Nilsson
Karolinska Institutet
Stockholm, Sweden
1 member likes this: judkil21
Joined: May 2022
Posts: 15
J
Forum Member
OP Offline
Forum Member
J
Joined: May 2022
Posts: 15
Thanks

Indeed dynami.1.coor is a pdb coordinate file, not a binary file, I will remove the obseleted structures

That's what I want to do is to compare the rmsd of the main set with the comparison set, I tried again by adding orient protein before the rms calculation

define heavy select (segid HETE .and. .not. type h*) end
define protein sele type ca end
coor orient rms sele protein end
coor rms select heavy end
set rmsd ?rms


I obtain


CHARMM> coor orient rms select heavy end
SELRPN> 20 atoms have been selected out of 32659
**WARNING** ALL SELECTED COORDINATES UNDEFINED

**** WARNING **** FOR THIS OPERATION, THERE WERE 20 MISSING COORDINATES

Last edited by judkil21; 09/21/22 10:06 AM.
Joined: Sep 2003
Posts: 4,860
Likes: 10
Forum Member
Online Content
Forum Member
Joined: Sep 2003
Posts: 4,860
Likes: 10
This is an outline of what you want to do:
read rtf ...
read param ...
read psf ...
read coor pdb name setup.pdb
read coor pdb name dynami.1.coor COMP
coor orie rms sele type ca end


Lennart Nilsson
Karolinska Institutet
Stockholm, Sweden
Joined: May 2022
Posts: 15
J
Forum Member
OP Offline
Forum Member
J
Joined: May 2022
Posts: 15
I have a system containing a ligand (an hetero molecule), a protein (containing amino acids) and a water box

I have realized a simulation of molecular dynamics

I would like to calculate the rsmd between, the coordinates of the ligand in the starting frame (setup.pdb, pdb format) and the coordinates of the ligand in the last frame (dynami.1.coor, pdb format)

Last edited by judkil21; 09/21/22 11:48 AM.
Joined: May 2022
Posts: 15
J
Forum Member
OP Offline
Forum Member
J
Joined: May 2022
Posts: 15
read rtf ...
read param ...
read psf ...

read psf card append name setup.psf
read coor pdb name setup.pdb
read coor pdb name dynami.1.coor COMP

define heavy select (segid HETE .and. .not. type h*) end ! For the ligands of the two structures that I want to align

coor orie rms sele heavy end


Prompt :

CHARMM> coor orie rms sele heavy end
SELRPN> 20 atoms have been selected out of 32599
**WARNING** ALL SELECTED COORDINATES UNDEFINED

**** WARNING **** FOR THIS OPERATION, THERE WERE 20 MISSING COORDINATES

Last edited by judkil21; 09/21/22 12:13 PM.
Page 1 of 3 1 2 3

Moderated by  BRBrooks, lennart, rmv 

Link Copied to Clipboard
Powered by UBB.threads™ PHP Forum Software 7.7.5
(Release build 20201027)
Responsive Width:

PHP: 7.3.31-1~deb10u1 Page Time: 0.010s Queries: 38 (0.006s) Memory: 0.7870 MB (Peak: 0.8821 MB) Data Comp: Off Server Time: 2022-11-27 13:53:02 UTC
Valid HTML 5 and Valid CSS