What is the original derivation of the Poisson distribution? I am learning the Poisson distribution. I understand it, but its probability mass function is not natural to me. I think its probability mass function seems to be derived from somewhere with more contexts (or from one application). Where is the original source?

I think it first appeared in Poisson's "Recherches sur la probabilité des jugements en matière criminelle et en matière civile'  (researches into the probabilities of judgements in criminal and civil cases).
He did not define it through the PMF but as a limit of a Binomial distribution (using the law of large numbers which he calls "Bernoulli's theorem")
At some point he obtains the formula
$$
P = \left(1 + \omega + \frac{\omega^2}{1 \cdot 2} + \dots + \frac{\omega^n}{1\cdot 2  \cdots n} \right) e^{-\omega}
$$
for the probability of obtaining at most $n$ successess out of $\mu$ trials in which there is a probability $\frac{\omega}{\mu}$ of success (for $\mu$ "very large").
This is the cumulative distribution function of a Poisson distribution with parameter $\omega$ evaluated at $n$.
You can see the formula above here page 206 (it's in French).
