def _sampleV_preoptimized(self,R,z,maxVT):
        """
        NAME:

           _sampleV_preoptimized

        PURPOSE:

           sample a radial, azimuthal, and vertical velocity at R,z;
           R,z can be an array of positions maxVT is already optimized

        INPUT:

           R - Galactocentric distance (can be Quantity)

           z - height (can be Quantity)
           
           maxVT - an array of pre-optimized maximum vT at corresponding R,z

        OUTPUT:

           a numpy array containing the sampled velocity, (vR, vT, vz),
           where each row correspond to the row of (R,z)

        HISTORY:

           2018-08-09 - Written - Samuel Wong (University of Toronto)

        """
        length = numpy.size(R)
        out= numpy.empty((length,3)) #Initialize output
        #Determine the maximum of the velocity distribution
        maxVR= numpy.zeros(length)
        maxVz= numpy.zeros(length)
        logmaxVD= self(R,maxVR,maxVT,z,maxVz,log=True,use_physical=False)
        #Now rejection-sample
        #Intiialize boolean index of position remaining to be sampled
        remain_indx = numpy.full(length,True)
        while numpy.any(remain_indx):
            nmore= numpy.sum(remain_indx)
            propvR= numpy.random.normal(size=nmore)*2.*self._sr
            propvT= numpy.random.normal(size=nmore)*2.*self._sr+maxVT[remain_indx]
            propvz= numpy.random.normal(size=nmore)*2.*self._sz
            VDatprop= self(R[remain_indx],propvR,propvT,z[remain_indx],propvz,
                           log=True, use_physical=False)-logmaxVD[remain_indx]
            VDatprop-= -0.5*(propvR**2./4./self._sr**2.+
                             propvz**2./4./self._sz**2.+
                             (propvT-maxVT[remain_indx])**2./4./self._sr**2.)
            accept_indx= (VDatprop > numpy.log(numpy.random.random(size=nmore)))
            vR_accept= propvR[accept_indx]
            vT_accept= propvT[accept_indx]
            vz_accept= propvz[accept_indx]
            #Get the indexing of rows of output array that need to be updated
            #with newly accepted velocity
            to_change= numpy.copy(remain_indx)
            to_change[remain_indx]= accept_indx
            out[to_change]= numpy.stack((vR_accept,vT_accept,vz_accept), axis = 1)
            #Removing accepted sampled from remain index
            remain_indx[remain_indx]= ~accept_indx
        return out
def _sampleV_preoptimized ( self , R , z , maxVT ) : length = numpy . size ( R ) out = numpy . empty ( ( length , 3 ) ) maxVR = numpy . zeros ( length ) maxVz = numpy . zeros ( length ) logmaxVD = self ( R , maxVR , maxVT , z , maxVz , log = True , use_physical = False ) remain_indx = numpy . full ( length , True ) while numpy . any ( remain_indx ) : nmore = numpy . sum ( remain_indx ) propvR = numpy . random . normal ( size = nmore ) * 2. * self . _sr propvT = numpy . random . normal ( size = nmore ) * 2. * self . _sr + maxVT [ remain_indx ] propvz = numpy . random . normal ( size = nmore ) * 2. * self . _sz VDatprop = self ( R [ remain_indx ] , propvR , propvT , z [ remain_indx ] , propvz , log = True , use_physical = False ) - logmaxVD [ remain_indx ] VDatprop -= - 0.5 * ( propvR ** 2. / 4. / self . _sr ** 2. + propvz ** 2. / 4. / self . _sz ** 2. + ( propvT - maxVT [ remain_indx ] ) ** 2. / 4. / self . _sr ** 2. ) accept_indx = ( VDatprop > numpy . log ( numpy . random . random ( size = nmore ) ) ) vR_accept = propvR [ accept_indx ] vT_accept = propvT [ accept_indx ] vz_accept = propvz [ accept_indx ] to_change = numpy . copy ( remain_indx ) to_change [ remain_indx ] = accept_indx out [ to_change ] = numpy . stack ( ( vR_accept , vT_accept , vz_accept ) , axis = 1 ) remain_indx [ remain_indx ] = ~ accept_indx return out
NAME:

           _sampleV_preoptimized

        PURPOSE:

           sample a radial, azimuthal, and vertical velocity at R,z;
           R,z can be an array of positions maxVT is already optimized

        INPUT:

           R - Galactocentric distance (can be Quantity)

           z - height (can be Quantity)
           
           maxVT - an array of pre-optimized maximum vT at corresponding R,z

        OUTPUT:

           a numpy array containing the sampled velocity, (vR, vT, vz),
           where each row correspond to the row of (R,z)

        HISTORY:

           2018-08-09 - Written - Samuel Wong (University of Toronto)