In this paper, we consider the following nonlinear Choquard equation with magnetic field \(\begin{aligned} \begin{aligned} \left\{ \begin{array}{l} \displaystyle \bigg (\frac{\varepsilon }{i}\nabla -A(x)\bigg )^{2}u+V(x)u=\varepsilon ^{\mu -N}\left( \,\,\int \limits _{{\mathbb {R}}^{N}}\frac{|u(y)|^{2_{\mu }^{*}}+F(|u(y)|^{2})}{|x-y|^{\mu }}\text {d}y\right) \left( |u|^{2_{\mu }^{*}-2}u+\frac{1}{2_{\mu }^{*}}f(|u|^{2})u\right) \hspace{1.14mm}\text{ in }\hspace{1mm} {\mathbb {R}}^{N},\\ \displaystyle u\in H^{1}({\mathbb {R}}^{N},{\mathbb {C}})\\ \end{array} \right. \end{aligned} \end{aligned}\) where \(\varepsilon >0\) is a small parameter, \(N\ge 3\) , \(0<\mu <N\) , \(2_{\mu }^{*}=\frac{2N-\mu }{N-2}\) , \(V(x):{\mathbb {R}}^{N}\rightarrow {\mathbb {R}}^{N}\) and \(A(x):{\mathbb {R}}^{N}\rightarrow {\mathbb {R}}^{N}\) is a continuous potential, f is a continuous subcritical term, and F is the primitive function of f. Under a local assumption on the potential V, by the variational methods, the penalization techniques and the Ljusternik–Schnirelmann theory, we prove the multiplicity and concentration properties of nontrivial solutions of the above problem for \(\varepsilon >0\) small enough.