I'm afraid this should really be a comment not an answer, since I am giving a reference rather than a description. But a full description is a bit long.

The description you want is indeed known, and is very nicely written up in the later chapters of the book `Néron models' by Bosch et al. For example, theorem 9.5.4 relates the Néron model to a certain quotient of the relative Picard. Theorem 9.6.1 tells you what the group of connected components is. The identity component itself is also discussed in 9.5.4. The assumption that the singularities be e.g. nodal is not needed, the results are very general. The main thing you want is something like the geometric multiplicities of the components of the special fibre having gcd = 1, which is for example satisfied in the nodal case, or more generally if there exists an étale quasi-section.

I'm afraid I don't quite follow your paragraph about loops, but I hope you can resolve things from the reference - please ask if not. Certainly there are close connections between the combinatorics of the graph and the component group in the nodal case; this can be extracted quite easily from the references given above, though perhaps there is somewhere else it is written more simply. Roughly, the toric rank is the rank of the homology of the graph (see e.g. example 9.2.8 of [loc.cit.]). The component group is the `critical group' of the graph. This goes by many names, you can extract a definition from lemma 9.5.9. I like to think of it in terms of electrical resistance in a circuit, but the people who referee my papers seem to disagree ;-).

In the special case of modular curves, you might also for example be interested in the paper `On Néron models, divisors and modular curves' by
Bas Edixhoven, which makes some parts more explicit.

Thanks to Owen Biesel for pointing this question out to me.

**Added later**

Now to try to allow for this generalised jacobian on the generic fibre. Some parts of the above description go through OK, but for others I get a bit stuck, as we shall see.

First, I find it helpful to relate the generalised jacobian to the relative Picard of another curve. You have the curve $\mathcal C/Z_p$, and a finite set $S$ of smooth sections. Let $\mathcal C_S$ be the curve obtained by `glueing these sections together' - more formally this should be constructed as a pushout, see for example [Ferrand, CONDUCTEUR, DESCENTE ET PINCEMENT] for a nice general theory. This $\mathcal C_S/Z_p$ has non-smooth locus consisting of exactly one section (with some multiplicity) and a bunch of points lying over thickenings of the closed point in $Z_p$ (the non-smooth locus is always finite unramified in the ordinary double point case). In particular, the generic fibre is irreducible, with exactly one node. Write $P$ for the total-degree-zero part of the relative Picard of $C_S/Z_p$. Then the generic fibre of $P$ is just the degree zero part of the relative Pic of the generic fibre of $C_S$, in particular it is a smooth connected group scheme over $Q_p$. In fact, it is exactly the generalised jacobian in the sense of Serre --- this can be deduced from chapter V of [Serre, Algebraic groups and class fields].

OK, so we have this group scheme P, whose generic fibre is the generalised jacobian. I want to relate $P$ to the N\'eron model. Note that $P$ is separated over $Z_p$ if and only if the special fibre of $C_S$ (or equivalently of $C$ is irreducible. So in general it is not separated, so is certainly not the N\'eron model. However, this is quite easy to fix - one just takes the largest separated quotient (equivalently, divide out by the closure of the unit section), to obtain a separated group scheme, let's call it $Q$. Is this $Q$ now the N\'eron model? Note that it is not obvious, since $C_S$ is not regular, so the usual argument about extending line bundles fails.

In general I do not know the answer. However, in the nodal case (where $S$ has cardinality 2) the answer is indeed yes, $Q$ is the N\'eron model *added: if one makes some suitable blowups of he special fibre first. *. This follows from a recent result of Giulio Orecchia, [Semi-factorial nodal curves and Néron models of jacobians, https://arxiv.org/abs/1602.03700 ]. The point is that all the edges in the dual graph have thickness 1 or infinity (since you started with a regular model $\mathcal C$), so Orecchia's `circuit coprime' condition is automatically satisfied once one makes some blowups to avoid degenerate loops. Note that of course $Q$ is actually a N\'eron lft-model, since it need not be of finite type.

The intersection of the closure of the unit section with the fibrewise connected component of identity will be trivial, so the identity component of the N\'eron model $Q$ will just be Pic^0 of the special fibre of $\mathcal C_S$, which has a nice explicit description as discussed above.

Even in the nodal case, this leaves the question of what the component group is. It should be possible to read this off from the dual graph of the special fibre together with a labelling of the edges by the thickness of the singularity there. Some variation of the constructions in [BLR, Neron models] should suffice for this, but I have never thought about how to do it.

It seems to me that the more difficult point is to relate this largest separated quotient $Q$ of $P$ to the N\'eron model. This is what is done by Orecchia in the case where $S$ has two elements. It seems quite possible his methods would generalise to this case, but I am not sure - maybe you should ask him...