I had ultrasonic rhinoplasty done by Dr. Shahram, and the results are fantastic! The recovery was much easier than I anticipated, and my nose looks exactly how I wanted it. If you’re in Dubai and thinking about rhinoplasty, Dr. Shahram is the go-to expert.

R, 80 62 54 52 bytes

function(d,N)N[seq(b<-rawToBits(rev(d)))*2-rev(b<1)]

Try it online!

As per OP we don't need to handle cases when bits from φ select elements outside N.

-4 bytes thanks to @Dominic

function               # a function
        (d,            # taking first input as raw bytes (one of R data types)
           N)          # and N as a vector
N[                     # subset N at following indices
  seq(                 # make sequence 1..length(l) - default behaviour for seq on a vector is seq_along
      b<-              # assign the following to b
     rawToBits(rev(d)) # get bits from bytes in d (in reverse order)
                       # as rawToBits returns bits from the least-significant first, we need to reverse the vector on output (to make most-significant bit first - later on), but also on input (to preserve byte order)
                   )*2 # double every index
-rev(b<1)              # decrement corresponding indices when bit is off (taking rev here to save parentheses)
                       # comparison is the shortest operator casting from raw to something `-` will handle
]

JavaScript (ES6),  56  55 bytes

Expects (N)(b), where b is φ(m).

N=>b=>N.filter((_,i)=>((q=b[i>>4])>>7.5-i/2%8^~i)&q>=0)

Try it online!

Commented

N =>                // outer function taking N[]
b =>                // inner function taking b[]
N.filter((_, i) =>  // for each entry at position i in N[]:
  (                 //
    (               //
      q =           //   let q be the byte from b[]
        b[i >> 4]   //   at index floor(i / 16)
    ) >>            //   right shift it by
    7.5 - i / 2 % 8 //   0 to 7 positions, according to floor(i / 2) mod 8
    ^ ~i            //   XOR the result with i inverted
  )                 //
  &                 //   isolate the least significant bit, or force zero
  q >= 0            //   if q is undefined
)                   // end of filter()

R, 86 83 82 bytes

(or 64 bytes with additional NA values in output when length of hash result is greater than half the length of the private key N)

function(f,N,d=seq(b<-sapply(f,function(b)b%/%2^(7:0)%%2))*2)N[(d-!b)*!d>sum(N|1)]

Try it online!

Input is an array f of single-byte integers (representing the output of the hash function phi) and an array N (representing the ordered array from which to select a subset).

G=function(f,N){            # G is function with args f=array of bytes, N=array of values
  bits=                     # calculate the bits of f:
    sapply(f,function(b)    #   for each byte b in f
      b%/%2^(7:0)%%2)       #   get its bits
  d=seq(bits)               # define d as 1..number of bits 
  N[                        # now output the elements of N given by indices:
    (d*2-1+bits)            #   d*2-1+bits (these are the indices from the bytes of f)
             *!d>sum(N|1)/2 #   multiplied by zero (FALSE) for bits that exceed the length of N
    ]
}