Re: [AMBER] leap cannot read MOL2 file that it has written before: small fix suggested

From: Jan-Philip Gehrcke <jgehrcke.googlemail.com>
Date: Thu, 10 Nov 2011 23:32:30 +0100

Hey,

yes, obviously the `sName` attribute of the container corresponding to x
is not set after loading x via loadpdb.

I've done a small test to illustrate what is missing. In mol2File.c,
around line 293, I've written:

VP0(("sName before ContainerSetName: %s\n", uUnit->cHeader.sName));
ContainerSetName(&uUnit->cHeader, "powername");
VP0(("lol cHeader: %s\n", uUnit->cHeader.sName));

The rest of the file is unchanged, especially this line is active again:

strcpy(sTemp, sContainerName((CONTAINER) uUnit));

The output is:

> Writing mol2 file: test.mol2
> sName before ContainerSetName:
> sName after ContainerSetName: powername

"powername" then appears in the MOL2 file, i.e.
`sContainerName((CONTAINER) uUnit)` returned a non-empty string.


Hence, we have to make sure that `ContainerSetName(&uUnit->cHeader,
name);` is called at some point after loading a PDB file to an object
with a proper name.

Jan-Philip


On 10.11.2011 22:50, M. L. Dodson wrote:
> Replying to myself again, it appears that when a molecule
> is read in from a pdb, there is no molecule name assigned
> (in the mol2-format sense). Maybe we should initialize the
> molecule name in the code before we even get started. If we
> need it for a molecule read from a pdb file, then the default
> is there. I suggest something like REPLACEME. If the name
> is assigned by some other way, then it overwrites the default.
>
> Bud
>
> On Nov 10, 2011, at 3:30 PM, M. L. Dodson wrote:
>
>> Replying to myself.
>>
>> If you write a pdb with x, then use it the way you describe,
>> you get the behavior you indicated. Something is going
>> wrong. What leap has after the combine is different from
>> what it has from reading the pdb file:
>>
>> source leaprc.ff99SB
>> x = combine {ACE THR NME}
>> savepdb x test.pdb
>> savemol2 x test.mol2 0
>> y = loadmol2 test.mol2
>> quit
>>
>> versus
>>
>> source leaprc.ff99SB
>> x = loadpdb test.pdb
>> savemol2 x test.mol2 0
>> y = loadmol2 test.mol2
>> quit
>>
>> One thing I noticed is that x from combine has TER cards in
>> the written pdb file, whereas the pdb file I used in the
>> second case had the TER cards removed. (But the behavior
>> you described happened with or without the TER cards.)
>>
>> Puzzled,
>> Bud
>>
>> On Nov 10, 2011, at 3:09 PM, M. L. Dodson wrote:
>>
>>> Excuse me for jumping in.
>>>
>>> Your statement is not correct in the general case. Maybe it
>>> is for your test.pdb.
>>>
>>> See attachments.
>>>
>>> Command line was: tleap -f test.in> test.log
>>>
>>> Bud Dodson
>>>
>>> On Nov 10, 2011, at 1:26 PM, Jan-Philip Gehrcke wrote:
>>>
>>>> Hello,
>>>>
>>>> in leap you can use any test.pdb file as input and try
>>>>
>>>>> source leaprc.ff99SB
>>>>> x = loadpdb test.pdb
>>>>> savemol2 x out.mol2 0
>>>>> y = loadmol2 out.mol2
>>>>
>>>> which will result in something like
>>>>
>>>>> Loading PDB file: ./test.pdb
>>>>> [...]
>>>>> Writing mol2 file: out.mol2
>>>>> [...]
>>>>> Loading Mol2 file: ./out.mol2
>>>>> Reading MOLECULE named [...]
>>>>> Fatal Error: last line read: [...]
>>>>
>>>> This is because the molecule name is not written by leap and then it is
>>>> expected during parsing. This is correct, since the molecule name is
>>>> required by the MOL2 format specification.
>>>>
>>>> Now, I've quickly looked into the code of mol2File.c.
>>>>
>>>> `sContainerName((CONTAINER) uUnit)` is used to get the object's name.
>>>> But it seems to be broken and returns an "empty string" which results in
>>>> the problem described above.
>>>>
>>>> It is really hard to figure out how to get the object's name back. After
>>>> some minutes I gave up. The code in mol2File.c itself is badly
>>>> structured and I don't know if it's worth investing time in it. Besides
>>>> fixing the molecule name issue, there is more to do.. for example it
>>>> prints "Truncating residue name for PDB format: CLEU -> CLE" (and then
>>>> nevertheless writes "CLEU" to the MOL2 file).
>>>>
>>>> To make a long story short.. a quick fix would be to just write a
>>>> default molecule name to any MOL2 file in order to at least stick to the
>>>> MOL2 format specs. This change worked for me:
>>>>
>>>>
>>>> --- src/leap/mol2File.c.old 2011-11-10 20:08:23.081930188 +0100
>>>> +++ src/leap/mol2File.c 2011-11-10 20:09:22.043721048 +0100
>>>> .. -288,7 +288,7 ..
>>>> iResidueCount = zUnitIOAmberOrderResidues( uUnit );
>>>>
>>>>
>>>> -strcpy( sTemp, sContainerName((CONTAINER) uUnit));
>>>> +//strcpy( sTemp, sContainerName((CONTAINER) uUnit));
>>>>
>>>> iCount = 0;
>>>> lTemp = lLoop((OBJEKT) uUnit, BONDS);
>>>> .. -302,7 +302,7 ..
>>>>
>>>> iBondCount = iCount;
>>>>
>>>> -fprintf(fOut, "%s\n", sTemp) ;
>>>> +fprintf(fOut, "molecule_name\n") ;
>>>> fprintf(fOut, "%5d %5d %5d 0 1 \n",
>>>> iAtomCount,iBondCount,iResidueCount) ;
>>>> fprintf(fOut, "SMALL\n");
>>>> fprintf(fOut, "USER_CHARGES\n");
>>>>
>>>>
>>>> -->
>>>>
>>>>> Loading PDB file: ./test.pdb
>>>>> [...]
>>>>> Writing mol2 file: out.mol2
>>>>> [...]
>>>>> Loading Mol2 file: ./out.mol2
>>>>> Reading MOLECULE named molecule_name
>>>>
>>>>
>>>>
>>>> Maybe it's useful.
>>>>
>>>> Jan-Philip
>>>>
>>>>
>>>> _______________________________________________
>>>> AMBER mailing list
>>>> AMBER.ambermd.org
>>>> http://lists.ambermd.org/mailman/listinfo/amber
>>>
>>> --
>>> M. L. Dodson
>>> Business email: activesitedynamics-at-gmail-dot-com
>>> Personal email: mldodson-at-comcast-dot-net
>>> Gmail: mlesterdodson-at-gmail-dot-com
>>> Phone: eight_three_two-five_63-386_one
>>>
>>>
>>>
>>> <test.mol2><test.log><test.in>
>>>
>>>
>>> _______________________________________________
>>> AMBER mailing list
>>> AMBER.ambermd.org
>>> http://lists.ambermd.org/mailman/listinfo/amber
>>
>> --
>> M. L. Dodson
>> Business email: activesitedynamics-at-gmail-dot-com
>> Personal email: mldodson-at-comcast-dot-net
>> Gmail: mlesterdodson-at-gmail-dot-com
>> Phone: eight_three_two-five_63-386_one
>>
>>
>>
>>
>>
>>
>>
>> _______________________________________________
>> AMBER mailing list
>> AMBER.ambermd.org
>> http://lists.ambermd.org/mailman/listinfo/amber
>


_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Thu Nov 10 2011 - 15:00:02 PST
Custom Search