fortran實現正餘弦變換
阿新 • • 發佈:2018-11-06
本來是想正餘弦變換都做的,根據前人研究:正弦變換的結果比餘弦更精確。所以本程式碼只做了正弦變化,如果讀者想實現餘弦變換,只需要把對應的正弦變換系數改為餘弦變換系數即可。在末尾本文會給出250點正餘弦變換系數。 程式碼如下: !.. 本程式碼將實現正餘弦變換 !.. 測試函式 !.. ∫exp(-w)sin(wt)dw !.. 正弦變換離散表示式 !.. f(t)=∑{K(w)*wsl} * coef !.. coef = 1/t*sqrt(pi/2) Module testSineTransform implicit none integer, parameter :: nt = 46 real(kind=8), parameter :: pi = acos(-1.d0) real(kind=8) :: t(nt), t1 = 5, t2 = 50 !.. 單位:s real(kind=8) :: wsl(-149:100), w, delta contains subroutine cal_SineTransform implicit none integer :: i, j, fileid real(kind=8) :: s1(nt), s2(nt) !.. s1為數值解, s2為解析解 do i = 1, nt t(i) = t1 + dble(i-1) * ( t2 - t1 ) / dble( nt - 1) s2(i) = func2( t(i) )!.. 解析表示式 end do open( newunit = fileid, file = "sin250.txt" ) do i = -149, 100 read( fileid,* ) wsl(i) end do close( fileid ) !.. 計算正弦變換 delta = log(10.d0) / 20.d0 !.. 理論推導而來,無需更改 s1 = 0.d0 open( newunit = fileid, file = "result.dat" ) do i = 1, nt do j = -149, 100 w = exp( dble(j)*delta ) / t(i) !.. 理論推導而來,無需更改 s1(i) = s1(i) + func1( w ) * wsl(j) end do s1(i) = s1(i) / t(i) * sqrt( pi/2.d0 ) !.. 1/t*sqrt(pi/2)為正弦變換出來的係數 write( fileid, '(*(g0,3x))' ) t(i), s1(i), s2(i), abs( s1(i)-s2(i) )/s2(i)*100 end do close( fileid ) end subroutine cal_SineTransform real(kind=8) function func1(w) !.. 核函式 implicit none real(kind=8) :: w func1 = exp(-1.d0*w) end function func1 real(kind=8) function func2(t) !.. 解析解 implicit none real(kind=8) :: t func2 = t / ( t*t + 1.d0 ) end function func2 End module testSineTransform Program main use testSineTransform implicit none call cal_SineTransform End program main 下面給出250點正弦變換系數 正弦變換系數: 1.156447054648637e-016 1.455880584491685e-016 1.832845064354326e-016 2.307415227397174e-016 2.904863665331116e-016 3.657006686082886e-016 4.603898648210914e-016 5.795965001557580e-016 7.296687626330448e-016 9.185985474711453e-016 1.156447054648634e-015 1.455880584491681e-015 1.832845064354321e-015 2.307415227397166e-015 2.904863665331103e-015 3.657006686082864e-015 4.603898648210901e-015 5.795965001557501e-015 7.296687626330329e-015 9.185985474711278e-015 1.156447054648612e-014 1.455880584491647e-014 1.832845064354265e-014 2.307415227397080e-014 2.904863665330965e-014 3.657006686082646e-014 4.603898648210535e-014 5.795965001556953e-014 7.296687626329462e-014 9.185985474709900e-014 1.156447054648389e-013 1.455880584491301e-013 1.832845064353717e-013 2.307415227396208e-013 2.904863665329599e-013 3.657006686080477e-013 4.603898648207092e-013 5.795965001551463e-013 7.296687626320796e-013 9.185985474696154e-013 1.156447054646209e-012 1.455880584487839e-012 1.832845064348231e-012 2.307415227387515e-012 2.904863665315808e-012 3.657006686058640e-012 4.603898648172482e-012 5.795965001496609e-012 7.296687626233827e-012 9.185985474558366e-012 1.156447054624371e-011 1.455880584453228e-011 1.832845064293375e-011 2.307415227300576e-011 2.904863665178019e-011 3.657006685840241e-011 4.603898647826370e-011 5.795965000948059e-011 7.296687625364435e-011 9.185985473180428e-011 1.156447054405990e-010 1.455880584107116e-010 1.832845063744824e-010 2.307415226431180e-010 2.904863663800132e-010 3.657006683656435e-010 4.603898644365244e-010 5.795964995462562e-010 7.296687616670500e-010 9.185985459401473e-010 1.156447052222165e-009 1.455880580645994e-009 1.832845058259313e-009 2.307415217737235e-009 2.904863650021138e-009 3.657006661818234e-009 4.603898609754006e-009 5.795964940607466e-009 7.296687529731007e-009 9.185985321611760e-009 1.156447030383951e-008 1.455880546034779e-008 1.832845003404190e-008 2.307415130797806e-008 2.904863512231297e-008 3.657006443436261e-008 4.603898263641547e-008 5.795964392056795e-008 7.296686660335758e-008 9.185983943714872e-008 1.156446812001731e-007 1.455880199922772e-007 1.832844454852834e-007 2.307414261403871e-007 2.904862134332549e-007 3.657004259617825e-007 4.603894802516352e-007 5.795958906554214e-007 7.296677966382738e-007 9.185970164759857e-007 1.156444628179741e-006 1.455876738807346e-006 1.832838969341425e-006 2.307405567480775e-006 2.904848355357452e-006 3.656982421491810e-006 4.603860191324607e-006 5.795904051742752e-006 7.296591027122000e-006 9.185832376013115e-006 1.156422790075531e-005 1.455842127958865e-005 1.832784114711636e-005 2.307318629427335e-005 2.904710567601305e-005 3.656764044807223e-005 4.603514087536846e-005 5.795355521527330e-005 7.295721667410295e-005 9.184454558911375e-005 1.156204422278983e-004 1.455496047235377e-004 1.832235621566062e-004 2.306449358701818e-004 2.903332902825535e-004 3.654580713105382e-004 4.600053900083832e-004 5.789871945981126e-004 7.287031465075219e-004 9.170683245122596e-004 1.154022098295481e-003 1.452037965023791e-003 1.826756086858394e-003 2.297767484717661e-003 2.889577754406488e-003 3.632790471922652e-003 4.565537633632412e-003 5.735207473181395e-003 7.200470157233521e-003 9.033651050110031e-003 1.132334380467910e-002 1.417727733464931e-002 1.772499346509631e-002 2.212023128150156e-002 2.754164273887034e-002 3.419148067351661e-002 4.228846080133485e-002 5.205416439533162e-002 6.368324102630288e-002 7.729769668970302e-002 9.286191525747088e-002 1.100669688846222e-001 1.281325720732841e-001 1.455759302997095e-001 1.598450602566018e-001 1.675994998118554e-001 1.619683745893242e-001 1.351613440916820e-001 8.220710085812771e-002 -5.378512256133189e-003 -1.293051042068928e-001 -2.723112991886841e-001 -3.978468461558793e-001 -4.394993669496010e-001 -3.176425577428793e-001 1.665756359292827e-002 4.630044483056220e-001 7.298213861736705e-001 3.957559489684660e-001 -4.902404478827365e-001 -1.020512348109581e+000 5.251407980261734e-002 1.253864451347815e+000 -1.317577223583377e-001 -1.549527153893942e+000 1.822344193602740e+000 -1.058204213789976e+000 2.629725097214586e-001 1.625578365785221e-001 -2.871269632839427e-001 2.688993168200287e-001 -2.085718537012279e-001 1.489333044554680e-001 -1.023510671993418e-001 6.920616476794303e-002 -4.650126535306088e-002 3.114481903045837e-002 -2.081379273802720e-002 1.390327480927853e-002 -9.298886523349949e-003 6.224588805575879e-003 -4.160927227996026e-003 2.775386859464604e-003 -1.851896921274413e-003 1.239790988028041e-003 -8.315445781273338e-004 5.558587783122442e-004 -3.696533331522879e-004 2.460602463017792e-004 -1.651618293595525e-004 1.114024919622099e-004 -7.389805836530149e-005 4.938931786488272e-005 -3.300905019751596e-005 2.206139792679832e-005 -1.474460111811994e-005 9.854464474452981e-006 -6.586171392444854e-006 4.401827590217889e-006 -2.941934696111318e-006 1.966224160030928e-006 -1.314113957933395e-006 8.782800707768828e-007 -5.869931432254631e-007 3.923132969292148e-007 -2.622002057839472e-007 1.752399127209358e-007 -1.171205297822945e-007 7.827679370240814e-008 -5.231581895777865e-008 3.496495939305279e-008 -2.336861793837973e-008 1.561827366110047e-008 -1.043837820431878e-008 6.976426582137068e-009 -4.662652272535390e-009 3.116255286086561e-009 -2.082730266047119e-009 1.391980105248554e-009 -9.303214367194074e-010 6.217746736151215e-010 -4.155593212088843e-010 2.777365447188800e-010 -1.856235303493902e-010 1.240603574666303e-010 -8.291498532420789e-011 5.541572611672626e-011 -3.703676348776440e-011 2.475329560347511e-011 -1.654371455635934e-011 1.105689099773430e-011 -7.389805845555023e-012 4.938931788889779e-012 -3.300905020390523e-012 2.206139792849932e-012 -1.474460111857239e-012 9.854464474573356e-013 -6.586171392476662e-013 4.401827590226633e-013 -2.941934696113489e-013 餘弦變換系數: 250點餘弦變換濾波係數(王華軍, 2004) n=-149~-100 n=-99~-50 n=-49~0 n=1~50 n=51~100 3.25931064E-09 1.03068452E-06 3.25929012E-04 4.47184906E-02 -1.61220740E-04 3.65700668E-09 1.15644705E-06 3.65697770E-04 3.54831893E-02 1.08092999E-04 4.10322899E-09 1.29755494E-06 4.10318805E-04 2.04513710E-02 -7.22432985E-05 4.60389865E-09 1.45588058E-06 4.60384082E-04 -2.05568659E-03 4.82833692E-05 5.16565924E-09 1.63352488E-06 5.16557757E-04 -3.36444840E-02 -3.22698961E-05 5.79596500E-09 1.83284506E-06 5.79584963E-04 -7.63113734E-02 2.15673888E-05 6.50317969E-09 2.05648598E-06 6.50301672E-04 -1.26643432E-01 -1.44144331E-05 7.29668762E-09 2.30741523E-06 7.29645743E-04 -1.86274938E-01 9.63379869E-06 8.18701817E-09 2.58896246E-06 8.18669301E-04 -2.46135314E-01 -6.43869076E-06 9.18598547E-09 2.90486366E-06 9.18552617E-04 -2.90298367E-01 4.30325981E-06 1.03068452E-08 3.25931064E-06 1.03061964E-03 -2.98807588E-01 -2.87605752E-06 1.15644705E-08 3.65700668E-06 1.15635541E-03 -2.44850614E-01 1.92219555E-06 1.29755494E-08 4.10322898E-06 1.29742549E-03 -9.95001792E-02 -1.28468770E-06 1.45588058E-08 4.60389864E-06 1.45569774E-03 1.35603930E-01 8.58613212E-07 1.63352488E-08 5.16565923E-06 1.63326661E-03 4.09054093E-01 -5.73848918E-07 1.83284506E-08 5.79596499E-06 1.83248024E-03 5.77439580E-01 3.83528434E-07 2.05648598E-08 6.50317967E-06 2.05597066E-03 4.57925568E-01 -2.56328896E-07 2.30741523E-08 7.29668760E-06 2.30668732E-03 -6.98345034E-02 1.71315859E-07 2.58896247E-08 8.18701814E-06 2.58793429E-03 -7.05728308E-01 -1.14497912E-07 2.90486366E-08 9.18598542E-06 2.90341135E-03 -7.84315028E-01 7.65239831E-08 3.25931064E-08 1.03068452E-05 3.25725924E-03 2.52561836E-01 -5.11443386E-08 3.65700668E-08 1.15644704E-05 3.65410908E-03 1.11765483E+00 3.41820077E-08 4.10322899E-08 1.29755492E-05 4.09913616E-03 7.96337157E-02 -2.28453369E-08 4.60389865E-08 1.45588057E-05 4.59811762E-03 -1.54929662E+00 1.52685420E-08 5.16565924E-08 1.63352486E-05 5.15749377E-03 9.69765096E-01 -1.02046372E-08 5.79596500E-08 1.83284503E-05 5.78443175E-03 5.22594045E-01 6.82020717E-09 6.50317969E-08 2.05648593E-05 6.48688994E-03 -1.38177462E+00 -4.55824398E-09 7.29668762E-08 2.30741515E-05 7.27368023E-03 1.42319844E+00 3.04647464E-09 8.18701817E-08 2.58896236E-05 8.15452380E-03 -1.10540419E+00 -2.03609279E-09 9.18598547E-08 2.90486352E-05 9.14009379E-03 7.66503550E-01 1.36081024E-09 1.03068452E-07 3.25931043E-05 1.02420355E-02 -5.12775503E-01 -9.09489244E-10 1.15644705E-07 3.65700639E-05 1.14729492E-02 3.41609163E-01 6.07851604E-10 1.29755494E-07 4.10322858E-05 1.28463164E-02 -2.28175614E-01 -4.06253922E-10 1.45588058E-07 4.60389807E-05 1.43763378E-02 1.52559552E-01 2.71517338E-10 1.63352488E-07 5.16565843E-05 1.60776460E-02 -1.01971829E-01 -1.81466961E-10 1.83284506E-07 5.79596384E-05 1.79648247E-02 6.81663603E-02 1.21282340E-10 2.05648598E-07 6.50317806E-05 2.00516669E-02 -4.55882998E-02 -8.10583142E-11 2.30741523E-07 7.29668532E-05 2.23500306E-02 3.04810683E-02 5.41748312E-11 2.58896247E-07 8.18701492E-05 2.48681715E-02 -2.03615795E-02 -3.62074189E-11 2.90486366E-07 9.18598088E-05 2.76082657E-02 1.35977270E-02 2.41990082E-11 3.25931064E-07 1.03068387E-04 3.05629294E-02 -9.08987662E-03 -1.61732600E-11 3.65700668E-07 1.15644614E-04 3.37101508E-02 6.08272853E-03 1.08092999E-11 4.10322899E-07 1.29755364E-04 3.70064047E-02 -4.06761271E-03 -7.22432985E-12 4.60389865E-07 1.45587876E-04 4.03767590E-02 2.71487245E-03 4.82833692E-12 5.16565924E-07 1.63352230E-04 4.37019801E-02 -1.81128418E-03 -3.22698961E-12 5.79596500E-07 1.83284142E-04 4.68002258E-02 1.21137650E-03 2.15673888E-12 6.50317969E-07 2.05648083E-04 4.94047004E-02 -8.12064069E-04 -1.44144331E-12 7.29668762E-07 2.30740795E-04 5.11324187E-02 5.43426304E-04 9.63379869E-13 8.18701817E-07 2.58895218E-04 5.14507096E-02 -3.61946442E-04 -6.43869076E-13 9.18598547E-07 2.90484914E-04 4.96316608E-02 2.40820925E-04 4.30325981E-13