Vjezba 2: Koeficijent apsorpcije

http://sail.zpf.fer.hr/laserlab/beamabs.nb

1. Unos podataka

Unose se parovi {(Subscript[x, i]), (Subscript[P, i])} gdje je (Subscript[x, i]) u cm, a (Subscript[P, i]) u mW:

data = {{2, 0.38}, {2.5, 0.355}, {5, 0.24}, {6, 0.20}, {8, 
          0.14}}

n = Length@data

datafig = {{Red, AbsolutePointSize[5], 
            Point /@ data}};
Show[Graphics@datafig, Axes -> True, 
        AxesLabel -> {"x [cm]", "P [mW]"}];

2. Linearizirani model: (Subscript[y, i]) = ln (Subscript[P, i]) = ln (Subscript[P, 0]) - α (Subscript[x, i]) = 
            a + b (Subscript[x, i])

Pogledati upute 2.4.1 i A.5

Funkcija koja racuna redak tablice:

redakfun[{x_, p_}] := 
      Module[{y = Log[p]}, {x, p, y, (x^2), x y, (y^2)}]

(tablica = redakfun /@ data) // TableForm

Zbroj vrijednosti u stupcima tablice:

{sx, tmp, sy, sx2, sxy, sy2} = Total /@ Transpose@tablica

Velicina Δ prema izrazu (A.19):

del = n sx2 - (sx^2)

Ocekivana vrijednost parametra a prema izrazu (A.17):

a = (sx2 sy - sx sxy)/del

Ocekivana vrijednost parametra b prema izrazu (A.18):

b = (n sxy - sx sy)/del

Standardna pogreska parametra b prema izrazu (A.23):

db = 
      Sqrt[(n/(del (n - 2))) (sy2 - a sy - b sxy)]

Koeficijent apsorpcije je, u lineariziranom modelu, konacno:

alphafitlin = -b ± db

Rel. pogreska:

100  Abs[db/b]

\(\(Plot[
        Exp[a + b\ x], {x, 
          0, \((3/2)\) 
            Max@\(First@\(Transpose@
                  data\)\)}, 
PlotRange -> {\(-0.05\) 
              Exp@a, 1.05 
               Exp@a}, 
AxesLabel -> {"x [cm]", "P [mW]"}, 
Epilog -> datafig];\)\)

3. Nelinearni model: (Subscript[P, i]) = (Subscript[P, 0]) Exp[-α (Subscript[x, i])]

Pogledati upute 2.4.2 i A.4

3.1 Prilagodbeni postupak (fit)

Za pretpostavljenu vrijednost α funkcija p0fun racuna ocekivanu vrijednost (Subscript[P, 0]) prema formuli (2.12) iz uputa:

p0fun[α_, data_] := 
      Module[{xi, pi, tmp}, 
{xi, pi} = 
          Transpose@data; 
tmp = 
          Exp[-α xi]; 
pi . tmp/
          tmp . tmp
]

Funkcija s2fun racuna (S^2) prema formuli (2.11) uz vrijednost (Subscript[P, 0]) dobivenu funkcijom p0fun:

s2fun[α_, data_] := 
      Module[{xi, pi, tmp}, 
{xi, pi} = 
          Transpose@data; 
tmp = 
          pi - p0fun[α, data] 
              Exp[-α xi]; 
tmp . 
          tmp
]

Crtamo (S^2) u ovisnosti o parametru α:

\(\(Plot[s2fun[α, data], {α, 0, 0.5}, 
        AxesLabel -> {"α", \*""\!\(S\^2\)""}, 
        PlotPoints -> 64];\)\)

Pronalazimo minimum (S^2) u odnosu na parametar α:

s2min = 
      FindMinimum[s2fun[α, data], {α, 0, 0.5}]

Ocekivane vrijednosti parametara:

{alphahat, p0hat} = {α, p0fun[α, data]} /. 
        Last[s2min]

\(\(Plot[
        p0hat\ Exp[\(-alphahat\)\ x], {x, 
          0, \((3/2)\) 
            Max@\(First@\(Transpose@
                  data\)\)}, 
PlotRange -> {\(-0.05\) 
              p0hat, 1.05  
              p0hat}, 
AxesLabel -> {"x [cm]", "P [mW]"}, 
Epilog -> datafig];\)\)

3.2 Procjena pogreske (nije opisano u uputama)

Napomena: postupak procjene pogreske nije opisan u uputama.

Definira se reducirani (χ^2):

chi2redfun[α_, 
          data_] = (n - 2) s2fun[α, data]/
            s2fun[alphahat, data];

Druga derivacija po α reduciranog (χ^2) u minimumu:

fdd =  
      D[chi2redfun[α, data], {α, 2}] /. {α -> 
            alphahat}

Standardna ((χ^2)-curvature) pogreska parametra α:

sigalph = 2/fdd // Sqrt

Kontrolni graf:

\(\(Plot[\((n - 2)\)\ s2fun[α, data]/First[s2min], {α, 
          alphahat - 1.5  sigalph, alphahat + 1.5  sigalph}, 
        FrameLabel -> {"α", \*""\!\(S\^2\)""}, 
        PlotRange -> {n - 2 - 0.5, n - 2 + 1.5}, Frame -> True, 
        Axes -> None, 
        GridLines -> {{alphahat, alphahat - sigalph, 
              alphahat + sigalph}, {n - 2, n - 2 + 1}}];\)\)

Konacno apsorpcijski koeficijent:

alphafitnonlin = alphahat ± sigalph

Rel. pogreska:

100  Abs[sigalph/alphahat]