Hi list,
I could not find how to Generate a density map from a structure with a *fixed* dimension :
When using IMP.em.SampledDensityMap(ps,resolution,apix) with a fixed size for apix, the grid size depends on the resolution. I would like to write a function that would allow me to generate maps with the same number of vertices (and the same space coverage) for different resolutions. I thus tried :
def pdbFile2dmap (pdbStream,densityFile,resolution,apix,sel,bboxCoverage=1.2):
'''
pdbStream : where to read the structure
densityFile : name of the file where to drop the density map
resolution : resolution of the density map
apix : voxel size
sel : selector for reading structure
bboxCoverage : approximate bbox coverage of the density map
'''
m= IMP.Model()
# read protein
mh = IMP.atom.read_pdb(pdbStream,m,sel)
# add radius info to each atom
IMP.atom.add_radii(mh)
# compute bbox, and map size in voxels
bbox = IMP.atom.get_bounding_box(mh)
diag = bbox.get_corner(0) - bbox.get_corner(1)
nx = int(bboxCoverage * diag[0] / apix)
ny = int(bboxCoverage * diag[1] / apix)
nz = int(bboxCoverage * diag[2] / apix)
# compute origin
origin = bbox.get_corner(0) + (1-bboxCoverage)/2 * diag
dmap = IMP.em.SampledDensityMap()
# dmap.CreateVoidMap(nx,ny,nz)
ps = IMP.Particles(IMP.core.get_leaves(mh))
dmap.set_particles(ps)
dmap.set_origin(origin)
dmap.get_header_writable().set_number_of_voxels(nx,ny,nz)
dmap.get_header_writable().set_resolution(resolution)
dmap.update_voxel_size(apix)
dmap.resample()
dmap.calcRMS() # computes statistic stuff about the map and insert it in the header
print dmap.get_header().show(),"\n"
IMP.em.write_map(dmap,densityFile,IMP.em.MRCReaderWriter())
that resulted in :
python(10802) malloc: *** mmap(size=18446744073709314048) failed (error code=12)
*** error: can't allocate region
*** set a breakpoint in malloc_error_break to debug
Traceback (most recent call last):
File "./pdb2density.py", line 181, in <module>
pdbFile2dmap(pdbStream,densityFile,resolution,apix,sel,1.2)
File "./pdb2density.py", line 75, in pdbFile2dmap
dmap.update_voxel_size(apix)
File "/usr/local/lib/python2.6/site-packages/IMP/em/__init__.py", line 428, in update_voxel_size
def update_voxel_size(self, *args): return _IMP_em.DensityMap_update_voxel_size(self, *args)
MemoryError: std::bad_alloc
Can you tell me if my approach is OK, and how to solve my allocation problem ?
I would also like to know if there are easy means in IMP to
1.) Strip a density map : out of an existing map, create a new map or restrict the existing one, so as to conserve only a central subset of voxels
2.) Re-sample a map : Conserve the density map "box" size, and change the number of voxels that cover it
--Ben