pro pdf2d, mass=massw, x=xd, y=yd
@common_blocks.inc
xyimage = [800,800]/2
mad = max(d, min=mid)
mav = max(v, min=miv)
binsize1 =  (mad-mid)/xyimage[0]
binsize2 =  (mav-miv)/xyimage[1]
if N_elements(binsize2) ne 1 then binsize2=.1
if not keyword_set(massw) then begin
   case 1 OF
      (keyword_set(xd) and not keyword_set(yd)): pd = hist2d(xd, v, vol, $
                        binsize1=binsize1, binsize2=binsize2, obin1=obin1,obin2=obin2) 
      (keyword_set(xd) and keyword_set(yd)): pd = hist2d(xd, yd, vol, $
                        binsize1=binsize1, binsize2=binsize2, obin1=obin1,obin2=obin2) 
      (not keyword_set(xd) and keyword_set(yd)): pd = hist2d(d, yd, vol, $
                        binsize1=binsize1, binsize2=binsize2, obin1=obin1,obin2=obin2) 
      else: pd = hist2d(d, v, vol, $
                        binsize1=binsize1, binsize2=binsize2, obin1=obin1,obin2=obin2) 
      endcase
   title = 'Volume weighted pdf. binsize:'+string(binsize1)+string(binsize2)
end  else begin
   case 1 OF
      keyword_set(xd) and not keyword_set(yd): pd = hist2d(xd, v, mass, $
                        binsize1=binsize1, binsize2=binsize2, obin1=obin1,obin2=obin2) 
      keyword_set(xd) and keyword_set(yd): pd = hist2d(xd, yd, mass, $
                        binsize1=binsize1, binsize2=binsize2, obin1=obin1,obin2=obin2) 
      not keyword_set(xd) and keyword_set(yd): pd = hist2d(d, yd, mass, $
                        binsize1=binsize1, binsize2=binsize2, obin1=obin1,obin2=obin2) 
      else: pd = hist2d(d, v, mass, $
                        binsize1=binsize1, binsize2=binsize2, obin1=obin1,obin2=obin2) 
   endcase
   title = 'Mass weighted pdf. binsizes:'+string(binsize1)+string(binsize2)
endelse
ind = where(pd gt 0.)   ; find second minimum
minsec = min(pd[ind])  
;pd[ind] = minsec

iimage, alog10(pd > minsec), TITLE=title, view_zoom=.75,  $
        image_dimensions=[max(obin1)-min(obin1), max(obin2)-min(obin2)], $
        image_location=[min(obin1), min(obin2)], /insert_colorbar, $
        /disable_splash_screen, Identifier=idimage

   idTool = ITGETCURRENT(TOOL=oTool)
   idPlot = oTool->FindIdentifiers('*X AXIS', /OPERATIONS)
   void = oTool -> DoAction(idPlot)
   idPlot = oTool->FindIdentifiers('*Y AXIS', /OPERATIONS)
   void = oTool -> DoAction(idPlot)

return
end