迭代反投影演算法

㈠ 年北京坐標系與西安坐標系的轉換方法

在礦業權實地核查准備工作階段，收集到的地質、測繪等相關資料、圖件和礦業權登記數據中，所涉及的地理數據可能是不同大地坐標系下的坐標數據。從實際情況來看，礦業權拐點坐標大多採用的是1954年北京坐標系，礦區已有的測量控制點和測量資料多數採用的也是1954年北京坐標系。本次礦業權實地核查測量工作採用的是1980西安坐標系，在實地測量和數據整理中涉及1954年北京坐標系與1980西安坐標系的轉換。下面簡要介紹二者之間轉換的理論與方法。

（一）高斯投影正算和反算

將大地坐標換算為平面直角坐標，叫做高斯投影正算，是在同一橢球中進行，不存在誤差。其常用量定義和公式如下：

a為橢球長半軸

b為橢球短半軸

f為橢球扁率

e為第一偏心率

e＇為第二偏心率

全國礦業權實地核查技術方法指南研究

B為緯度，單位為弧度

全國礦業權實地核查技術方法指南研究

M為子午圈曲率半徑

N為卯酉圈曲率半徑

子午線弧長X

設有子午線上兩點p1和p2,p1在赤道上，p2的緯度為B,p1、p2間的子午線弧長X計算公式：

全國礦業權實地核查技術方法指南研究

例如，1980西安坐標系a=6378140,e²=0.006694385,A'＝1.005052506,B'＝0.002531556209,C'＝2.656901555E-06,D'＝3.470075599E-09,E'＝4.916542167E-12,F '＝7.263137253E-15,G'＝1.074009912E-17以B=30°弧度值0.5235987756為例，在Y=0時算得X=3320114.946。

當Y≠0,l≠0時則需要採用下列積分和逐次趨近的方法。

（1）高斯正算公式（利用點的經緯度計算XY坐標）

全國礦業權實地核查技術方法指南研究

（2）高斯反算公式（利用點的XY坐標計算經緯度）

全國礦業權實地核查技術方法指南研究

（3）底點緯度B_f迭代公式

全國礦業權實地核查技術方法指南研究

直到B_i-1-B_i小於某一個指定數值，即可停止迭代。

式中

全國礦業權實地核查技術方法指南研究

國家測繪局經過改進，將7個系數改為5個算出各橢球的值，採用公式如下：

（1）高斯投影正算（B,L→x,y）

全國礦業權實地核查技術方法指南研究

式中：X₀=C₀B-cosB（C₁sinB+C₂sin²B+C₂sin⁵B+C₄sin⁷B）

m₀=lcosB

l=L－中央子午線經度值（弧度）

L,B為該點的經緯度值。

全國礦業權實地核查技術方法指南研究

式中：t=tanB，η²=e'²cos²B,

C,C₀,C₁,C₂,C₃,C₄,e²為橢球常數

（2）高斯投影反算（x,y→B,L）

全國礦業權實地核查技術方法指南研究

式中：t=tanB_f，η²=e'²cos²B_f，

，K₁K₂,K₃,K₄為橢球常數。

各坐標系橢球常數如表4-1。

表4-1 各大地坐標系橢球常數

國家測繪局採用的公式編程更加容易，高斯投影的正算、反算因為是在同一橢球下進行，公式是嚴密的，不存在誤差，電算操作非常方便。現在網上很多軟體有這種功能。度、分、秒輸入使用小數形式，小數點前面是度，小數點後前兩位為分，後兩位為秒，再後面為秒的十進制小數。如25.23451124其值為25°23′45.1124″，正反算已經成了非常簡單的事。高斯正算、反算必須考慮到橢球參數，橢球不同結果是不同的。必須考慮到中央子午線位置。因為各帶中都有重復點，本次實地核查要求使用3度帶，所有Y坐標必須帶有3°帶的帶號，不允許使用獨立坐標系或假定坐標系。

（二）參心坐標與空間直角坐標的關系

空間直角坐標X、Y、Z與大地坐標B、L、H間的關系表示如下：

全國礦業權實地核查技術方法指南研究

大地坐標B、L、H 與空間直角坐標X、Y、Z間的關系表示如下：

全國礦業權實地核查技術方法指南研究

式中

在轉換中對於不知道橢球高的控制點可將控制點的大地高置為0，放在橢球面上計算，三維就變成二維，其效果更好。

（三）坐標系統轉換

1954年北京坐標系與1980西安坐標系的轉換通常有兩種方法：四參數轉換法和七參數轉換法。

1.四參數轉換法

所謂四參數轉換是兩個平移參數，一個旋轉參數，一個尺度比。不考慮什麼橢球，在小范圍內按平面坐標直接平移、旋轉、縮放。最少條件是兩個公共點，多公共點時可以使用最小二乘法，刪除殘差大的點。這在區域面積小的情況下是可以的，一般不宜超過40平方千米。四參數轉換模型如下：

x₂＝Δx+x₁（1+m）cosa-y（1+m）sina

y₂＝Δx+x₁（1+m）sina-y（1+m）cosa

2.七參數轉換法

該方法適用於橢球間的坐標轉換。其實質是原橢球空間直角坐標（X₁,Y₁,Z₁）與新橢球空間直角坐標（X₂,Y₂,Z₂）間的轉換。橢球間的坐標轉換至少需要3個公共點，解算七參數。轉換公式採用的是布爾莎公式，法方程的解算採用高斯消元法。高斯消元法，是線性代數中的一個演算法，可用來為線性方程組求解，求出矩陣的秩，以及求出可逆方陣的逆矩陣。當用於一個矩陣時，高斯消元法會產生出一個「行梯陣式」。高斯消元法可以用在電腦中來解決數千條等式及未知數。迭代法較消元法的殘差大。

橢球間的坐標轉換適用基於橢球的參心（地心）坐標系間的轉換，而不適用於基於平面的獨立坐標系間以及獨立坐標系和參心（地心）坐標系間的轉換。基於橢球的坐標轉換中（七參數），橢球→橢球的轉換實際上是在空間直角坐標系中完成的。完整的變換過程如下（以「平面→平面」為例）：（x₁,y₁,H₁）→（B₁,L₁,H₁）→（X₁,Y₁,Z₁）→（X₂,Y₂,Z₂）→（B₂,L₂,H₂）→（x₂,y₂,H₂）。首先把直角坐標系下的直角坐標，原公共點中的1954年北京坐標轉換成2000國家大地經緯度坐標，再轉換為1954年北京坐標系的參心坐標，公共點的1980西安坐標做同樣轉換。利用兩個橢球的參心（地心）坐標求得轉換參數，利用該參數直接將1954年北京坐標系下的坐標轉換成1980西安坐標系下的坐標。在上述過程中，高程H₁、H₂是大地高（橢球高）。大地高＝正常高＋測區高程異常。如果不需要轉換高程的話，可以將高程和高程異常全部置為0。不可將1954年北京坐標系坐標所帶的正常高直接代入。

七參數的轉換模型如下：

（1）七參數轉換模型

全國礦業權實地核查技術方法指南研究

式中：ΔB，ΔL為同一點位在兩個坐標系下的緯度差、經度差（弧度）；

a，Δf為橢球長半軸差（米）、扁率差（無量綱）；

X，ΔY，ΔZ為平移參數（米）；

ε_x，ε_y，ε_z為旋轉參數（弧度）；

m為尺度參數（無量綱）。

最少3個公共點可以解求出七個參數。

（2）三維七參數轉換模型

全國礦業權實地核查技術方法指南研究

全國礦業權實地核查技術方法指南研究

式中：ΔB，ΔL，ΔH為同一點位在兩個坐標系下的緯度差（弧度）、經度差（弧度）、大地高差（米）；

ρ為一個弧度的秒值，180×3600/π弧度/秒；

a為橢球長半軸差（米）；

f為扁率差（無量綱）；

X，ΔY，ΔZ為平移參數（米）；

εx，ε_y，ε_z為旋轉參數（弧度）；

m為尺度參數（無量綱）。

最少3個公共點可以解求出七個參數。

七參數適用於整個測區的轉換，面積小於2000平方千米的可以一次轉換完成，面積大的可以分區轉換，各分區之間應選公共點，以保證數據的接邊精度。關於殘差，國家規定以1∶2000圖為例，殘差為圖上0.1毫米即實地20厘米，超過3倍中誤差的點刪除。為了保證礦業權礦界拐點轉化的精度，本次礦業權實地核查規定殘差超過實地0.1米一般不宜使用，實際上比國家規定的精度嚴，相當於國家規定的1/6。

（四）利用坐標轉換軟體進行坐標轉換

以上介紹了1954年北京坐標系和1980西安坐標系轉換的理論，在實際轉換時可以採用相關的軟體來完成。目前，市場上有多種坐標轉換軟體可供選擇。在選擇軟體時，應注意部分軟體轉換的精度可能達不到本次礦業權實地核查的要求。下面以經天測繪技術公司開發的測量計算工具包軟體V4.05為例，介紹坐標轉換方法。

該軟體界面如圖4-3。該軟體可以進行高斯正算、高斯反算、坐標換帶、橢球間的轉換，可以批量導入，可以保存數據、保存公共點，包括了坐標轉換所需的相關計算功能。另外，該軟體還能實現2000國家大地坐標系與1954年北京坐標系、1980西安坐標系、WGS-84坐標系以及獨立坐標系的轉換。

圖4-3 經天測繪技術公司開發的測量計算工具包軟體界面

坐標系統變換，可以採用平面坐標轉換中的多公共點相似變換和橢球坐標轉換。小面積可以採用多公共點相似變換。限制在400平方千米左右，不超過1 幅1∶50000圖。它與中央子午線無關、高程需要置為0，計算參數的輸入文件為文本文件，格式為：

點號，原X 坐標，原Y坐標，新X 坐標，新Y坐標

需要轉換的輸入文件格式為：

點號，原X 坐標，原Y坐標

參數計算點數不超過30個，文件可以導入，公共點可以保存，參數也可以保存。轉換坐標可以導入，轉換後的坐標可以保存。需要注意的是，轉換坐標的位數與計算參數的坐標位數應一致。計算參數不使用帶號，轉換後坐標也沒有帶號。圖4-4中的算例X捨去前4位，Y捨去前3位。

圖4-4 多公共點平面相似變換窗口

面積較大的測區應使用7參數轉換。在橢球間坐標轉換開關下，有平面－平面、大地－平面、平面－大地、大地－大地4個子開關。對於采礦權，可使用平面－平面；對於探礦權，使用大地－大地，小數後位數較多，根據需要可將尾部刪去。輸入文件的格式與上述相同，需要輸入中央子午線，Y坐標不加帶號，在不知道1954年北京坐標、1980西安坐標的橢球高的情況下，可在高程欄輸入0，測區高程異常輸入0，探礦權是大地坐標格式，小數點前3位為°，後2位為′，3、4位為″，後面為十進制的秒的小數，如108°33′15″8563，輸入108.33158563，由於控制點坐標是X、Y格式，可用高斯投影反算將控制點變為大地坐標格式。或是使用高斯坐標正算把探礦權登記坐標轉換為直角坐標，計算完成後再使用高斯坐標反算將1980西安坐標轉換為2000國家坐標。圖4-5表示一個縣的采礦權轉換過程，Y坐標略去了前3位數。

圖4-5 橢球間平面坐標轉換窗口

需要注意的是，該軟體沒有採用軟體狗加密，但需要注冊才能用，採用機器碼注冊，一個軟體只能裝一台計算機專用。

㈡ lbp是什麼意思

lbp是Local Binary Patterns的縮寫，中文名字叫做「線性反投影演算法」，指局部二值模式，最初功能為輔助圖像局部對比度，並不是一個完整的特徵描述子。

lbp線性反投影演算法又稱累加法，是最早使用的一種簡單成像演算法。它將通過某點的所有投影射線進行累加，再反向估算出該點的密度值，從成像觀點分析，它是不完全的雷登逆變換。

lbp演算法計算步驟：

1、基於均質靈敏度信息，利用線性反投影演算法獲得初始圖像。

2、利用已獲得的介電常數分布，求解正問題，得到一組模擬電容值，將該值與測量電容值進行比較，若誤差已達到滿意值，演算法結束，否則進行下步。

3、修正靈敏度信息。

4、根據上步已經修正後的靈敏度信息，利用測量電壓重新進行線性反投影；返回第二步，並進行循環迭代，直到獲得滿意的結果為止，迭代結束。

㈢ 迭代法的演算法

迭代是數值分析中通過從一個初始估計出發尋找一系列近似解來解決問題（一般是解方程或者方程組）的過程，為實現這一過程所使用的方法統稱為迭代法（Iterative Method）。
一般可以做如下定義：對於給定的線性方程組（這里的x、B、f同為矩陣，任意線性方程組都可以變換成此形式），用公式 （代表迭代k次得到的x，初始時k=0)逐步帶入求近似解的方法稱為迭代法（或稱一階定常迭代法）。如果存在，記為x*，稱此迭代法收斂。顯然x*就是此方程組的解，否則稱為迭代法發散。
跟迭代法相對應的是直接法（或者稱為一次解法），即一次性的快速解決問題，例如通過開方解決方程x +3= 4。一般如果可能，直接解法總是優先考慮的。但當遇到復雜問題時，特別是在未知量很多，方程為非線性時，我們無法找到直接解法（例如五次以及更高次的代數方程沒有解析解，參見阿貝耳定理），這時候或許可以通過迭代法尋求方程（組）的近似解。
最常見的迭代法是牛頓法。其他還包括最速下降法、共軛迭代法、變尺度迭代法、最小二乘法、線性規劃、非線性規劃、單純型法、懲罰函數法、斜率投影法、遺傳演算法、模擬退火等等。
利用迭代演算法解決問題，需要做好以下三個方面的工作： 例 1 ：一個飼養場引進一隻剛出生的新品種兔子，這種兔子從出生的下一個月開始，每月新生一隻兔子，新生的兔子也如此繁殖。如果所有的兔子都不死去，問到第 12 個月時，該飼養場共有兔子多少只？
分析：這是一個典型的遞推問題。我們不妨假設第 1 個月時兔子的只數為 u 1 ，第 2 個月時兔子的只數為 u 2 ，第 3 個月時兔子的只數為 u 3 ，……根據題意，「這種兔子從出生的下一個月開始，每月新生一隻兔子」，則有
u 1 = 1 ， u 2 = u 1 + u 1 × 1 = 2 ， u 3 = u 2 + u 2 × 1 = 4 ，……
根據這個規律，可以歸納出下面的遞推公式：
u n = u(n － 1）× 2 (n ≥ 2）
對應 u n 和 u(n － 1），定義兩個迭代變數 y 和 x ，可將上面的遞推公式轉換成如下迭代關系：
y=x*2
x=y
讓計算機對這個迭代關系重復執行 11 次，就可以算出第 12 個月時的兔子數。參考程序如下：
cls
x=1
for i=2 to 12
y=x*2
x=y
next i
print y
end
例 2 ：阿米巴用簡單分裂的方式繁殖，它每分裂一次要用 3 分鍾。將若干個阿米巴放在一個盛滿營養參液的容器內， 45 分鍾後容器內充滿了阿米巴。已知容器最多可以裝阿米巴 220,220個。試問，開始的時候往容器內放了多少個阿米巴？請編程序算出。
分析：根據題意，阿米巴每 3 分鍾分裂一次，那麼從開始的時候將阿米巴放入容器裡面，到 45 分鍾後充滿容器，需要分裂 45/3=15 次。而「容器最多可以裝阿米巴2^ 20 個」，即阿米巴分裂 15 次以後得到的個數是 2^20。題目要求我們計算分裂之前的阿米巴數，不妨使用倒推的方法，從第 15 次分裂之後的 2^20 個，倒推出第 15 次分裂之前（即第 14 次分裂之後）的個數，再進一步倒推出第 13 次分裂之後、第 12 次分裂之後、……第 1 次分裂之前的個數。
設第 1 次分裂之前的個數為 x 0 、第 1 次分裂之後的個數為 x 1 、第 2 次分裂之後的個數為 x 2 、……第 15 次分裂之後的個數為 x 15 ，則有
x 14 =x 15 /2 、 x 13 =x 14 /2 、…… x n-1 =x n /2 (n ≥ 1）
因為第 15 次分裂之後的個數 x 15 是已知的，如果定義迭代變數為 x ，則可以將上面的倒推公式轉換成如下的迭代公式：
x=x/2 （x 的初值為第 15 次分裂之後的個數 2^20）
讓這個迭代公式重復執行 15 次，就可以倒推出第 1 次分裂之前的阿米巴個數。因為所需的迭代次數是個確定的值，我們可以使用一個固定次數的循環來實現對迭代過程的控制。參考程序如下：
cls
x=2^20
for i=1 to 15
x=x/2
next i
print x
end
ps:java中冪的演算法是Math.pow（2,20）；返回double，稍微注意一下
例 3 ：驗證谷角猜想。日本數學家谷角靜夫在研究自然數時發現了一個奇怪現象：對於任意一個自然數 n ，若 n 為偶數，則將其除以 2 ；若 n 為奇數，則將其乘以 3 ，然後再加 1。如此經過有限次運算後，總可以得到自然數 1。人們把谷角靜夫的這一發現叫做「谷角猜想」。
要求：編寫一個程序，由鍵盤輸入一個自然數 n ，把 n 經過有限次運算後，最終變成自然數 1 的全過程列印出來。
分析：定義迭代變數為 n ，按照谷角猜想的內容，可以得到兩種情況下的迭代關系式：當 n 為偶數時， n=n/2 ；當 n 為奇數時， n=n*3+1。用 QBASIC 語言把它描述出來就是：
if n 為偶數 then
n=n/2
else
n=n*3+1
end if
這就是需要計算機重復執行的迭代過程。這個迭代過程需要重復執行多少次，才能使迭代變數 n 最終變成自然數 1 ，這是我們無法計算出來的。因此，還需進一步確定用來結束迭代過程的條件。仔細分析題目要求，不難看出，對任意給定的一個自然數 n ，只要經過有限次運算後，能夠得到自然數 1 ，就已經完成了驗證工作。因此，用來結束迭代過程的條件可以定義為：n=1。參考程序如下：
cls
input Please input n=;n
do until n=1
if n mod 2=0 then
rem 如果 n 為偶數，則調用迭代公式 n=n/2
n=n/2
print —;n;
else
n=n*3+1
print —;n;
end if
loop
end
迭代法開平方：
#include<stdio.h>
#include<math.h>
void main()
{
double a,x0,x1;
printf(Input a: );
scanf(%lf,&a);//因為a是double型數據，所以要用%lf，而不是%f
if(a<0)
printf(Error! );
else
{
x0=a/2;
x1=(x0+a/x0)/2;
do
{
x0=x1;
x1=(x0+a/x0)/2;
}while(fabs(x0-x1）>=1e-6）；
}
printf(Result: );
printf(sqrt(%g)=%g ,a,x1）；
}
求平方根的迭代公式：x1=1/2*(x0+a/x0）。
演算法：1.先自定一個初值x0，作為a的平方根值，在我們的程序中取a/2作為a的初值；利用迭代公式求出一個x1。此值與真正的a的平方根值相比，誤差很大。
⒉把新求得的x1代入x0中，准備用此新的x0再去求出一個新的x1.
⒊利用迭代公式再求出一個新的x1的值，也就是用新的x0又求出一個新的平方根值x1，此值將更趨近於真正的平方根值。
⒋比較前後兩次求得的平方根值x0和x1，如果它們的差值小於我們指定的值，即達到我們要求的精度，則認為x1就是a的平方根值，去執行步驟5；否則執行步驟2，即循環進行迭代。
迭代法是用於求方程或方程組近似根的一種常用的演算法設計方法。設方程為f(x)=0，用某種數學方法導出等價的形式x=g(x），然後按以下步驟執行：
⑴ 選一個方程的近似根，賦給變數x0；
⑵ 將x0的值保存於變數x1，然後計算g(x1），並將結果存於變數x0；
⑶ 當x0與x1的差的絕對值還小於指定的精度要求時，重復步驟⑵的計算。
若方程有根，並且用上述方法計算出來的近似根序列收斂，則按上述方法求得的x0就認為是方程的根。上述演算法用C程序的形式表示為：
【演算法】迭代法求方程的根
{ x0=初始近似根；
do {
x1=x0；
x0=g(x1）； /*按特定的方程計算新的近似根*/
} while (fabs(x0-x1）>Epsilon）；
printf（「方程的近似根是%f 」，x0）；
}
迭代演算法也常用於求方程組的根，令
X=（x0，x1，…，xn-1）
設方程組為：
xi=gi(X) (I=0，1，…，n-1）
則求方程組根的迭代演算法可描述如下：
【演算法】迭代法求方程組的根
{ for (i=0;i
x=初始近似根；
do {
for (i=0;i
y=x;
for (i=0;i
x=gi(X);
for (delta=0.0,i=0;i
if (fabs(y-x)>delta) delta=fabs(y-x）；
} while (delta>Epsilon）；
for (i=0;i
printf（「變數x[%d]的近似根是 %f」，I，x）；
printf（「 」）；
}
具體使用迭代法求根時應注意以下兩種可能發生的情況：
⑴ 如果方程無解，演算法求出的近似根序列就不會收斂，迭代過程會變成死循環，因此在使用迭代演算法前應先考察方程是否有解，並在程序中對迭代的次數給予限制；
⑵ 方程雖然有解，但迭代公式選擇不當，或迭代的初始近似根選擇不合理，也會導致迭代失敗。
遞歸
遞歸是設計和描述演算法的一種有力的工具，由於它在復雜演算法的描述中被經常採用，為此在進一步介紹其他演算法設計方法之前先討論它。
能採用遞歸描述的演算法通常有這樣的特徵：為求解規模為N的問題，設法將它分解成規模較小的問題，然後從這些小問題的解方便地構造出大問題的解，並且這些規模較小的問題也能採用同樣的分解和綜合方法，分解成規模更小的問題，並從這些更小問題的解構造出規模較大問題的解。特別地，當規模N=1時，能直接得解。
【問題】 編寫計算斐波那契（Fibonacci）數列的第n項函數fib（n）。
斐波那契數列為：0、1、1、2、3、……，即：
fib(0)=0;
fib⑴=1;
fib(n)=fib(n-1）+fib(n-2） （當n>1時）。
寫成遞歸函數有：
int fib(int n)
{ if (n==0) return 0;
if (n==1） return 1;
if (n>1） return fib(n-1）+fib(n-2）；
}
遞歸演算法的執行過程分遞推和回歸兩個階段。在遞推階段，把較復雜的問題（規模為n）的求解推到比原問題簡單一些的問題（規模小於n）的求解。例如上例中，求解fib(n），把它推到求解fib(n-1）和fib(n-2）。也就是說，為計算fib(n），必須先計算fib(n-1）和fib(n- 2），而計算fib(n-1）和fib(n-2），又必須先計算fib(n-3）和fib(n-4）。依次類推，直至計算fib⑴和fib(0），分別能立即得到結果1和0。在遞推階段，必須要有終止遞歸的情況。例如在函數fib中，當n為1和0的情況。
在回歸階段，當獲得最簡單情況的解後，逐級返回，依次得到稍復雜問題的解，例如得到fib⑴和fib(0）後，返回得到fib⑵的結果，……，在得到了fib(n-1）和fib(n-2）的結果後，返回得到fib(n）的結果。
在編寫遞歸函數時要注意，函數中的局部變數和參數知識局限於當前調用層，當遞推進入「簡單問題」層時，原來層次上的參數和局部變數便被隱蔽起來。在一系列「簡單問題」層，它們各有自己的參數和局部變數。
由於遞歸引起一系列的函數調用，並且可能會有一系列的重復計算，遞歸演算法的執行效率相對較低。當某個遞歸演算法能較方便地轉換成遞推演算法時，通常按遞推演算法編寫程序。例如上例計算斐波那契數列的第n項的函數fib(n）應採用遞推演算法，即從斐波那契數列的前兩項出發，逐次由前兩項計算出下一項，直至計算出要求的第n項。
【問題】 組合問題
問題描述：找出從自然數1、2、……、n中任取r個數的所有組合。例如n=5，r=3的所有組合為：⑴5、4、3 ⑵5、4、2 ⑶5、4、1
⑷5、3、2 ⑸5、3、1 ⑹5、2、1
⑺4、3、2 ⑻4、3、1 ⑼4、2、1
⑽3、2、1
分析所列的10個組合，可以採用這樣的遞歸思想來考慮求組合函數的演算法。設函數為void comb(int m,int k）為找出從自然數1、2、……、m中任取k個數的所有組合。當組合的第一個數字選定時，其後的數字是從餘下的m-1個數中取k-1數的組合。這就將求m 個數中取k個數的組合問題轉化成求m-1個數中取k-1個數的組合問題。設函數引入工作數組a[ ]存放求出的組合的數字，約定函數將確定的k個數字組合的第一個數字放在a[k]中，當一個組合求出後，才將a[ ]中的一個組合輸出。第一個數可以是m、m-1、……、k，函數將確定組合的第一個數字放入數組後，有兩種可能的選擇，因還未去頂組合的其餘元素，繼續遞歸去確定；或因已確定了組合的全部元素，輸出這個組合。細節見以下程序中的函數comb。
【程序】
# include
# define MAXN 100
int a[MAXN];
void comb(int m,int k)
{ int i,j;
for (i=m;i>=k;i--)
{ a[k]=i;
if (k>1）
comb(i-1,k-1）；
else
{ for (j=a[0];j>0;j--)
printf（「%4d」，a[j]);
printf（「 」）；
}
}
}
void main()
{ a[0]=3;
comb（5,3）；
}
【問題】 背包問題
問題描述：有不同價值、不同重量的物品n件，求從這n件物品中選取一部分物品的選擇方案，使選中物品的總重量不超過指定的限制重量，但選中物品的價值之和最大。
設n 件物品的重量分別為w0、w1、…、wn-1，物品的價值分別為v0、v1、…、vn-1。採用遞歸尋找物品的選擇方案。設前面已有了多種選擇的方案，並保留了其中總價值最大的方案於數組option[ ]，該方案的總價值存於變數maxv。當前正在考察新方案，其物品選擇情況保存於數組cop[ ]。假定當前方案已考慮了前i-1件物品，現在要考慮第i件物品；當前方案已包含的物品的重量之和為tw；至此，若其餘物品都選擇是可能的話，本方案能達到的總價值的期望值為tv。演算法引入tv是當一旦當前方案的總價值的期望值也小於前面方案的總價值maxv時，繼續考察當前方案變成無意義的工作，應終止當前方案，立即去考察下一個方案。因為當方案的總價值不比maxv大時，該方案不會被再考察，這同時保證函數後找到的方案一定會比前面的方案更好。
對於第i件物品的選擇考慮有兩種可能：
⑴ 考慮物品i被選擇，這種可能性僅當包含它不會超過方案總重量限制時才是可行的。選中後，繼續遞歸去考慮其餘物品的選擇。
⑵ 考慮物品i不被選擇，這種可能性僅當不包含物品i也有可能會找到價值更大的方案的情況。
按以上思想寫出遞歸演算法如下：
try（物品i，當前選擇已達到的重量和，本方案可能達到的總價值tv)
{ /*考慮物品i包含在當前方案中的可能性*/
if（包含物品i是可以接受的）
{ 將物品i包含在當前方案中；
if (i
try(i+1,tw+物品i的重量，tv);
else
/*又一個完整方案，因為它比前面的方案好，以它作為最佳方案*/
以當前方案作為臨時最佳方案保存；
恢復物品i不包含狀態；
}
/*考慮物品i不包含在當前方案中的可能性*/
if （不包含物品i僅是可男考慮的）
if (i
try(i+1,tw,tv-物品i的價值）；
else
/*又一個完整方案，因它比前面的方案好，以它作為最佳方案*/
以當前方案作為臨時最佳方案保存；
}
為了理解上述演算法，特舉以下實例。設有4件物品，它們的重量和價值見表：
物品 0 1 2 3
重量 5 3 2 1
價值 4 4 3 1
並設限制重量為7。則按以上演算法，下圖表示找解過程。由圖知，一旦找到一個解，演算法就進一步找更好的佳。如能判定某個查找分支不會找到更好的解，演算法不會在該分支繼續查找，而是立即終止該分支，並去考察下一個分支。
按上述演算法編寫函數和程序如下：
【程序】
# include
# define N 100
double limitW,totV,maxV;
int option[N],cop[N];
struct { double weight;
double value;
}a[N];
int n;
void find(int i,double tw,double tv)
{ int k;
/*考慮物品i包含在當前方案中的可能性*/
if (tw+a.weight<=limitW)
{ cop=1;
if (i
else
{ for (k=0;k
option[k]=cop[k];
maxv=tv;
}
cop=0;
}
/*考慮物品i不包含在當前方案中的可能性*/
if (tv-a.value>maxV)
if (i
else
{ for (k=0;k
option[k]=cop[k];
maxv=tv-a.value;
}
}
void main()
{ int k;
double w,v;
printf（「輸入物品種數 」）；
scanf（（「%d」，&n);
printf（「輸入各物品的重量和價值 」）；
for (totv=0.0,k=0;k
{ scanf（「%1f%1f」，&w,&v);
a[k].weight=w;
a[k].value=v;
totV+=V;
}
printf（「輸入限制重量 」）；
scanf（「%1f」，&limitV);
maxv=0.0;
for (k=0;k find(0,0.0,totV);
for (k=0;k
if (option[k]) printf（「%4d」，k+1）；
printf（「 總價值為%.2f 」，maxv);
}
作為對比，下面以同樣的解題思想，考慮非遞歸的程序解。為了提高找解速度，程序不是簡單地逐一生成所有候選解，而是從每個物品對候選解的影響來形成值得進一步考慮的候選解，一個候選解是通過依次考察每個物品形成的。對物品i的考察有這樣幾種情況：當該物品被包含在候選解中依舊滿足解的總重量的限制，該物品被包含在候選解中是應該繼續考慮的；反之，該物品不應該包括在當前正在形成的候選解中。同樣地，僅當物品不被包括在候選解中，還是有可能找到比目前臨時最佳解更好的候選解時，才去考慮該物品不被包括在候選解中；反之，該物品不包括在當前候選解中的方案也不應繼續考慮。對於任一值得繼續考慮的方案，程序就去進一步考慮下一個物品。
【程序】
# include
# define N 100
double limitW;
int cop[N];
struct ele { double weight;
double value;
} a[N];
int k,n;
struct { int ;
double tw;
double tv;
}twv[N];
void next(int i,double tw,double tv)
{ twv.=1;
twv tw=tw;
twv tv=tv;
}
double find(struct ele *a,int n)
{ int i,k,f;
double maxv,tw,tv,totv;
maxv=0;
for (totv=0.0,k=0;k
totv+=a[k].value;
next(0,0.0,totv);
i=0;
While (i>=0)
{ f=twv.;
tw=twv tw;
tv=twv tv;
switch(f)
{ case 1: twv.++;
if (tw+a.weight<=limitW)
if (i
{ next(i+1,tw+a.weight,tv);
i++;
}
else
{ maxv=tv;
for (k=0;k
cop[k]=twv[k].!=0;
}
break;
case 0: i--;
break;
default: twv.=0;
if (tv-a.value>maxv)
if (i
{ next(i+1,tw,tv-a.value);
i++;
}
else
{ maxv=tv-a.value;
for (k=0;k
cop[k]=twv[k].!=0;
}
break;
}
}
return maxv;
}
void main()
{ double maxv;
printf（「輸入物品種數 」）；
scanf（（「%d」，&n);
printf（「輸入限制重量 」）；
scanf（「%1f」，&limitW);
printf（「輸入各物品的重量和價值 」）；
for (k=0;k
scanf（「%1f%1f」，&a[k].weight,&a[k].value);
maxv=find(a,n);
printf（「 選中的物品為 」）；
for (k=0;k
if (option[k]) printf（「%4d」，k+1）；
printf（「 總價值為%.2f 」，maxv);
}

㈣ lms演算法是什麼

LMS（Least mean square）演算法，即最小均方誤差演算法。

lms演算法的特點

根據小均方誤差准則以及均方誤差曲面，自然的我們會想到沿每一時刻均方誤差 的陡下降在權向量面上的投影方向更新，也就是通過目標函數的反梯度向量來反 復迭代更新。由於均方誤差性能曲面只有一個唯一的極小值，只要收斂步長選擇恰當， 不管初始權向量在哪，後都可以收斂到誤差曲面的小點，或者是在它的一個鄰域內。

㈤ 誰能幫忙說下CT原理和反投影重建演算法是神馬書上內容太詭異了，希望用自己的經驗總結簡單一點說明。

把採集到的圖象用仿射變換配准，
為了加快運行速度可以先進行展開。
配准這一步可以在空間域，
也可在頻率域進行
然後按配准結果將這些圖象插合成一幅圖象，
再用最小二乘法求解線性方程組即可。

注意，
最好使用超鬆弛迭代法求解，
但是遇到0的時候結果可能有較大出入，
解決辦法中的一種是圖象矩陣所有元素全部加上1，
計算完成後再全部減去1，
然後再512級灰度量化

這是最簡單的重構方法，
沒有考慮圖象的模糊效應。

此外，如果有矩陣維度問題，
有兩種解決辦法，
一是將插合圖象變成正方形圖象，
一是將各插合行，列按權值累加，
反向映射，
後一種速度快些，
也不必直接求解方程，
但是不具有通用性。

㈥ camshift演算法可以和什麼結合

camshift利用目標的顏色直方圖模型將圖像轉換為顏色概率分布圖，初始化一個搜索窗的大小和位置，並根據上一幀得到的結果自適應調整搜索窗口的位置和大小，從而定位出當前圖像中目標的中心位置。

分為三個部分：
1--色彩投影圖（反向投影）：
(1).RGB顏色空間對光照亮度變化較為敏感，為了減少此變化對跟蹤效果的影響，首先將圖像從RGB空間轉換到HSV空間。(2).然後對其中的H分量作直方圖，在直方圖中代表了不同H分量值出現的概率或者像素個數，就是說可以查找出H分量大小為h的概率或者像素個數，即得到了顏色概率查找表。(3).將圖像中每個像素的值用其顏色出現的概率對替換，就得到了顏色概率分布圖。這個過程就叫反向投影，顏色概率分布圖是一個灰度圖像。

2--meanshift
meanshift演算法是一種密度函數梯度估計的非參數方法，通過迭代尋優找到概率分布的極值來定位目標。
演算法過程為：
(1).在顏色概率分布圖中選取搜索窗W
(2).計算零階距：

計算一階距：

計算搜索窗的質心：

(3).調整搜索窗大小
寬度為；長度為1.2s；
(4).移動搜索窗的中心到質心，如果移動距離大於預設的固定閾值，則重復2)3)4)，直到搜索窗的中心與質心間的移動距離小於預設的固定閾值，或者循環運算的次數達到某一最大值，停止計算。關於meanshift的收斂性證明可以google相關文獻。

3--camshift
將meanshift演算法擴展到連續圖像序列，就是camshift演算法。它將的所有幀做meanshift運算，並將上一幀的結果，即搜索窗的大小和中心，作為下一幀meanshift演算法搜索窗的初始值。如此迭代下去，就可以實現對目標的跟蹤。
演算法過程為：
(1).初始化搜索窗
(2).計算搜索窗的顏色概率分布（反向投影）
(3).運行meanshift演算法，獲得搜索窗新的大小和位置。
(4).在下一幀圖像中用(3)中的值重新初始化搜索窗的大小和位置，再跳轉到(2)繼續進行。

㈦ 地震波層析技術

地學層析成像技術（Geotomography，簡稱GT）是地球物理學與層析成像技術相結合的產物，它是利用地震波或電磁波在地下介質中的傳播特性，通過反演演算法來重建地下介質的分布圖像。目前，GT技術已廣泛用於地球物理學的多個研究領域，從地球科學研究中地球內部精細結構的探測，石油勘探開發中岩性圈閉的尋找，金屬礦勘探中深部盲礦的勘查，到工程、災害和環境地質中查明地下構造、裂隙等的分布規律，它已發展成為研究地下精細結構的一種有效手段。

6.3.2.1 基本原理

層析成像技術是根據對物體外部獲取的某種物理量的測定值（或稱投影）進行處理以重建物體內部圖像的一種技術。所謂層析成像，是對物體進行逐層剖析成像，相當於把物體切成片。用波去穿透物體，讓波帶出關於物體內部的信息，通過對這些信息的處理來重建物體的內部圖像。

從物體內部圖像重建的角度看，一物體切片的圖像是兩個空間變數（x，y）的函數，稱為圖像函數，記為f（x，y）。用不同方向的入射波穿過物體，觀測到的波場信息至少是入射波方向θ和觀測點位置ρ兩個變數的函數，稱為投影函數，可記為u（θ，ρ）。1917年奧地利數學家Radon證明，已知所有入射角的投影函數u（θ，ρ）可以恢復惟一的圖像函數f（x，y），人們稱之為Radon變換，它成為層析成像的理論基礎。

從物理角度來看，各種波動（地震波、電磁波等）在時間和空間上都是連續的。當它們在介質中傳播時，由於介質物性上的差異（如密度、速度、介電常數、電導率等）使得它們的傳播速度的大小及方向發生變化，能量吸收也因介質而異。因此，當波穿過某一物體時，必定會把物體內部的物性參數的信息攜帶到物體外部來，只要在物體外部測得波穿過物體後的有關參數（即投影函數），通過適當的反演方法便可重建物體內部圖像。

因此，從本質而言，層析成像屬於反演問題，因為它是通過參數的觀測結果來求解參數的空間分布。一般來說，層析成像方法適用於能以數據形式獲得某種參數（例如地震波傳播時間等）的條件，但入射波必須穿透物體並且參數的變化應該能測得出來。

圖6.3.4 跨孔地震層析成像觀測系統

6.3.2.2 數據採集

在地球物理層析成像中，一般都是將待成像的區域劃分為許多單元組成的網格，單元被稱為「像元」。假定每個像元為常數，它代表了圖像函數f（x，y）。像元的大小取決於諸多因素，如震源和接收器的間距等。由發射源發射出的具有一定方向性的射束被稱為「射線」，射線穿過成像區域到達接收點，接收器記錄到變化了的信號。為了獲得多方向的射線覆蓋，數據採集常在物體的一側布設線性發射源排列，而在另一側沿平行發射源排列布設接收器排列。這實際上是跨孔排列（圖6.3.4所示），它是數據採集方法中最典型的一種。

6.3.2.3 成像方法

地震層析成像方法大致可分為兩種類型：一種是基於射線理論的圖像重建技術，在數學上也就是由一個函數的線積分反求這個函數的問題，它可以化為相似的代數方程組。走時反演成像和振幅反演成像均屬這種方法。當射線為直線時，這類方法比較成熟。在實際應用中，由於客觀條件的限制，經常會遇到數據不全的問題。另外，當必須考慮射線彎曲的情況時，理論上還有不少困難。另一種是基於波動方程反演的散射（或衍射）層析成像方法，在數學上它屬於解偏微分方程反問題，波場成像屬此類方法。

走時反演成像方法是目前使用較為成熟的方法，主要包括反投影法（BPT）、代數重建法（ART）、聯合迭代重建法（SIRT）、共軛梯度最小二乘法（CGLS）和正交分解最小二乘法（LSQR）等。在工程、環境調查中，通常在現場需要對數據進行處理和解釋，故一般選用較為簡單的成像方法。

6.3.2.4 在工程和環境調查中的應用

跨孔地震層析成像在岩溶勘察中的應用，是基於完整灰岩與岩溶（包括充填物）、溶蝕裂隙及上覆土層之間存在明顯的彈性縱波波速差異。一般來說，完整石灰岩的彈性縱波速度大於4500 m/s，而溶蝕裂隙發育灰岩的彈性縱波波速則在2800～4500 m/s之間，岩溶充填物及上覆土層的彈性縱波波速小於2800 m/s。因此，在岩溶發育地區開展地震層析成像勘察具有良好的地球物理前提。

例如，在位於廣州市雅崗與佛山南海市和順之間的廣合大橋，設計有4個主橋墩，32根沖孔灌注樁。該橋基岩為石炭紀灰岩，基岩面埋深約為19～36 m，上部覆蓋層為第四紀沖淤積、殘積淤泥、砂及粘土，基岩面起伏變化大，岩溶裂隙非常發育。在主橋墩施工中採用跨孔地震層析成像方法對4個主橋墩進行勘察。勘察要求查明主橋墩樁位及樁周基岩面埋深，岩溶裂隙發育及分布情況。勘察深度范圍在岩面以上，-5～-50 m高程。要求分辨線性尺度1.0 m以上的溶洞及溶蝕裂隙發育帶。

勘察中共布置了26個鑽孔，可組合成46對跨孔地震層析成像剖面。跨孔距最大為16.02 m，最小為5.95 m。每對剖面測試范圍為基岩面以上1/2跨孔距（且不小於5.0 m）至孔底，激發接收點距為1.0 m。震源主頻高於500 Hz，感測器頻響范圍為5～4000 Hz，接收儀器頻響為10～4000 Hz，采樣間隔為31.25 μs。

對46對跨孔地震層析成像剖面走時數據進行處理、反演，得到了46對跨孔地震層析成像縱波速度色譜圖。其中L26剖面跨孔地震層析成像（縱波）速度色譜圖（如圖6.3.5所示）。從圖中可以看出，速度色譜圖從上至下大致可以分為3個速度帶。根據縱波速度與土層、溶洞、溶蝕裂隙發育及完整（或基本完整）灰岩的相關關系並結合兩側鑽孔資料，可將L26剖面跨孔地震CT速度色譜圖從上至下解釋為3個地質單元：覆蓋層、岩溶及溶蝕裂隙發育帶（局部含完整岩塊）、完整（或基本完整）灰岩（局部偶有小溶洞或溶蝕裂隙）。同理對其他45對CT剖面分別作了地質解釋，推薦了各根樁的樁端持力層高程，同時指出了各根樁施工中可能碰到的不良地質體情況。

圖6.3.5 大橋L26測線地震CT速度色譜圖及地質解釋剖面

㈧ 求一份詳盡的基於迭代反投影演算法的超解析度重建的源代碼

這個太高深了