(* Independent measurements, all with the same error deviation *)
(Subscript[t, real])=23;
a=Table[1,{i,3},{j,1}]
(Subscript[b, real])=Table[(Subscript[t, real]),{i,3},{j,1}]
V=DiagonalMatrix[{1,1,1}]

errors=RandomVariate[NormalDistribution[0,1],{3,1}]

bhat=(Subscript[b, real])+errors
P=Inverse[Transpose[a].V.a] (* Kovariancni matice t-(Subscript[x, guess]) *)
moorePenrose[x_]:=P.Transpose[a].V.x
(Subscript[x, guess])=moorePenrose[bhat][[1,1]]

(* Probability distribution of the calculated t, i.e. normal distribution with expected value (Subscript[x, guess]) and variance P. . Red line marks the real value for comparison. *)
Plot[PDF[NormalDistribution[(Subscript[x, guess]),P[[1,1]]],x]//Evaluate,{x,21,24},Filling->Axis,PlotRange->All,GridLines->{{{(Subscript[t, real]),{Red,Thick}}},None}]

(* Independent measurements, all with the first two receivers having the same error variance and the third one being more faulty *)
(Subscript[t, real])=23;
a=Table[1,{i,3},{j,1}];
(Subscript[b, real])=Table[(Subscript[t, real]),{i,3},{j,1}];
V=DiagonalMatrix[{1,1,1/4}];

errors=Join[RandomVariate[NormalDistribution[0,1],2],RandomVariate[NormalDistribution[0,4],1]]

bhat=(Subscript[b, real])+errors
P=Inverse[Transpose[a].V.a] (* Kovariancni matice t-(Subscript[x, guess]) *)
moorePenrose[x_]:=P.Transpose[a].V.x
(Subscript[x, guess])=moorePenrose[bhat][[1,1]]

(* Probability distribution of the calculated t, i.e. normal distribution with expected value (Subscript[x, guess]) and variance P. Red line marks the real value for comparison. *)
Plot[PDF[NormalDistribution[(Subscript[x, guess]),P[[1,1]]],x]//Evaluate,{x,21,25},Filling->Axis,PlotRange->All,GridLines->{{{(Subscript[t, real]),{Red,Thick}}},None}]

(* Even the faulty measurements count *)
(Subscript[t, real])=23;
a=Table[1,{i,3},{j,1}];
(Subscript[b, real])=Table[(Subscript[t, real]),{i,3},{j,1}];
V=DiagonalMatrix[{1,1,1/25}];

errors=Join[RandomVariate[NormalDistribution[0,1],2],RandomVariate[NormalDistribution[0,5],1]]

bhat=(Subscript[b, real])+errors
P=Inverse[Transpose[a].V.a] (* Kovariancni matice t-(Subscript[x, guess]) *)
moorePenrose[x_]:=P.Transpose[a].V.x
(Subscript[x, guess])=moorePenrose[bhat][[1,1]]

bhat2=((Subscript[b, real])+errors)[[1;;2]]
P2=Inverse[Transpose[a[[1;;2]]].IdentityMatrix[2].a[[1;;2]]] (* Kovariancni matice t-(Subscript[x, guess]) *)
moorePenrose[x_]:=P.Transpose[a[[1;;2]]].IdentityMatrix[2].x
(Subscript[x2, guess])=moorePenrose[bhat2][[1,1]]

(* Blue distribution is estimated from all the measurements (i.e. even the faulty one), the green one uses low-error measurements only. Red line marks the real value for comparison. *)
g1=Plot[PDF[NormalDistribution[(Subscript[x, guess]),P[[1,1]]],y]//Evaluate,{y,21,24},Filling->Axis,PlotRange->All,PlotStyle->Blue,GridLines->{{{(Subscript[t, real]),{Red,Thick}}},None}];
g2=Plot[PDF[NormalDistribution[(Subscript[x2, guess]),P2[[1,1]]],y]//Evaluate,{y,21,24},Filling->Axis,PlotRange->All,PlotStyle->Green];
Show[g1,g2]

ListPlot