How to implement a M-H step in a Gibbs sampling

Question

I am having trouble implementing a Metropolis Hastings step in a Gibbs sampling problem. The following code was taken from https://www.stat.colostate.edu/computationalstatistics/

Details: It is a capture recapture study, with seven total draws from the population. ci are the number drawn at time i, mi is the number of unique individuals drawn at time i.

The goal is to estimate the total population, N , using gibbs and or MH sampling.

As in, from the original problem , a likelihood is used is an "M(t)" model , that is,

$L(N, \alpha |c,r) \alpha \frac{N!}{(N-r)!} \Pi_{i=1}^{I=7}\alpha_{i}^{c_{i}}(1-\alpha_{i})^{N-c_{i}}$

From this it can be derived that N has a scaled negative binomial type distribution.

Thus Gibbs sampling could be implemented as follows:

ci = c(30,22,29,26,31,32,35)  
mi = c(30,8,17,7,9,8,5)
r = sum(mi)
I=7
    num.its=100000

     alpha= matrix(0,num.its,I)
     N=rep(0,num.its)
     N[1]=sample(84:500,1) 
     set.seed(4)


    for (i in 2:num.its) {
       alpha[i,]=rbeta(I,ci+.5,N[i-1]-ci+.5)
        N[i]=rnbinom(1,r+1,1-prod(1-alpha[i,]))+r
  }

My question is, I want to understand how I could use a Metropolis-Hastings step instead if I did not know the form of N being a negative binomial. For example, using a proposal density of Poisson from the previous value.

How could this be implemented?

So could I use the above as my 'f' as the full likelihood in MH? and my 'g' could be Poisson of the previous N?

What I have tried, but is not working is

alpha[i,]=rbeta(I,ci+.5,N[i-1]-ci+.5)


      xo<-N[i-1]
      xn<-rpois(1,lambda=xo)

     s1<-factorial(xn)/factorial(xo)




   a=s1*(prod(1-alpha[i,]))^(xn-xo)                  




   u<-runif(1,0,1)

      if(u<=a){
       N[i]=xn
     }
     else if(u>a){
     N[i]=xo
   }

Because I am having trouble using the whole factorial terms since the numbers are so large I get errors. In the way I have it it only decreases quickly to zero instead of evening out the correct final value of around 90.

Another issue I have is that, when we subtract 84 from the factorial, I get errors because it is possible that the new N is less than 84.