I would recommend to use the DOS analyzer in VNL, then you don't need to worry about the Python code... Sure, it only gives tetrahedron spectrum, but why, exactly, do you need the Gaussian version of the DOS? The tetrahedron method is generally better for bulk materials unless you already have a very high k-point sampling (in which case it's also good).
Still, as long as you have the DOS object you can plot it using either tetrahedron or Gaussian. You only need to compute the DOS once. But you can't just take code from two scripts and put them together, you need to understand how the code works. Of course you don't need any molecule configuration if you are working with a bulk configuration - it's just that a molecule is used as an example of how to use the ProjectionList class.
What you need is to take the Gaussian spectrum example from the DOS manual page and just add a ProjectionList defining what you want to project on. Like
dos_g = dos.gaussianSpectrum(energies, broadening = 0.2*eV, projection_list=ProjectionList(elements=[Silicon]))
or use atom indices, or angular momenta, depending on what you want to project on.