Using decays to ϕ-meson pairs, the inclusive production of charmonium states in b-hadron decays is studied with pp collision data corresponding to an integrated luminosity of 3.0fb-1, collected by the LHCb experiment at centre-of-mass energies of 7 and 8 TeV. Denoting by BC≡B(b→CX)×B(C) the inclusive branching fraction of a b hadron to a charmonium state C that decays into a pair of ϕ mesons, ratios RC2C1≡BC1/BC2 are determined as Rηc(1S)χc0=0.147±0.023±0.011, Rηc(1S)χc1=0.073±0.016±0.006, Rηc(1S)χc2=0.081±0.013±0.005, Rχc0χc1=0.50±0.11±0.01, Rχc0χc2=0.56±0.10±0.01 and Rηc(1S)ηc(2S)=0.040±0.011±0.004. Here and below the first uncertainties are statistical and the second systematic. Upper limits at 90% confidence level for the inclusive production of X(3872), X(3915) and χ c 2(2 P) states are obtained as Rχc1X(3872)<0.34, Rχc0X(3915)<0.12 and Rχc2χc2(2P)<0.16. Differential cross-sections as a function of transverse momentum are measured for the ηc(1 S) and χc states. The branching fraction of the decay Bs0 is measured for the first time, B(Bs0)=(2.15±0.54±0.28±0.21B)×10-6. Here the third uncertainty is due to the branching fraction of the decay Bs0, which is used for normalization. No evidence for intermediate resonances is seen. A preferentially transverse ϕ polarization is observed. The measurements allow the determination of the ratio of the branching fractions for the ηc(1 S) decays to ϕϕ and pp¯ as B(ηc(1S))/B(ηc(1S)→pp¯)=1.79±0.14±0.32.