mpfit

게시판 IDL Q&A mpfit

태그: , ,

4 답변 글타래를 보이고 있습니다
  • 글쓴이
    • #1144 Reply
      chobo
      회원

      안녕하세요. ^^
      mpfit을 이용하여 lognormal fitting을 하려고 합니다.
      chi-square 계산이 infininte값이 나와 에러가 발생합니다.
      이런 경우에는 어떻게 파라미터를 해줘야 하는지요?
      매번 도움만 받아서 죄송합니다 *^^*

      idl code는 다음과 같습니다.

      x = findgen(100)*0.5 + 0.25
      y = [0.0000000, 0.0000000, 0.0000000, 0.0000000, 0.0142621, 0.0153541, 0.0604877, 0.0610517, 0.0606418, 0.0524113, 0.0450828, 0.0393483, 0.0335146, $ 0.0299068,0.0269496,0.0244353,0.0227490,0.0206914,0.0193604,0.0179925,0.0167907,0.0154459,0.0146801,0.0139811,0.0130469,0.0117275,0.0111692,0.0103157, $
      0.0098059,0.0094023,0.0087956,0.0080920,0.0077276,0.0072754,0.0070978,0.0065996,0.0061705,0.0059837,0.0057553,0.0054508,0.0049526,0.0048165,0.0044843, $ 0.0043159,0.0040160,0.0037992,0.0035247,0.0035109,0.0032294,0.0030126,0.0029065,0.0027358,0.0027266,0.0024820,0.0023875,0.0022329,0.0020392,0.0019953, $
      0.0019307,0.0018615,0.0016562,0.0016424,0.0014925,0.0013771,0.0013771,0.0014210,0.0012133,0.0011534,0.0010980,0.0010150,0.0010403,0.0009642,0.0009665,$ 0.0008420,0.0008627,0.0007866,0.0008004,0.0007474,0.0006320,0.0006966,0.0006390,0.0005628,0.0005075,0.0005628,0.0005236,0.0004959,0.0004613,0.0004291, $
      0.0003668,0.0004360,0.0004775,0.0004291,0.0003345,0.0003368,0.0003506,0.0003345,0.0002699,0.0003045,0.0003183,0.0003137]

      window, 1
      device, decomposed=0
      loadct,39
      plot, x, y, psym=2, color=255, yrange=[0, 0.1]

      ; —- Lognormal fitting from origin s/w
      ; the parameter from origin
      P =fltarr(4)
      P[0] = 0.00228 & P[1] = 0.30216 & p[2] = 0.48267 & p[3] = 5.59144
      y_origin = P[0] + P[1] / (x*P[2]*sqrt(2.*!pi)) * exp(-1 * (alog(x/P[3])) * (alog(x/P[3])) /(2.*P[2]*P[2]))
      oplot, x, y_origin, color = 50, thick = 2

      ;– Lognormal fitting from mfit —————————–
      expr = ‘P[0] + P[1] / (x*P[2]*sqrt(2.*!pi)) * exp(-1 * (alog(x/P[3])) * (alog(x/P[3])) /(2.*P[2]*P[2]))’
      start=[0.001, 0.1, 0.1, 1.0] ; input starting point of fitting
      weights=replicate(1.0, n_elements(x)) ; No weighting
      result=mpfitexpr(expr, x, y, start, weights=weights)
      print, result

      y_fitted = result[0] + result[1] / (x*result[2]*sqrt(2.*!pi)) * exp(-1* (alog(x/result[3]))* (alog(x/result[3])) /(2.*result[2]*result[2]))
      oplot, x, y_fitted,color=100,thick=2

      end

      ; === IDL 결과 ====
      % MPFITEXPR: Number of parameters: 4 (initialized to zero)
      Iter 1 CHI-SQUARE = -NaN DOF = 96
      P(0) = 0.000000
      P(1) = 0.000000
      P(2) = 0.000000
      P(3) = 0.000000
      % MPFITEXPR: ERROR: parameter or function value(s) have become infinite; check model function for over- and underflow
      0.00000000 0.00000000 0.00000000 0.00000000

    • #1145 Reply
      Jonghyuk
      회원

      핵심적인 부분은,

      result=mpfitexpr(expr, x, y, start, weights=weights)

      이 부분인데요, 원래 문법이 이렇습니다.

      result=mpfitexpr(expr, x, y, err, start, weights=, …)

      즉, err에 해당하는 부분이 빠졌고, 그래서 mpfitexpr은 start 변수를 err라고 간주하고, 시작값은 주어지지 않았으므로 모두 0이라고 간주하고 피팅을 시작하게 됩니다.
      weights= 키워드가 사용되면 err는 무시되기 때문에, start값이 err로 잘못 입력된 것은 큰 문제 없이 넘어가지만, 시작값이 [0,0,0,0] 디폴트값으로 되는 것은 첫판에 문제가 생깁니다. 분모에 0이 들어가는 일이 발생하는 것이죠.

      start=[1, 1, 1, 1]
      result=mpfitexpr(expr, x, y, sqrt(y), start, weights=weights)
      와 같이 수정하여 실행하시면 만족스러운 계산 결과를 얻으실 것 같습니다. 물론, Positional Parameter의 특성상 순서를 지켜서 기입해야 하기 때문에, 실제 계산에서는 무시될 err 값을 sqrt(y)로 넣어 주었습니다. 그래야, start 변수를 시작값으로 인식할 테니까요.

      다음과 같이 결과가 나옵니다.
      ….
      Iter 15 CHI-SQUARE = 0.00070653355 DOF = 51
      P(0) = 0.00185790
      P(1) = 0.203406
      P(2) = 0.306719
      P(3) = 4.50020
      0.00185790 0.203406 0.306719 4.50020

      그리고, mpfitexpr 함수에는 yfit 키워드가 있어서, 이를 사용하시면, y_fitted 값을 따로 계산하실 필요가 없습니다.

      result=mpfitexpr(expr, x, y, sqrt(y), start, weights=weights, yfit=y_fitted)

      이렇게 쓰시면 되지요.

      parinfo 키워드를 사용하시면 가장 상세한 피팅 설정을 하실 수 있습니다. 특정 파라메터에 대해서 고정을 시킨다든지, 상한값, 하한값의 범위를 정해 준다든지(발산하지 못하도록) 지정할 수 있습니다. parinfo 안에 시작값을 정할 수도 있습니다(파라메터 고정 설정까지 하면, 시작값이 곧 최종값으로 나오게 됩니다).

      이렇게 설정합니다.

      parinfo = replicate({value:0.D, fixed:0, limited:[0,0], limits:[0.D,0]}, 4)
      parinfo[3].limited[1]=1
      parinfo[3].limits[1]=10.
      parinfo[*].value=[1, 1, 1, 1]

      첫줄은 parinfo 구조체를 정의하고 이런 것을 4개 만드는 것입니다(파라메터 개수만큼). 위 예제의 구조체 필드가 가지는 값은 디폴트 값이라고 할 수 있습니다. value는 초기값을 의미하고, fixed (0 or 1)은 고정할 것인가를 설정합니다. limited[하한,상한] (0 or 1)은 하한값을 지정할 것인지, 상한값을 지정할 것인지 결정합니다. 그리고 limits 필드를 이용하여 상한값, 하한값을 설정합니다.
      두번째 줄의 의미는 네번째 파라메터의 상한값 (limited[0]은 하한설정, limited[1]은 상한 설정)을 설정하겠다는 것이고,
      세번째 줄의 의미는 네번째 파라메터의 상한값을 10으로 하겠다는 것입니다. 이보다 커질 수가 없게 됩니다.
      네번째 줄은 모든 파라메터의 시작값을 1로 하겠다는 의미가 됩니다.
      사용은, 다음과 같이 할 수 있습니다.

      result=mpfitexpr(expr, x, y, sqrt(y), start, weights=weights, yfit=y_fitted, parinfo=parinfo)

    • #1149 Reply
      chobo
      회원

      매번 감사합니다.
      글을 올리시는 시간이 범상치 않습니다.
      얼른 IDL 고수가 되어서 조금이라도 덜 귀찮게 해드려야되는데요….*^^*

    • #2005 Reply
      chobo
      회원

      안녕하세요..

      mpfit으로 gauss fitting을 하고자 합니다.
      코드는 아래와 같이 작성했는데,
      자꾸 에러가 발생합니다.
      syntax error라는데… ㅜㅜ
      오후 내내 삽질하고 있습니다.
      어느 부분이 문제일까요?

      pro mpfit_test
      x = findgen(21, start = 51.)
      y = [0.259140, 0.320440, 0.352395, 0.415606, 0.466519, 0.618763, 0.702596, 0.842151, 0.980824, 1.01404, $
      1.12439, 0.975067, 0.943691, 0.847766, 0.708441, 0.573130, 0.483977, 0.393586, 0.360299, 0.332075, 0.296387]

      help, x
      help, y

      window, xs = 500, ys = 400
      plot, x, y, psym=2, color = 250

      ; — Gasuss fitting —
      expr = ‘P(0) + GAUSS1(X, P(1:3))’
      err = replicate(0., n_elements(y)) ; no error
      weights = replicate(1., n_elements(y)) ; no weighting = same weighting
      start = [0.2, 1., 1., 1.]
      result = MPFITEXPR(expr, x, y, err, start, weights=weights, yfit = y_fitted)
      oplot, x, y_fitted, color = 200

      stop
      end

      IDL> .RESET_SESSION
      IDL> .compile -v ‘D:\SpaceSienceTeam\SPECRED\codes\mpfit_test.pro’
      % Compiled module: MPFIT_TEST.
      IDL> mpfit_test
      % Compiled module: MPFIT_TEST.
      X FLOAT = Array[21]
      Y FLOAT = Array[21]
      % Compiled module: MPFITEXPR.
      % Compiled module: MPFIT.

      _f = P(0) + GAUSS1(X, P(1:3))
      ^
      % Syntax error.
      % MPFITEXPR: ERROR: execution of “P(0) + GAUSS1(X, P(1:3))” failed.
      % MPFITEXPR: check syntax and parameter usage
      % Illegal subscript range: RESULT.
      % Execution halted at: MPFIT_TEST 21 D:\SpaceSienceTeam\SPECRED\codes\mpfit_test.pro
      % $MAIN$
      IDL>

    • #2009 Reply
      Jonghyuk
      회원

      GAUSS1 함수는 IDL 내장함수가 아닙니다. 보시고 계신 문서를 보시면 앞에 GAUSS1 함수를 정의하는 코드가 있을 거예요. 이를 작성하시거 컴파일 해 놓으면 실행에 문제가 없을 것 같습니다.
      이 문서는 어쨌든 사용자가 직접 만든 함수도 얼마든지 fitting 할 수 있다는 것을 보여 주기 위해서 만든 예제입니다. GAUSS 함수 Fitting이 목적이라면 IDL 내장의 GAUSSFIT 함수를 사용하시는 것이 더 간편할 수 있습니다. 물론 이는 MPFIT과 무관한 함수입니다.

4 답변 글타래를 보이고 있습니다
'mpfit'에 답변달기
글쓴이 정보:




취소