Mutational signatures analysis is a powerful technique for uncovering the mutational processes involved in cancer. Current approaches are based on non-negative matrix factorization (NMF), however, this ignores the non-homogeneous occurrence of mutations across the genome. We introduce a flexible new Bayesian method using Poisson point processes to model the activities of mutational signatures as they vary across the genome. Using covariate-dependent factorized intensity functions, our Poisson process factorization (PPF) model generalizes the standard NMF model to include regression coefficients that capture the effect of genomic features on the mutation rates from each latent process. Furthermore, our method employs sparsity-inducing hierarchical priors to automatically infer the number of active latent factors in the data. We present algorithms to obtain maximum a posteriori estimates and uncertainty quantification via Markov chain Monte Carlo. We demonstrate the method on simulated data and on real data from breast cancer.