Normally I wouldn't ask such question before even read the material from the Wolfram Documentation about what I'm trying to do but this time I have a pretty tigh dead end. I apologize for that. So, what I need to do is to extract a data from a NDSolve integration to use it in different ways. I will show my whole code and make a precise description for what I need. The code is the following:

AU = 1.49597871*10^13;
n = 0.5*10^15; 
cg = 6.67384/10^8; 
kb = 1.3806488/10^16; 
mh = 1.0078250321*(1.660538921/10^24); 
mf = 1.98892*10^33; 
Te[ϕ_] := 2200 + 440*Cos[ϕ]; 
Ra[ϕ_] := 1.9 + 0.1*Cos[ϕ + Pi]; 
t[r_, ϕ_] := Te[ϕ]*(1 - Sqrt[1 - (Ra[ϕ]/r)^2])^(1/4);

k1[r_, ϕ_] := (6.99/10^14)*(t[r, ϕ]/300)^2.8*Exp[-1950/t[r, ϕ]]; 
k2[r_, ϕ_] := (1.59/10^11)*(t[r, ϕ]/300)^1.2*Exp[-9610/t[r, ϕ]]; 
k17[r_, ϕ_] := (3.14/10^13)*(t[r, ϕ]/300)^2.7*Exp[-3150/t[r, ϕ]]; 
k18[r_, ϕ_] := (2.05/10^12)*(t[r, ϕ]/300)^1.52*Exp[-1736/t[r, ϕ]]; 
k62[r_, ϕ_] := (1.77/10^11)*Exp[178/t[r, ϕ]]; 
k63[r_, ϕ_] := (1.85/10^11)*(t[r, ϕ]/300)^0.95*Exp[-8571/t[r, ϕ]]; 
k94[r_, ϕ_] := (1.65/10^12)*(t[r, ϕ]/300)^1.14*Exp[-50/t[r, ϕ]]; 
k138[r_, ϕ_] := (5.94/10^17)*(t[r, ϕ]/300)^0.17*Exp[65.9/t[r, ϕ]]; 
k141[r_, ϕ_] := (1/10^11)*Exp[-4800/t[r, ϕ]]; 
k143[r_, ϕ_] := (1/10^11)*Exp[-4800/t[r, ϕ]]; 
k144[r_, ϕ_] := 1/10^9; 

r1 = k1[r, ϕ]*nH[r]*nOH[r]; 
r2 = k2[r, ϕ]*nH[r]*nH2O[r]; 
r17 = k17[r, ϕ]*nH2[r]*nO[r]; 
r18 = k18[r, ϕ]*nH2[r]*nOH[r]; 
r62 = k62[r, ϕ]*nO[r]*nOH[r]; 
r63 = k63[r, ϕ]*nO[r]*nH2O[r]; 
r94 = k94[r, ϕ]*nOH[r]*nOH[r]; 
r138 = k138[r, ϕ]*nTi[r]*nO[r]; 
r141 = k141[r, ϕ]*nTi[r]*nOH[r]; 
r143 = k143[r, ϕ]*nTi[r]*nH2O[r]; 
r144 = k144[r, ϕ]*nTiO[r]*nH2O[r];

dens[r_, ϕ_] := n*Exp[(-((2*cg*mf*mh)/(kb*Te[ϕ]*Ra[ϕ]*AU)))*
           (1 - (Ra[ϕ]/r)^(1/2))]; 
dldens[r_, ϕ_] := D[dens[r, ϕ], r]/dens[r, ϕ];

eqns = {
   (dldens[r, ϕ]*nH[r] + Derivative[1][nH][r])/dens[r, ϕ] == 
     -r1 + r17 + r18 - r2 + r62, 
   (dldens[r, ϕ]*nOH[r] + Derivative[1][nOH][r])/dens[r, ϕ] == 
     -r1 - r141 + r17 - r18 + r2 - r62 + r63 - r94, 
   (dldens[r, ϕ]*nO[r] + Derivative[1][nO][r])/dens[r, ϕ] == 
      r1 - r138 - r62 - r63 + r94, 
   (dldens[r, ϕ]*nH2[r] + Derivative[1][nH2][r])/dens[r, ϕ] == 
      r1 - r17 - r18 + r2, 
   (dldens[r, ϕ]*nH2O[r] + Derivative[1][nH2O][r])/dens[r, ϕ] == 
     -r143 - r144 + r18 - r2 + r94, 
   (dldens[r, ϕ]*nTi[r] + Derivative[1][nTi][r])/dens[r, ϕ] == 
     -r138 - r141 - r143 - r144, 
   (dldens[r, ϕ]*nTiO[r] + Derivative[1][nTiO][r])/dens[r, ϕ] == 
      r138 + r141 + r143 - r144}; 

ic = {nH[2] == 2/10^7, nOH[2] == 9.8/10^7, nO[2] == 2/10^5, nH2[2] == 5/10, 
      nH2O[2] == 8/10^6, nTi[2] == 1.96/10^7, nTiO[2] == 2/10^8}; 

sol = Table[
   NDSolveValue[{eqns, ic} /. ϕ -> ϕi, {nH, nOH, nO, nH2, nH2O, nTi, nTiO},            
      {r, 2, 6}, MaxStepSize -> 0.001, MaxSteps -> Infinity,
                 AccuracyGoal -> 0, PrecisionGoal -> 0, Method -> {"BDF"}], 
      {ϕi, 0, 2*Pi, Pi/5}];
 

Question: how could I extract the data from the solution to a especific value of r and ϕ?

Let's say, for example, that I want to know nOH[6] for all values of ϕ, so

nOH[6] for ϕ=0
nOH[6] for ϕ=Pi/5
nOH[6] for ϕ=Pi

and so on. Is there some way to do it?

Moreover, is there some way to export the data from the solution where I can obtain an arrangement like the following?

Thank you very much.