Ben/Daniel - 
I agree that initializing a map from a bounding-box makes sense.
Please see revision 1644 that now contains an example on how to do that 
em/examples/generate_density_map_of_fixed_dimension.py

I have added a function that generates a header from a bounding box and show how to use it in the example.

Keren.
On Aug 2, 2010, at 7:57 AM, Daniel Russel wrote:

On Aug 1, 2010, at 10:07 PM, Keren Lasker wrote:

hi Ben - 
You can initiate a SampledDensityMap from a DensityHeader, which you can set to any dimension you want.
Does that make sense ?
What fields of the density header must be set in order for construction to work right?
set_number_of_voxels()
set_resolution()
set_{x,y,z}origin()
Objectpixelsize_ (there doesn't seem to be a set function, right?)
Anything else?


Our Ben - I am still for some reason not getting these messages, only after forwarding from Daniel :)
This message was old, I just revived it.


On Aug 1, 2010, at 12:31 PM, Daniel Russel wrote:

Keren? I have comments below.

On Jul 22, 2010, at 6:51 AM, Benjamin SCHWARZ wrote:

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 ?
There DenistyMap/DensityHeader is still a bit finicky (that is, there are protocols which have to be followed, and, even worse, they don't tend to be documented). Keren is working on a cleaner version. My guess it the problem is due to not setting the "upper right" value or the length of the edge of each voxel. I don't see the function to do that offhand. Keren?


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
There is get_transformed_into() which could be used. This may be a post 1.0 addition.

2.) Re-sample a map : Conserve the density map "box" size, and change the number of voxels that cover it
There is get_resampled(). Also may be a post 1.0 addition.

 
   --Ben
_______________________________________________
IMP-users mailing list
IMP-users@salilab.org
https://salilab.org/mailman/listinfo/imp-users

_______________________________________________
IMP-users mailing list
IMP-users@salilab.org
https://salilab.org/mailman/listinfo/imp-users

_______________________________________________
IMP-users mailing list
IMP-users@salilab.org
https://salilab.org/mailman/listinfo/imp-users

_______________________________________________
IMP-users mailing list
IMP-users@salilab.org
https://salilab.org/mailman/listinfo/imp-users