tag:blogger.com,1999:blog-13668189902893893282024-02-08T09:34:32.410-08:00Dr.B. Witchayangkoon: GPS & GIS for Civil & Transportation EngineeringGPS Precise Point Positioningดร.บุญทรัพย์ วิชญางกูรhttp://www.blogger.com/profile/09210204627753495612noreply@blogger.comBlogger7125tag:blogger.com,1999:blog-1366818990289389328.post-58665468401649890252009-05-06T04:13:00.000-07:002012-12-07T02:50:38.584-08:00Dr.Boonsap Witchayangkoon @ Faculty of Engineering @ Thammasat University<span style="font-size: 180%;">รองศาสตราจารย์ ดร.บุญทรัพย์ วิชญางกูร<br />Boonsap Witchayangkoon<br />(PhD, Associate Professor)</span><br />
<br />
ภาควิชาวิศวกรรมโยธา<br />
<span style="font-size: 130%;"><a href="http://www.engr.tu.ac.th/" target="_blank">คณะวิศวกรรมศาสตร์ มหาวิทยาลัยธรรมศาสตร์<br />Engineering Building, Thammasat University</a><br />Pathumtani, THAILAND<br />Tel : (66)02-564-3001-9 ext. 3101, 3039<br />wboon @ engr.tu.ac.th</span><br />
<br />
<span style="font-size: 180%;"><span style="font-weight: bold;">Education</span></span><br />
<ul style="font-family: georgia;">
<li><span style="font-size: 130%;">PhD <a href="http://spatial.umaine.edu/" target="_blank">(Spatial Information Science and Engineering) 2000, University of Maine, USA</a></span></li>
<li><span style="font-size: 130%;">B.Eng. (Honors) (Civil Engineering) 1992, <a href="http://www.kmutt.ac.th/" target="_blank">King Mongkut's Institute (now University) of Technology Thonburi, Thailand</a></span></li>
<li><span style="font-size: 130%;">LL.B. 2007, <a _blank="" href="http://www.stou.ac.th%20target=/">Sukhothai Thammathirat Open University.</a></span></li>
</ul>
<br />
<br />
<span style="font-size: 180%; font-weight: bold;">Committee Member </span><span style="font-size: 130%;">(partial list):</span><br />
<ul style="font-family: georgia;"><span style="font-size: 130%;"><br />
<li>Committee Member of the Thammasat University Board of Trusstee 2005-2007 กรรมการสภามหาวิทยาลัยธรรมศาสตร์ 2548-50</li>
<br />
<li>Committee Member of special engineering task of the Thailand's Office of The Civil Service Commission for promotion consideration to escalate position of the officer: 2009. คณะกรรมการเฉพาะกิจ สาขาวิศวกรรมศาสตร์ 1 สำนักงาน ก.พ. (พิจารณาผลงานซี 10) ปี 2552-ปัจจุบัน. </li>
<br />
<li>Committee Member of the Faculty Senate of Thammasat Univ 2002-now. กรรมการสภาอาจารย์มหาวิทยาลัยธรรมศาสตร์ 2545-ปัจจุบัน</li>
</span></ul>
<br />
กรรมการประเมินวิทยานิพนธ์ดีเ่ด่น มหาวิทยาลัยธรรมศาสตร์ พฤศจิกายน 2555 <br />
<br />
ผู้ประเมินบทความวิจัยระดับบัณฑิตศึกษา มหาวิทยาลัยเชียงใหม่ ครั้งที่ 3 (2555)<br />
<br />
กรรมการประเมินตำแหน่งทางวิชาการ ผู้ช่วยศาสตราจารย์ จุฬาลงกรณ์มหาวิทยาลัย 2554<br />
กรรมการประเมินตำแหน่งทางวิชาการ ผู้ช่วยศาสตราจารย์ มหาวิทยาลัยเทคโนโลยีราชามงคลตะวันออก 2555<br />
<br />
กรรมการตรวจสอบดุษฎีนิพนธ์ (doctoral thesis examiner) ของ Valeriu Dragan มหาวิทยาลัย Politehnica University of Bucharest (ประเทศ โรมาเนีย) 2012<br />
<br />
<br />
<span style="font-size: 180%; font-weight: bold;">International Invited Committee</span><br />
<ul style="font-family: georgia; font-size: 130%;"><br />
<li>2009 - The panel (including Assoc.Prof.Sukum Attavavutichai, Dr.Kanlaya Reuksuppasompon, Dr.Chalintorn N.Burian and Assoc.Prof.Dr.Boonsap Witchayangkoon) for the presentation of students from the Hong Kong universities participating in the Dragon Foundation's Global Citizenship Programme 2009 in Bangkok, Thailand, July 3, 2009.</li>
<br />
<li>2008 - Reviewer of applications for the Saudi Arabia's King Abdullah University of Science and Technology's Discovery Scholarship Award (the first batch of KAUST graduate students)</li>
</ul>
<br />
<span style="font-size: 180%; font-weight: bold;">Scholar Visitor</span><br />
<ul style="font-family: georgia;">
<li><span style="font-size: 130%;">Visiting Faculty (UMAP) at the <a href="http://www.ipc.shizuoka.ac.jp/%7Esemsato/satomurae.htm" target="_blank">Geoscience Institute, Shizuoka University</a>, Japan: March 15- June 10, 2006</span></li>
</ul>
<br />
<br />
<br />
<br />
<span style="font-size: 180%; font-weight: bold;">Fields of Interest</span><br />
Spatial Information Engineering - วิศวกรรมภูมิสารสนเทศ<br />
Informational Construction Management - สารสนเทศในการจัดการงานก่อสร้าง<br />
GPS/GIS & analysis -<br />
Transportation & Information of safety & accidents analysis<br />
Positioning of military weapons & cannon explosions<br />
Operations research in civil engineering<br />
Basic Research in Civil Engineering & related fields<br />
GPS-Aided Air Traffic Management<br />
GPS application PPP for studying of Troposphere and Ionopshere<br />
Engineering & Law - law-aided technology, environmental law, ..<br />
<br />
<span style="font-size: 180%; font-weight: bold;">On-Going Research Work</span><br />
<ul>
<li>GPS Supplementary Artificial Rainmaking Via Precipitable Water Vapor Information from Tropospheric Analysis</li>
<br />
<li>Traffic Calming System on Thammasat University - e.g. Analysis of Speed Bump & ..</li>
</ul>
<br />
<br />
<span style="font-size: 180%; font-weight: bold;">Graduate Students Research Group</span><br />
<ul>
<li>สัญญา นามี Sanya Namee, PhD Candidate</li>
<br />
<li>อภิรัฐ ฉัตรานุสรณ์ Master Candidate</li>
<br />
<li>Congratulations just graduate: P.Sritipun ภาคภูมิ ศรีตริพันธุ์ Master Candidate</li>
</ul>
<br />
<span style="font-size: 180%; font-weight: bold;">Publications</span><br />
<b>Book</b><br />
<ul>
<li>คู่มือ Easy Google SketchUp 8 (2554) </li>
<li>คู่มือ <a href="http://www.infopress.co.th/product_detail.php?pr_id=0000000616">Easy Google SketchUp 7</a> หรือ <a href="http://www.se-ed.com/eShop/%28A%28A8n6VrM4ywEkAAAANWE2MmRmYjctNjRmNy00MWNhLWFlNGQtZmVlMmJmYWZjMDVm_Nh3HWioSN96G95BhD49VIPadFo1%29%29/Products/Detail.aspx?CategoryId=0&No=9786162000782">ร้านซีเอ็ด</a></li>
</ul>
<br />
<b>Papers</b><br />
<ul><br />
<li>Namee, S., Witchayangkoon, B., and Karoonsoontawong, A. 2012. Fuzzy Logic Modeling Approach for Risk Area Assessment for Hazardous Materials Transportation. American Transactions on Engineering & Applied Sciences, 1(2): 127-142.</li>
<li>Chatranusorn, Aprirath., and Witchayangkoon, B. 2012. GPS Supplementary Artificial Raninmaking via Precipitable Water Vapor Information from Tropospheric Analysis. Thammasat Journal of Science & Technology, 20(2): 99-106. </li>
<li>Namee, S. and Witchayangkoon, B. 2011. Crossroads Vertical Speed Control Devices: Suggestion from Observation. International Transaction J of Engineering, Management, and Applied Sciences and Technologies, 2(2), 161-171.</li>
<li>Sirimontree, S., B.Witchayangkoon, N. Khaosri, and J.Teerawong. 2011. Shear Strength of Reinforced Concrete Bream Strengthened by Transverse External Post-tension. American J of Engineering and Applied Science, 4(1), 108-115.</li>
<li>Witchayangkoon, B., S.Namee, V.Temrungsri, and S.Intaratad. 2011. Traffic Calming at Crossroads Using Speed Bumps/Humps: A Case Study in Thammasat University, Rangsit Campus, Thailand. J of Science and Technology, 19(1) (Jan-Mar), 60-71. (in Thai language)</li>
<li>Witchayangkoon, B. and P. Dumrongchai. 2009. Computational Linearized Least Square for 2D Positioning Analysis Using MathCad Programming. <a href="http://www.kmutt.ac.th/rippc/journal48.htm"> </a><a href="http://www.kmutt.ac.th/rippc/journal48.htm" target="_blank">KMUTT Research and Development Journal.</a></li>
<li>Sritipun, P. and B. Witchayangkoon. 2009. GIS Database of Wind Speed for the central region of Thailand. <span style="font-style: italic;">Proc. the 14th National Convention on Civil Engineering</span>, Section of Survey and Geographic Information System Engineering (SGI), Nakorn-Ratchasrima, THAILAND, May 13-15, 2009.</li>
<li>Namee, S., and B. Witchayangkoon. 2008. An Application of Geographic Information System for Assessment Risk Area from Hazardous Materials Transport. <span style="font-style: italic;">Proc. the 5th National Transport Conference</span> (section Transport System & Sustainable Development Strategies), Thailand, Dec 19, 2008: 76-83.</li>
<li><span style="font-size: 100%;">Witchayangkoon, B. 2004. A GIS Database of Pile Dynamic Load Tests for the Bangkok Region. <span style="font-style: italic;">Proc. International Conference on Civil and Environmental Engineering (ICCEE)</span>, Higashi-Hiroshima, Japan.</span></li>
<li><span style="font-size: 100%;">Witchayangkoon, B. 2003. Estimation of Zenith Troposphere Using Dual-Frequency GPS Precise Point Position. <span style="font-style: italic;">Proc. the 4th Regional Symposium on Infrastructure Development in Civil Enigineering</span>, Bangkok, Thailand.</span></li>
<li><span style="font-size: 100%;">Wattanakul, N. and B.Witchayangkoon. 2002. "A survey and analysis of energy consumption of industrial factories within the central region of Thailand." <span style="font-style: italic;">Research and Development Journal of the Engineering Institute of Thailand</span>, Vol. 13(4): 11-15.</span></li>
<li><span style="font-size: 100%;">Witchayangkoon, B. "GPS-based vehicular velocity determination on the Chalermmahanakhon Expressway." <span style="font-style: italic;">Proc. the First Conference on Civil and Environmental Engineering (ICCEE-2002)</span>. Higashi-Hiroshima, Japan. Oct 30-31, 2002: 215-221.</span></li>
<li><span style="font-size: 100%;">Witchayangkoon, B. 2002. Roles of GIS-aided construction management. <span style="font-style: italic;">Proc. the 8th National Convention on Civil Engineering</span>. Hotel Sofitel Raja Orchid Khon Kaen, Thailand. Oct 23-25, 2002: CEM129-CEM134.</span></li>
<li><span style="font-size: 100%;">Tounnawarat, N., B.Witchayangkoon, and K. Chompooming. 2002. "GPS/GIS application for historical monuments around the Rattanagosin Island." <span style="font-style: italic;">Proc. the 8th National Convention on Civil Engineering</span>. Hotel Sofitel Raja Orchid Khon Kaen, Thailand. Oct 23-25.</span></li>
</ul>
<br />
<span style="font-size: 180%; font-weight: bold;">Courses Taught @ Faculty of Engineering @ Thammasat University</span><br />
Surveying<br />
Engineering Mechanics<br />
Strength of Material<br />
Geographic Information System (GIS)<br />
Global Positioning System (GPS)<br />
Research Methodology<br />
Applied information system<br />
Spatial information analysis<br />
<br />
<span style="font-size: 180%; font-weight: bold;">Professional Services</span><br />
<ul><br />
<li>Teachings and research advising for national & international institutes and universities</li>
<br />
<li>Consultancy assistance for many academics and professional projects.</li>
<br />
<li>Devoted advisor for students projects for building or repairing schools in the very rural areas of Thailand since 2001.</li>
<br />
<li>Professional Engineer (PE)</li>
<br />
<li>Contributing reader of various national and international journals</li>
<li>Energy & Environmental Conservation</li>
<br />
<li><span style="font-size: 130%;">อาจารย์หลักสูตร วิศวกรรมศาสตร์สองสถาบัน (TEP) คณะวิศวกรรมศาสตร์ ม.ธรรมศาสตร์</span></li>
</ul>
<br />
<a href="http://www.ce.engr.tu.ac.th/ce/eng/staff_view.php?company_code=145" target="_blank">This is a mirror homepage of Dr.Boonsap Witchayangkoon official web site @ Faculty of Engineering Thammasat University @ THAILAND<br /><br />เวปนี้ เป็นเวปสำรองของ รองศาสตราจารย์ ดร.บุญทรัพย์ วิชญางกูร ภาควิชาวิศวกรรมโยธา คณะวิศวกรรมศาสตร์ มหาวิทยาลัยธรรมศาสตร์<br /><span style="font-size: 85%;">http://www.ce.engr.tu.ac.th/ce/eng/staff_view.php?company_code=145</span></a><br />
<br />
<br />
<span style="font-size: 180%;"><span style="font-weight: bold;">THAILAND : We are here..</span></span><br />
<span style="font-family: verdana; font-size: 100%;">For the popular international engineering education, we offer the attractive<br /><a href="http://www.tep.engr.tu.ac.th/" target="_blank">Twining Engineering Program (TEP) and Thammasat English Program of Engineering (TEPE)</a>. Our programs are worldwide recognized to international standard of top-rank Thailand engineering education.<br /></span>ดร.บุญทรัพย์ วิชญางกูรhttp://www.blogger.com/profile/09210204627753495612noreply@blogger.comtag:blogger.com,1999:blog-1366818990289389328.post-65037783512583444052009-05-06T04:12:00.000-07:002009-05-11T21:44:57.863-07:00Extended Kalman Filter for GPS PPPfunction [Result,Residual,Relativity,Pminus,CLKphi,NumSV,W,Pkalman] ... <br /> = EKF(DOY,rTime,Ph1,Ph2,Pr1,Pr2,SP3X,SP3Y,SP3Z,SP3clk,Vx,Vy,Vz,x,Elev,trop,type) <br /> <br />% Run Extend Kalman Filtering <br />% <br />% %%%%%%%%% HELP %%%%%%%%%% <br />% <br />% [Result,Residual,Relativity,Pminus,CLKphi,NumSV,W,Pkalman] ... <br />% = EKF(DOY,rTime,Ph1,Ph2,Pr1,Pr2,SP3X,SP3Y,SP3Z,SP3clk,Vx,Vy,Vz,x,Elev,trop,type) <br />% <br />% INPUT <br />% <br />% DOY = Day of year [1 x 1] <br />% rTime = GPS second of RINEX file [GPSsec], [numep x 1] <br />% Ph1 = Carrier phase L1,cycles [numep x numsat] <br />% Ph2 = Carrier phase L2,cycles [numep x numsat] <br />% Pr1 = Pseudorange C/A code, metres [numep x numsat] <br />% Pr2 = Pseudorange P code , metres [numep x numsat] <br />% SP3X, SP3Y, SP3Z = Satellite coordinate (WGS84), metres [numep x numsat] <br />% SP3clk = Precise clock correction, microseconds [numep x numsat] <br />% Vx = Satellite velocity [numep x numsat], m/s <br />% Vy = Satellite velocity [numep x numsat], m/s <br />% Vz = Satellite velocity [numep x numsat], m/s <br />% x = Aproximate station coordinates, metres [3 x 1] <br />% Elev = Elevation angle, radians [numep x numsat] <br />% trop = 1 Estimate tropospheric correction from EKF <br />% = 2 Approximate tropospheric correction from Saastamoinen model <br />% = 3 Approximate tropospheric correction from Saastamoinen model & Neill mapping function <br />% type = 1 Normal <br />% = 2 transition matrix (Phi) is identity <br />% <br />% OUTPUT <br />% <br />% Result = [x y z Ni dT Tropo PDOP TDOP GDOP] <br />% Residual <br />% Relativity <br />% Pminus <br />% CLKphi <br />% NumSV <br />% W <br />% Pkalman <br />% <br />% Modified by Phakphong Homniam <br />% December 25, 2002 <br />% Orginal source code in Mathcad by Boonsap Witchayangkoon, 2000 <br /> <br />%%%%%%%%%% CONSTANT DECLARATION %%%%%%%%%% <br /> <br />% Global constants of the relevant parameters in constant_global.m file <br />constant_global <br /> <br />global NumSat NumPara gamma maskangle <br /> <br />% Satellite mask angle <br />maskangle = 15*deg; <br /> <br />% Kalman initialization constants <br />SDx = 0.0003; % (m) used in Q (dynamic model) <br />SDclk = 30; % (m) usde in Q (dynamic model) <br />SDph = 0.02; % (m) STDS of carrier phase (R) <br />SDpr = 0.2; % (m) STD of pseudorange (R) <br />SDcoorinit = 1; % (m) STD of initial parameters (P) <br />SDambinit = 10^10; % (m) STD of initial ambiguities (P) <br />SDclkinit = 500; % (m) STD for initial clock (P) <br />SDtropinit = 0.1; % (m) STD for initial tropo (P) <br /> <br />% Specification when estimating troposphere <br />% Tropospheric Correlation Time <br />thor = 18000; % seconds <br /> <br />% Number of satellite <br />NumSat = size(Ph1,2); <br /> <br />gamma = (f1/f2)^2; <br />tep = length(rTime); <br />TimeInterval = rTime(3)-rTime(2); <br /> <br />%%%%%%%%%% BEGIN %%%%%%%%%% <br /> <br />Result = []; <br />Residual = []; <br />Relativity = []; <br />Pminus = []; <br />CLKphi = []; <br />NumSV = []; <br />W = []; <br />Pkalman = []; <br /> <br />% Estimate receiver clock error for the first epoch for EKF <br />for i = 1:NumSat <br /> if (Ph1(1,i) ~= NA) & (Pr1(1,i) ~= NA) <br /> if SP3clk(1,i) ~= NAJ <br /> rho = sqrt((SP3X(1,i)-x(1))^2+(SP3Y(1,i)-x(2))^2 ... <br /> +(SP3Z(1,i)-x(3))^2); <br /> h = rho-c*SP3clk(1,i)*10^-6; <br /> FirstClkEst = Pr1(1,i)-h; <br /> break <br /> end <br /> end <br />end <br /> <br />switch trop <br />case 1 <br /> % Number of unknown parameters <br /> % (x,y,z)+(ambiguities)+Receiver Clock+Tropo <br /> NumPara = 3+NumSat+2; <br /> <br /> % Set martix Xm <br /> zenith = 2.3; <br /> Xstart = [x(:);zeros(NumSat,1);FirstClkEst;zenith]; <br /> Xm = Xstart; <br /> Xminus = Xm; <br /> <br /> % Set matrix Q <br /> qSD(1:3) = SDx; <br /> qSD(NumSat+4) = SDclk; <br /> qSD(NumSat+5) = 0.01/3600; % may be 0.1/60 or 0.01/3600 : i don't sure <br /> Q = diag(qSD.^2); <br /> <br /> % Set matrix Pminus <br /> pmSD(1:3) = SDcoorinit; <br /> pmSD(4:NumSat+3) = SDambinit; <br /> pmSD(NumSat+4) = SDclkinit; <br /> pmSD(NumSat+5) = exp(-TimeInterval/thor); <br /> Pm = diag(pmSD.^2); <br /> Pminus = Pm; <br /> <br />case {2,3} <br /> % Number of unknown parameters <br /> % (x,y,z)+(ambiguities)+Receiver Clock <br /> NumPara = 3+NumSat+1; <br /> <br /> % Set martix Xm <br /> Xstart = [x(:);zeros(NumSat,1);FirstClkEst]; <br /> Xm = Xstart; <br /> Xminus = Xm; <br /> <br /> % Set matrix Q <br /> qSD(1:3) = SDx; <br /> qSD(NumSat+4) = SDclk; <br /> Q = diag(qSD.^2); <br /> <br /> % Set matrix Pminus <br /> pmSD(1:3) = SDcoorinit; <br /> pmSD(4:NumSat+3) = SDambinit; <br /> pmSD(NumSat+4) = SDclkinit; <br /> Pm = diag(pmSD.^2); <br /> Pminus = Pm; <br />end <br /> <br />% Set matrix Phi <br />switch type <br />case 1 <br /> Phi = eye(NumPara); <br /> % Set Phi for the receiver clock <br /> % = 0 pure white noise <br /> % = 1 random walk <br /> Phi(NumSat+4,NumSat+4) = 0; <br />case 2 <br /> Phi = eye(NumPara); <br />end <br />[PR,Phase] = typeobs(Ph1,Ph2,Pr1,Pr2,tep); <br />ep = 1; <br /> <br />for i = 1:(tep-1) <br /> phase = Phase(i,:); <br /> pr = PR(i,:); prp1 = PR(i+1,:); <br /> sp3x = SP3X(i,:); sp3xp1 = SP3X(i+1,:); vx = Vx(i,:); <br /> sp3y = SP3Y(i,:); sp3yp1 = SP3Y(i+1,:); vy = Vy(i,:); <br /> sp3z = SP3Z(i,:); sp3zp1 = SP3Z(i+1,:); vz = Vz(i,:); <br /> sp3clk = SP3clk(i,:); sp3clkp1 = SP3clk(i+1,:); elev = Elev(i,:); <br /> <br /> % set matrix H <br /> [HhZ1,HhZ2,HhZ3,HhZ4] = MatrixHhZ(Xm,phase,pr,sp3x,sp3y,sp3z,sp3clk,elev,x,DOY,vx,vy,vz,trop); <br /> H = HhZ1; <br /> <br /> % Set matrix R <br /> rSD(1:NumSat) = SDph; <br /> rSD(NumSat+1:2*NumSat) = SDpr; <br /> R = diag(rSD.^2); % Observation eorror <br /> % Reweight observation at low elevation angles <br />% for j = 1:NumSat <br />% if (elev(j) ~= NA) & (elev(j) < 30*deg) <br />% if elev(j) < 25*deg <br />% lv = elev(j)+5*deg; <br />% else <br />% lv = elev(j); <br />% end <br />% R(j,j) = (SDph/sin(lv))^2; <br />% R(j+NumSat,j+NumSat) = (SDpr/sin(lv))^2; <br />% end <br />% end <br /> <br /> K = Pm*H'*inv(H*Pm*H'+R); % Calculate Kalman Gain <br /> w = K*HhZ2; <br /> W = [W,w]; <br /> Xestimate = Xm+w; % Using measurement to update estimate <br /> P = Pm-K*H*Pm; % Update the error covariance <br /> Pkalman = [Pkalman;P]; <br /> % Phi clk here <br /> Phi(NumSat+4,NumSat+4) ... % Xestimate or x??? <br /> = clkphiall(Xestimate,pr,sp3x,sp3y,sp3z,sp3clk,prp1,sp3xp1,sp3yp1,sp3zp1,sp3clkp1); <br /> Xm = Phi*Xestimate; % Project the state ahead <br /> Xminus = [Xminus,Xm]; <br /> <br /> switch type <br /> case 1 <br /> Pm = PhiPPhiT(Phi,P,trop)+Q; % Project the covariance ahead <br /> case 2 <br /> Pm = P+Q; <br /> end <br /> <br /> Pminus = [Pminus;Pm]; <br /> Result(1:NumPara,ep) = Xestimate; <br /> A = H4DOPs(phase,pr,sp3x,sp3y,sp3z,elev,x); <br /> ATA = A'*A; <br /> detATA = det(ATA); <br /> if detATA == 0 <br /> PDOP = NA; <br /> TDOP = NA; <br /> GDOP = NA; <br /> else <br /> Pvar = inv(ATA); <br /> PDOP = sqrt(Pvar(1,1)+Pvar(2,2)+Pvar(3,3)); <br /> TDOP = sqrt(Pvar(4,4)); <br /> GDOP = sqrt(Pvar(1,1)+Pvar(2,2)+Pvar(3,3)+Pvar(4,4)); <br /> end <br /> Result(NumPara+1,ep) = PDOP; <br /> Result(NumPara+2,ep) = TDOP; <br /> Result(NumPara+3,ep) = GDOP; <br /> Residual = [Residual,HhZ2]; <br /> NumSV(ep) = HhZ3; <br /> Relativity = [Relativity,HhZ4(:)]; <br /> CLKphi(ep) = Phi(NumSat+4,NumSat+4); <br /> ep = ep+1; <br />end <br /> <br />%%%%%%%%%% RESULT %%%%%%%%%% <br /> <br />Result = Result'; <br />Residual = Residual'; <br />Relativity = Relativity'; <br />Pminus; <br />CLKphi = CLKphi(:); <br />NumSV = NumSV(:); <br />W = W'; <br />Pkalman; <br />%%%%%%%%%% END %%%%%%%%%%ดร.บุญทรัพย์ วิชญางกูรhttp://www.blogger.com/profile/09210204627753495612noreply@blogger.comtag:blogger.com,1999:blog-1366818990289389328.post-8420856809521697612009-05-04T21:24:00.000-07:002009-05-22T09:55:20.515-07:00Neill Mapping Function Source Code for GPS/ GPS PPP/ GNSSfunction hmf = NeillMF(DOY,lat,h,elev) <br /> <br />% Neill Mapping Function <br />% <br />% %%%%%%%%%% HELP %%%%%%%%%% <br />% <br />%hmf = NeillMF(DOY,lat,h,elev) <br />% <br />% Input Data <br />% DOY is day of year. <br />% lat is latitude, degrees. <br />% long is longitude, degrees. <br />% h is height, meters. <br />% elev is elevation angle, degrees. <br />% <br />% Transwritten by Phakphong Homniam <br />% September 20, 2002 <br />% Original Mathcad source code by Boonsap Witchayangkoon, 2000 <br /> <br />% Coefficient for changing DOY to radian <br />day2rad = 2*pi/365.25; <br /> <br />% Coefficient for changing degrees to radians <br />deg2rad = pi/180; <br /> <br />% Latitude array for the hydrostatic mapping function coefficients <br />lat_hmf = [15,30,45,60,75]'; <br /> <br />% Average of coefficients a, b and c corresponding to the given latitude. <br />% See Table 5.4 <br />abc_avg = [1.2769934 2.9153695 62.610505 <br /> 1.2683230 2.9152299 62.837393 <br /> 1.2465397 2.9288445 63.721774 <br /> 1.2196049 2.9022565 63.824265 <br /> 1.2045996 2.9024912 64.258455]*10^-3; <br /> <br />% Amplitude of coefficients a, b and c corresponding to the given latitude. <br />abc_amp = [0 0 0 <br /> 1.2709626 2.1414979 9.0128400 <br /> 2.6523662 3.0160779 4.3497037 <br /> 3.4000452 7.2562722 84.795348 <br /> 4.1202191 11.723375 170.37206]*10^-5; <br /> <br />% Height correction <br />a_ht = 2.53*10^-5; <br />b_ht = 5.49*10^-3; <br />c_ht = 1.14*10^-3; <br /> <br />% Wet Mapping Function <br />% See Table5.5 <br />lat_wmf = lat_hmf; <br />abc_w2po = [5.8021897*10^-4 1.4275268*10^-3 4.3472961*10^-2 <br /> 5.6794847*10^-4 1.5138625*10^-3 4.6729510*10^-2 <br /> 5.8118019*10^-4 1.4572752*10^-3 4.3908931*10^-2 <br /> 5.9727542*10^-4 1.5007428*10^-3 4.4626982*10^-2 <br /> 6.1641693*10^-4 1.7599082*10^-3 5.4736038*10^-2]; <br /> <br />hs_km = h/1000; <br />DOY_atm = DOY-28; <br />if lat < 0 <br /> DOY_atm = DOY_atm+365.25/2; <br /> lat = abs(lat); <br />end <br />DOYr_atm = DOY_atm*day2rad; <br />cost = cos(DOYr_atm); <br />if lat <= lat_hmf(1,1) <br /> a = abc_avg(1,1); <br /> b = abc_avg(1,2); <br /> c = abc_avg(1,3); <br />end <br />for i = 1:4 <br /> if (lat >= lat_hmf(i,1)) & (lat <= lat_hmf(i+1,1)) <br /> dlat = (lat-lat_hmf(i))/(lat_hmf(i+1)-lat_hmf(i)); <br /> daavg = abc_avg(i+1,1)-abc_avg(i,1); <br /> daamp = abc_amp(i+1,1)-abc_amp(i,1); <br /> aavg = abc_avg(i,1)+dlat*daavg; <br /> aamp = abc_amp(i,1)+dlat*daamp; <br /> a = aavg+aamp*cost; % + or - i don't sure <br /> dbavg = abc_avg(i+1,2)-abc_avg(i,2); <br /> dbamp = abc_amp(i+1,2)-abc_amp(i,2); <br /> bavg = abc_avg(i,2)+dlat*dbavg; <br /> bamp = abc_amp(i,2)+dlat*dbamp; <br /> b = bavg+bamp*cost; % + or - i don't sure <br /> dcavg = abc_avg(i+1,3)-abc_avg(i,3); <br /> dcamp = abc_amp(i+1,3)-abc_amp(i,3); <br /> cavg = abc_avg(i,3)+dlat*dcavg; <br /> camp = abc_amp(i,3)+dlat*dcamp; <br /> c = aavg+aamp*cost; % + or - i don't sure <br /> end <br />end <br />if lat >= lat_hmf(5,1) <br /> a = abc_avg(5,1)+abc_amp(5,1)*cost; % + or - i don't sure <br /> b = abc_avg(5,2)+abc_amp(5,2)*cost; % + or - i don't sure <br /> c = abc_avg(5,3)+abc_amp(5,3)*cost; % + or - i don't sure <br />end <br />sine = sin(elev*deg2rad); <br />cose = cos(elev*deg2rad); <br /> <br />beta = b/(sine+c); <br />gamma = a/(sine+beta); <br />topcon = 1+(a/(1+b/(1+c))); <br />hmf(1,1) = topcon/(sine+gamma); <br />hmf(2,1) = -topcon*cose/((sine+gamma)^2*... <br /> (1-a/((sine+beta)^2*(1-b/(sine*c)^2)))); <br /> <br />beta = b_ht/(sine+c_ht); <br />gamma = a_ht/(sine+beta); <br />topcon = 1+a_ht/(1+b_ht/(1+c_ht)); <br />ht_corr_coef = 1/sine-topcon/(sine+gamma); <br />ht_corr = ht_corr_coef*hs_km; <br />hmf(1,1) = hmf(1,1)+ht_corr; <br />dhcc_del = -cose/sine^2+topcon*cose/((sine+gamma)^2*... <br /> (1-a_ht/((sine+beta)^2*(1-b_ht/(sine*c_ht)^2)))); <br />dht_corr_del = dhcc_del*hs_km; <br />hmf(2,1) = hmf(2,1)+dht_corr_del; <br />dht_corr_del = dhcc_del*hs_km; <br />hmf(2,1) = hmf(2,1)+dht_corr_del; <br />if lat <= lat_wmf(1,1) <br /> alat = abc_w2po(1,1); <br /> blat = abc_w2po(2,1); <br /> clat = abc_w2po(3,1); <br />end <br />for i = 1:4 <br /> if (lat >= lat_wmf(i,1)) & (lat <= lat_wmf(i+1,1)) <br /> dll = (lat-lat_wmf(i,1))/(lat_wmf(i+1,1)-lat_wmf(i,1)); <br /> da = abc_w2po(i+1,1)-abc_w2po(i,1); <br /> alat = abc_w2po(i,1)+dll*da; <br /> db = abc_w2po(i+1,2)-abc_w2po(i,2); <br /> blat = abc_w2po(i,2)+dll*db; <br /> dc = abc_w2po(i+1,3)-abc_w2po(i,3); <br /> clat = abc_w2po(i,3)+dll*dc; <br /> end <br />end <br />if lat >= lat_wmf(5,1) <br /> alat = abc_w2po(5,1); <br /> blat = abc_w2po(5,2); <br /> clat = abc_w2po(5,3); <br />end <br />sinelat = sin(elev*deg2rad); <br />coselat = cos(elev*deg2rad); <br />betalat = blat/(sinelat+clat); <br />gammalat = alat/(sinelat+betalat); <br />topconlat = 1+alat/(1+blat/(1+clat)); <br />wmf(1,1) = topconlat/(sinelat+gammalat); <br />wmf(2,1) = -topconlat/((sinelat+gammalat)^2*(coselat-alat/... <br /> ((sinelat+betalat)^2*coselat*(1-blat/(sinelat*clat)^2)))); <br />hmf(3,1) = wmf(1,1); <br />hmf(4,1) = wmf(2,1); <br /> <br />hmf; <br /> <br />%%%%%%%%%%END%%%%%%%%%%ดร.บุญทรัพย์ วิชญางกูรhttp://www.blogger.com/profile/09210204627753495612noreply@blogger.comtag:blogger.com,1999:blog-1366818990289389328.post-35441935529496230992009-05-03T21:48:00.000-07:002009-05-11T21:54:36.448-07:00Lagrange for SP3 Satellite Orbit and Clock Correction Interpolationfunction [xsp3,ysp3,zsp3,clksat,P1,P2,sp3Vx,sp3Vy,sp3Vz ...<br /> ,RV1,RV2,RVIO,saas,Elev,AZIMUTH,PH1RV,PH2RV,PHRVIO] ...<br /> = LG_SP3(DOY,rT,PH1,PH2,PR1,PR2,SP3time,SP3X,SP3Y,SP3Z,SATclk,x,SunX,SunY,SunZ,PTRH)<br /><br />% Lagrange:satellite orbit and satellite clock interpolation <br />% called by EKF. Computes PHI (dynamic model) for receiver clock<br />% <br />% %%%%%%%%%% HELP %%%%%%%%%%<br />% <br />% [xsp3,ysp3,zsp3,clksat,P1,P2,sp3Vx,sp3Vy,sp3Vz ...<br />% ,RV1,RV2,RVIO,saas,Elev,AZIMUTH,PH1RV,PH2RV,PHRVIO] ...<br />% = LG_SP3(DOY,rT,PH1,PH2,PR1,PR2,SP3time,SP3X,SP3Y,SP3Z,SATclk,x,SunX,SunY,SunZ,PTRH)<br />% <br />% Input data<br />% DOY = day of year<br />% rT = nominal (RINEX) time for epoch ep [GPSweek GPSsec], [numep x 2]<br />% PH1 = Carrier phase L1,cycles [numep+1 x numsat]<br />% PH2 = Carrier phase L2,cycles [numep+1 x numsat]<br />% PR1 = Pseudorange C/A code, metres [numep+1 x numsat]<br />% PR2 = Pseudorange P code , metres [numep+1 x numsat]<br />% SP3time = GPS time of SP3 file [GPSweek GPSsec], [numep x 2]<br />% SP3X, SP3Y, SP3Z = Satellite coordinate (WGS84), metres [numep+1 x numsat]<br />% x = Aproximate station coordinates, metres [3 x 1]<br />% SATclk = Precise clock correction, microseconds [numep+1 x numsat]<br />% SunX, SunY, SunZ = Sun coordinate (WGS84), metres [numep+1 x numsat]<br />% PTRH = [P TC RH]<br />%<br />% OUTPUT<br />% xsp3<br />% ysp3<br />% zsp3<br />% clksat<br />% P1<br />% P2<br />% sp3Vx<br />% sp3Vy<br />% sp3Vz<br />% RV1<br />% RV2<br />% RVIO<br />% saas<br />% Elev<br />% AZIMUTH<br />% PH1RV<br />% PH2RV<br />% PHRVIO<br />%<br />% Modified by Phakphong Homniam<br />% October 13, 2002<br />% Original Mathcad source code by Boonsap Witchayangkoon, 2000<br /><br />%%%%%%%%%% CONSTANT DECLARATION %%%%%%%%%%<br /><br />% Global constants of the relevant parameters in constant_global.m file<br />constant_global<br /><br />global maskangle antenna_offset IIR <br /><br />% Satellite mask angle<br />maskangle = 15*deg;<br /><br />% Satellite Antenna Offset<br />antenna_offset = [0.279 0 1.023]';<br /><br />% GPS Block IIR satellite<br />IIR = 13;<br />% For GPS Block IIR satellite, IGS has agreed upon to<br />% the satellite phase center offset [0 0 0]<br /><br />% Langrange InterPolation Constants<br />% Number of points used in LG orbit interpolation<br />nXYZ = 9;<br />% Number of points used in LG clock interpolation<br />nCLK = 5;<br /><br />% Environmental data for this station<br />P = PTRH(1); % Pressure P in milibar (mb)<br />TC = PTRH(2); % Temperature TC in degree CELCIUS<br />RH = PTRH(3); % Humedity H(%)<br />T = TC+273.2; % Temperature T in Kelvin<br /><br />%%%%%%%%%% BEGIN %%%%%%%%%%<br /><br />PRN = PH1(1,:);<br />PH1(1,:) = []; SP3X(1,:) = []; <br />PH2(1,:) = []; SP3Y(1,:) = []; <br />PR1(1,:) = []; SP3Z(1,:) = []; <br />PR2(1,:) = []; SATclk(1,:) = []; <br />rT = rT(:,2); SP3time = SP3time(:,2);<br /><br />x = x(:);<br />llh = ecef2lla(x',1);<br />wetzenith = SWZD(P,T,RH);<br />dryzenith = SHZD(P,T,RH);<br />[tep,NumSat] = size(PR1);<br />for j = 1:NumSat<br /> fprintf('\nProcessing satellite number %s\n',num2str(PRN(j)));<br /> CLKsat(1:tep,1) = NAJ;<br /> sp3x(1:tep,1) = NAJ;<br /> sp3y(1:tep,1) = NAJ;<br /> sp3z(1:tep,1) = NAJ;<br /> Prange1(1:tep,1) = NAJ;<br /> Prange2(1:tep,1) = NAJ;<br /> Vx(1:tep,1) = NAJ;<br /> Vy(1:tep,1) = NAJ;<br /> Vz(1:tep,1) = NAJ;<br /> R1(1:tep,1) = NA;<br /> R2(1:tep,1) = NA;<br /> RIO(1:tep,1) = NA;<br /> Tropo(1:tep,1) = NA;<br /> ELEV(1:tep,1) = NA;<br /> azimuth(1:tep,1) = NA;<br /> PH1R(1:tep,1) = NA;<br /> PH2R(1:tep,1) = NA;<br /> PHRIO(1:tep,1) = NA;<br /> <br /> r = size(IIR,1);<br /> % Assume no satellite offset<br /> NoOffset = 1;<br />% NoOffset = 0;<br />% for h = 1:r<br />% if PRN(j) == IIR(h)<br />% NoOffset = 1;<br />% end<br />% end<br /> for i = 1:tep<br />% i<br /> rT(i);<br /> if (PR1(i,j) ~= NA) & (PR2(i,j) ~= NA) & (PH1(i,j) ~= NA) & (PH2(i,j) ~= NA)<br /> CLKsat(i) = SVCLK(rT(i),SP3time,SATclk(:,j) ...<br /> ,PR1(i,j),nCLK,NAJ);<br /> if CLKsat(i) ~= NAJ<br /> [X,V] = SatPo(rT(i),PR1(i,j),CLKsat(i) ...<br /> ,SP3time,SP3X(:,j),SP3Y(:,j),SP3Z(:,j) ...<br /> ,x,nXYZ,[SunX(i),SunY(i),SunZ(i)]',NoOffset,NAJ);<br /> rho = sqrt((X(1)-x(1))^2+(X(2)-x(2))^2+(X(3)-x(3))^2);<br /> elev = elevation(X,x);<br /> if elev < maskangle<br /> continue<br /> end<br /> azimuth(i) = Az(llh(1),llh(2),X-x);<br /> ELEV(i) = elev;<br /> sp3x(i) = X(1);<br /> sp3y(i) = X(2);<br /> sp3z(i) = X(3);<br /> Vx(i) = V(1);<br /> Vy(i) = V(2);<br /> Vz(i) = V(3);<br /> nmf = NeillMF(DOY,llh(1)/deg,llh(3),elev/deg);<br /> Tropo(i) = dryzenith*nmf(1)+wetzenith*nmf(2);<br /> Relativity = -2*dot(X(:),V(:))/c;<br /> h = rho-c*CLKsat(i)*10^-6+Tropo(i)-Relativity;<br /> R1(i) = PR1(i,j)-h;<br /> R2(i) = PR2(i,j)-h;<br /> % See eq.9.46<br /> RIO(i) = cf1*PR1(i,j)-cf2*PR2(i,j)-h;<br /> PH1R(i) = lamda1*PH1(i,j)-h;<br /> PH2R(i) = lamda2*PH2(i,j)-h;<br /> % See eq.9.51<br /> PHRIO(i) = lamda1*(cf1*PH1(i,j)-cf1f2*PH2(i,j))-h;<br />% if (abs(R1(i)) < rclock) & (abs(R2(i)) < rclock)<br /> Prange1(i) = PR1(i,j);<br /> Prange2(i) = PR2(i,j);<br />% end<br /> end<br /> end<br /> end<br /> xsp3(1:tep,j) = sp3x(:);<br /> ysp3(1:tep,j) = sp3y(:);<br /> zsp3(1:tep,j) = sp3z(:);<br /> clksat(1:tep,j) = CLKsat(:);<br /> P1(1:tep,j) = Prange1(:);<br /> P2(1:tep,j) = Prange2(:);<br /> sp3Vx(1:tep,j) = Vx(:);<br /> sp3Vy(1:tep,j) = Vy(:);<br /> sp3Vz(1:tep,j) = Vz(:);<br /> RV1(1:tep,j) = R1(:);<br /> RV2(1:tep,j) = R2(:);<br /> RVIO(1:tep,j) = RIO(:);<br /> saas(1:tep,j) = Tropo(:);<br /> Elev(1:tep,j) = ELEV(:);<br /> AZIMUTH(1:tep,j) = azimuth(:);<br /> PH1RV(1:tep,j) = PH1R(:);<br /> PH2RV(1:tep,j) = PH2R(:);<br /> PHRVIO(1:tep,j) = PHRIO(:);<br />end<br /><br />%%%%%%%% RESULT %%%%%%%%<br /><br />xsp3 = [PRN;xsp3];<br />ysp3 = [PRN;ysp3];<br />zsp3 = [PRN;zsp3];<br />clksat = [PRN;clksat];<br />P1 = [PRN;P1];<br />P2 = [PRN;P2];<br />sp3Vx = [PRN;sp3Vx];<br />sp3Vy = [PRN;sp3Vy];<br />sp3Vz = [PRN;sp3Vz];<br />RV1 = [PRN;RV1];<br />RV2 = [PRN;RV2];<br />RVIO = [PRN;RVIO];<br />saas = [PRN;saas];<br />Elev = [PRN;Elev];<br />AZIMUTH = [PRN;AZIMUTH];<br />PH1RV = [PRN;PH1RV];<br />PH2RV = [PRN;PH2RV];<br />PHRVIO = [PRN;PHRVIO];<br /><br />%%%%%%%%%%END%%%%%%%%%%ดร.บุญทรัพย์ วิชญางกูรhttp://www.blogger.com/profile/09210204627753495612noreply@blogger.comtag:blogger.com,1999:blog-1366818990289389328.post-88167818212147395512008-12-19T07:50:00.000-08:002008-12-19T08:01:25.333-08:00QZSS: A new asset for PPP on Pacific Ocean of AsiaJapanese GPS-augmented navigation satellites on Pacific ocean of Asia called "Quasi-Zenith Satellite System (QZSS)" will give a new challenge to gain positioning accuracy within shorter time. This is true not only single point positioning, differential positioning, but also precise point positioning (PPP). Similar to GPS PPP, precise orbit and clock of QZSS must be made available for the augmentation analysis.<br /><br />We are very much looking forward to this new positioning space asset.ดร.บุญทรัพย์ วิชญางกูรhttp://www.blogger.com/profile/09210204627753495612noreply@blogger.comtag:blogger.com,1999:blog-1366818990289389328.post-58943221349556608592008-12-18T06:11:00.000-08:002009-05-03T19:55:39.582-07:00ELEMENTS OF GPS PRECISE POINT POSITIONING<span style="font-size:180%;"><em>ELEMENTS OF GPS PRECISE POINT POSITIONING</em></span><br />PhD Dissertation By Boonsap Witchayangkoon.<br />PhD Thesis Adviser: <a href="http://www.gnss.umaine.edu/9_contact.htm" target="_blank">Professor Dr. Alfred Leick</a>.<br />You can download freely from<br /><a href="http://www.cctechnol.com/uploads/167_0.pdf" target="_blank">http://</a><cite><a href="http://www.cctechnol.com/uploads/167_0.pdf">www.cctechnol.com/uploads/167_0.pdf</a><br /><br /><span style="font-size:150%;"><em>*** My dissertation work has been mentioned in the book<br /><a href="http://books.google.co.th/books%3Fid%3D4qE6xYjYSHgC%26pg%3DPA255%26lpg%3DPA255%26dq%3Dwitchayangkoon%26source%3Dbl%26ots%3Dfp9Ex2H4TF%26sig%3DiBz-fHlIKiucg48-EG9S5HfEtC0%26hl%3Dth&ei=XHr6SeTCOtSUkAXo9b2CBQ&sig2=giTV0r6HWti3VRcS6VQNIA&usg=AFQjCNENz_E2jt1UuyZqWZZJijA9zuJchw&ei=XHr6SeTCOtSUkAXo9b2CBQ&oi=book_result&ct=result&resnum=6"></a><a href="http://books.google.co.th/books?id=4qE6xYjYSHgC&printsec=frontcover">GPS Satellite Surveying (see page 254-256).</a><br /></em></span><br /><br /><br /><br /></cite>ดร.บุญทรัพย์ วิชญางกูรhttp://www.blogger.com/profile/09210204627753495612noreply@blogger.comtag:blogger.com,1999:blog-1366818990289389328.post-22551452461481551592008-12-18T05:57:00.000-08:002008-12-22T00:38:18.803-08:00GPS Precise Point Positioning (PPP)<span style="font-size:130%;"><span style="font-family: georgia;font-family:webdings;" >There are many more researcher paid attention to PPP. There are many claimed result reached only at a meter level for static positioning. This is very much misleading. This is mainly because PPP should attain centimeters accuracy, not a meter.</span><br /></span>ดร.บุญทรัพย์ วิชญางกูรhttp://www.blogger.com/profile/09210204627753495612noreply@blogger.com