当前位置:首页 >> 天文/地理 >>

p2


INTEGRATION OF INSAR AND GPS FOR VERTICAL DEFORMATION MONITORING : A CASE STUDY IN FAIAL AND PICO ISLANDS
(1)

Jo?o Catal?o(1), Giovanni Nico(1), Ramon Hanssen(2) and Cristina Catita(1) IDL, LATTEX, Univ. Lisboa, , 1749-016 Lisboa, Portugal, Email:jcfernandes@fc.ul.pt (2) Delft University of Technology, Kluyverweg 1, 2629 HS, Delft, PO Box 5058, 2600 GB Delft, Netherlands, Email: R.F.Hanssen@tudelft.nl network is not enough dense to strudy how the deformation is distributed. For that, a set of ASAR images (ascending and descending passes) were acquired between 2006 and 2009 in the scope of ESA Cat-1 (n. 3149) project and a map of the line-of-sight deformation was derived using the Persistent Scatterer approach. However, Faial and Pico are high vegetated islands with very unstable atmospheric conditions influencing negatively the interferometric processing. Also, the manmade structures are sparse reducing the density of Persistent Scatterers. Because of different geometry acquisitions, ascending and descending passes produce different sets of Persistent Scatterers, with complementary spatial coverage, and different line-ofsight velocities. Besides, the estimated velocities are relative to the master image (different from ascending and descending) and must be referred to an absolute velocity (in the sense of referred to a geodetic reference frame).

ABSTRACT In this work we investigate the integration of repeated GPS geodetic measurements and interferometric observations of synthetic aperture radar images from ascending and descending passes, aimed at the determination of a high resolution vertical deformation map. For that, a set of 57 ASAR-ENVISAT images (ascending and descending passes) were acquired between 2006 and 2009 over Faial and Pico islands in the scope of ESA Cat-1 (n. 3149) project. Two sets of interferograms were computed relative to ascending and descending masters SAR images. Each interferogram was corrected for atmospheric artifacts using a Weather Forecasting model. Initially, the horizontal component of the velocity (east and north) is assigned to each PS from interpolation of available GPS observations. Then, the vertical component of the velocity is determined, for each PS, from the line-of-sight velocity and the horizontal component of the velocity. Latter the vertical offsets are numerically determined by comparison between GPS (ITRF velocities) and the two sets of PS (ascending and descending). These values are then used to create the vertical deformation map of Faial and Pico islands with considerably better resolution and accuracy than single set of observations. The vertical deformation map has identified a large continuous area of subsidence on the west of Faial Island, on the flank of Capelinhos eruption cone, with a maximum subsidence range of 10 mm/yr. It has also revealed the subsidence of the summit crater of Pico Island (9 mm/yr) and a large area of subsidence on the west of the island, corresponding mostly to creep movement. 1. INTRODUCTION

Faial and Pico are two islands of Azores archipelago. Faial island has been, in the last century, affected by intense volcanic and seismic activity. On 1957 a volcanic eruption started on the west of the island increasing the area of Faial island lasting one year (Capelinhos eruption, 1957-58) [1] and on 1997 a Mw 6.2 earthquake destroyed partially the north-east of Faial and the north-west of Pico islands [2, 3]. Recently, a dense GPS network was installed on Faial and Pico islands with the purpose of monitoring the deformation within these islands. This network was surveyed 4 times (2001, 2003, 2004 and 2006) and has shown that some stations are subsiding by 9 mm/yr. However, the
_____________________________________________________ Proc. ‘Fringe 2009 Workshop’, Frascati, Italy, 30 November – 4 December 2009 (ESA SP-677, March 2010)

Figure 1. Azores Archipelago, Atlantic North. Faial and Pico islands are depicted. Large boxes highlight Synthetic Aperture Radar (SAR) coverage from ENVISAT-ASAR. Gloria fault (GF), East Azores Fracture Zone (EAFZ) and Mid-Atlantic Ridge (MID) are the major geological structures. The strategy used to overcome aforementioned problems is based on the combination of sparse GPS 3D-velocities with two sets of Persistent Scatterers determined from ascending and descending passes. This approach has been already proposed by [4] and [5,6] by the minimization of a energy function that combines GPS and INSAR data. In that solution a continuous deformation map given by INSAR is required for the

Markov Random Field regularization and further simulating annealing optimization. Both solutions requires for its implementation a continuous deformation map derived from INSAR (an interferograms) and a dense grid of GPS stations. In this paper we present the first results obtained by adapting a simple approach to merge GPS and INSAR where we have only few GPS sites (at least 4 sites, covering all the area) and two sets of Persistent Scatterers from ascending and descending passes. 2. METHODOLOGY

= ,

(5)

The main advantage of this decomposition is that the horizontal velocity is inherently smooth and can be interpolated to compute the horizontal velocity at the PS positions. This is generally verified on stable plate tectonics or inter-seismic periods. It is also worth to noting that the horizontal velocity has a small impact on the line-of-sight velocity measure in the area. Inserting the interpolated horizontal velocity in equation (2), the vertical velocity for a given PS is given by:
= 1 ( ? ? )

The outcome of PS-InSAR and GPS processing is the estimation of deformation velocity in a set of sparse terrain points. Generally, the density of PS points is higher than that of GPS stations. In this section a method is described to merge PS-InSAR and GPS measurements to get a map of the vertical component of deformation velocity. The method is suited to be applied to study deformations in regions characterized by stable plate movements or inter-seismic periods. The PS-InSAR estimates the velocity VPS which is the component along the radar line-of-sight (LOS) of the terrain deformation. The LOS is direction between the satellite and the point on the ground. The velocity of a Persistent Scatterer on the line of sight can be related to the velocity on a local reference system by:
= , , . ( , , ) = + +

(6)

The estimated vertical velocities are biased relatively to the local reference system and are affected by long wavelength residual un-modeled errors (unwrapping errors, orbital errors or atmospheric). The GPS observations are not expected to have significant systematic errors and can be used to remove bias and tilts from INSAR data. For that, the INSAR vertical velocities (given by equation (6)) are compared with GPS vertical velocity. Persistent scatterers falling within a circular area around each GPS station are identified and median value of their vertical velocities is computed and compared to the vertical velocity at the corresponding GPS station.

(1) (2)

The constant bias D between INSAR and GPS vertical velocity is computed by minimizing the difference: =
=1 ( )() + ? () 2

(7)

where VPS is the known velocity of a PS on the line-of sight, , , are the unknown components of the velocity vector relative to a local (geodetic) reference system and , , are the components of the unit vector pointing from the ground toward the satellite. The unit vector pointing from the ground to the satellite is given by ([7]):
= ? ? sin ? ? ? cos ? ? ? (3)

where NGPS is the number of GPS stations, ( ) is the median velocity value of PS’ belonging to the 200 m circular area around the i-th GPS station, and D is a 4parameter function representing a bias and a tilt [8].

The estimated bias and tilt are then removed from PS results so obtaining a map of the vertical velocity merging both GPS and ascending and descending INSAR data. 3. GPS MEASUREMENTS

The GPS velocities, known only at a few sparse sites are represented by:
= ( , , )

(4)

Where , are the east, north and vertical velocity components of a GPS site on a regional or global geodetic reference system (p.e. ETRS or ITRS systems). The 3D velocity is usually separated into the horizontal and vertical components. The horizontal velocity component is given by:

On 2001, in the scope of SARAZORES project, a dense geodetic network was installed on the central group of Azores archipelago (Terceira, Faial, Pico, Graciosa and S. Jorge islands) for geodynamic studies. The geodetic network consists on 64 stations (all network) distributed mostly on the borders of the islands. Some of these stations were already installed on 1999 in the aim of TANGO project. There are 20 stations on Faial and 12 in Pico with metallic benchmarks installed in outcrops of solid rock. In the scope of SARAZORES and KARMA projects, four GPS campaigns were conducted

five years apart from 2001 to 2006 [9]. Each site was surveyed for 12 to 24 hours per day over an average of 3 consecutive days. All the surveys were performed using dual-frequency GPS receivers, collecting data every 30 sec. The GPS data analysis has been conducted using the GAMIT and GLOBK software [10] in a three-step approach, described by [11]. For the first stage of the processing, we estimate station coordinates, the zenith delay of the atmosphere at each station, and orbital and Earth orientation parameters (EOPs), using doubly differenced GPS phase observations. All of these parameters are loosely constrained in this step. At this stage in the processing, we included GPS data from a few International GPS Service (IGS) and European Reference Frame (EUREF) sites, in order to serve as ties with the ITRF05 [12]. We used precise satellite orbits computed at the Scripps Orbit and Permanent Array Center at SIO [13] from data collected by the permanent tracking stations of the International GPS Service for Geodynamics [14]. In our GAMIT implementation, we performed bias-fixing and leastsquares adjustments using double-difference combinations of the ionosphere-free linear combination of the L1 and L2 phase observables recorded at each station. In a second step processing, we used the GLOBK software, a Kalman filtering network adjustment, in order to estimate site positions and velocities in the ITRF05 reference frame. At this stage of processing, we combined our GAMIT daily solutions with the loosely constrained solutions (performed by Scripps Orbit and Permanent Array Center) of the global IGS network into a single set for each survey to produce a single, free-network solution for the site coordinates and a single, coherent estimate of the satellite orbits. We then impose the reference frame using this combined solution by minimizing the position and velocity deviations of IGS core stations with respect to the ITRF05 while estimating an orientation, translation, and scale transformation. This estimation is based on 26 globally distributed stations, which we chose as the ones with the lowest uncertainties in position and velocity in the ITRF05 definition. The one sigma uncertainties for the GPS velocities of our sites were derived by scaling the formal errors by the square root of chi-square per degree of freedom of the final adjustment. The chi-square per degree of freedom for this solution was 1.05, indicating a reasonable constraints and a robust estimation of

parameters. About 50% of the sites have horizontal velocity uncertainties lower than 1 mm/yr.

Figure 2. GPS horizontal velocity (a) and GPS vertical velocity (b) for Faial and Pico islands.

4.

INSAR MEASUREMENTS

4.1 Dataset The interferometric dataset consisted of 57 ENVISATASAR images covering Faial, Pico and S. Jorge islands, from 2006 March to 2009 January. Twenty-nine images were acquired along the descending pass and 28 along the ascending one. The perpendicular baseline range from -756m to 944m with a mean value of 126 m. The master image was chosen to minimize both the temporal and perpendicular baseline. 4.2 Atmospheric Artefacts Mitigation Tropospheric delay maps were obtained using a methodology based on the Weather Research and Forecasting (WRF) model of meteorological relevant parameters and a raytracing procedure to estimate tropospheric delay on the SAR grid. Tropospheric delays at the date of SAR acquisitions were computed and used to derive the synthetic atmospheric phase screen of each interferogram of the InSAR time series [15].

Figure 3. Estimated Persistent Scatters for ascending (a) and descending (b) passes. The line-ofsight velocity is shown in mm.yr-1. 4.3 Processing The DORIS software (Delft University of Technology) was used for interferometric processing. The major difficulty in the interferometric processing was the coregistration step due to the minimal amount of land area and low coherence values related to the lush vegetation. A mask was applied to the original data improving considerably the co registration and interferometric processing. On Azores islands the coherence of interferograms is considerably low and the Persistent Scatterers approach seems to be the only reliable technique to extract useful information from interferograms. STAMPS software [16] was used to determine the Permanent Scatterers using the stack of interferograms already processed for the ascending and descending passes. The resulting maps of line-of-sight velocity estimated from PS-INSAR, ascending and descending passes, are displayed in figure 3. The number of PS is very low and irregularly distributed, mostly concentrated on the border of the islands and on few urban areas. STAMPS parameters were fine tuned and a final solution with 5 PS/Km2 on Faial and on Pico were determined.

5

DATA MERGING AND QUALITY ASSESSMENT

GPS and INSAR estimations of deformation are not directly comparable. In fact, while GPS processing gives absolute velocities referred to a regional or global geodetic reference system (e.g. ETRS or ITRS), InSAR technique can just measure relative velocities of PS’ with respect to reference point in the scene. Moreover, ascending and descending passes estimation of deformation are not directly comparable. To deal with this problem, GPS and PS-INSAR estimations of velocity were processed by means of the methodology described in section 2. The horizontal component of the velocity was assumed constant over the Pico and Faial Islands. In fact, GPS measurements are characterized by a small dispersion both in east and north component (1.2 mm/yr and 1.1 mm/yr,

Figure 5. Vertical velocity map resulting from the integration of ascending and descending PS’ and GPS data. The velocity is shown in mm.yr-1. respectively). This GPS estimate of the constant horizontal component of deformation velocity was used to derive the vertical component of velocity at each Persistent Scatterer. After this step, GPS absolute measurements of vertical velocity were used to calibrate InSAR relative measurements. A circular area with a radius of 200 m was studied around each GPS station. Know tectonic phenomena in these islands do not vary within areas having such a size. PS’ falling within each of those areas were identified and median value of their vertical velocities was computed and compared to the vertical velocity at the corresponding GPS station. This operation was repeated for each one of GPS and for both ascending and descending SAR orbits resulting in a number of PS’ varying from 8 to 41 PS. The constant bias D between InSAR and GPS estimation of vertical velocity was computed using equation (7) and removed from PS results. The result is a map of the vertical velocity merging GPS and ascending and descending INSAR data (figure 5). In this figure it is observed the result of the integration of both ascending and descending persistent scatterers with GPS velocity data. The vertical velocity map is reliable and reveals a vast area of vertical deformation on the west and east of Pico island and on the west of Faial close to Capelinhos eruption place. It is also observed a vertical subsidence on the top summit of Pico (with an altitude of 2400 m). Further analysis and interpretation of this vertical deformation map is foreseen. The accuracy of the vertical velocity map was determined by the comparison with GPS stations. For each GPS station a set of PS’ within a circle of 500m was identified. Of each set of PS’s it was computed the inter-quartile distance as a measure of dispersion of the vertical velocity around its median value. PS’ velocity doesn’t follow a normal distribution and may contain gross errors. Because of these, the inter-quartile range is more representative than the standard deviation as an estimate of the spread of data, since changes in the upper and lower 25% of the data do not affect it. In fig. 6, the inter-quartile range and median of the selected set of PS is represented. In this box-plot the median, upper and lower quartiles and minimum and maximum values are represented. For FAIM GPS station, on the SE of Faial island, 41 PS were selected and a interquartile range of 1.3 mm.yr-1 was determined. For all GPS stations the interquartile range of the PS’ vertical velocity is lower than the estimated GPS vertical velocity. The highest interquartlie range is 3.0 and 2.9 mmyr-1 for FCBR (SW of Faial). and FFAJ stations and are a consequence of the reduced number of PS’ close to these stations (only 4 and 6 , respectively).
8

4

0

mm/yr

-4

-8

-12
FAIM PMAD FVUL PRIB FCBR FVUN FAER FVAR PDRO FCDR FFAJ PSMT PCPI PPIL

GPS sites

Figure 6. Box-plot of PS’ vertical velocity selected within a circular area around each GPS station. The same selection set of PS’ associated with a GPS station, was used to compute the dispersion of the the inter-quartile range of the residuals. The residuals were computed as the difference between the median value of the PS’ falling within the circular area and the GPS station velocity. The residuals median and interquartile range were computed for four different data sets,

subdivided into: before merging (standard processing or atmospheric correction) and after merging (standard processing or atmospheric correction). The results are presented in Table 1 for Faial and Pico islands. It is seen that for both islands there are a bias between the PS velocity and the GPS velocity. Also, the interquartile range is higher than expected, particularly in Pico island with a score of 7 mmyr-1. It is important to notice that the vertical velocity is within the millimeter range. The merging strategy have produced a better deformation map removing the bias and reducing the dispersion to about 2 mmyr-1 for both islands. This means that most of PS’ estimated vertical velocities are statistical meaningful. The use of a WRF model to correct for path delay has also a positive impact on the final deformation map. The atmospheric correction improves the result, mostly before fusion. Table 1. Comparison between GPS and InSAR recovered vertical velocity. Values ? = Q(075)-Q(025) Faial Pico in (N=9) (N=6) -1
mm yr

Foundation for Science and technology ( FCT). We acknowledge ESA for providing the ENVISAT-ASAR images to research project n. 3149. 8. REFERENCES

Before Fusion

Standard PS Processing With ATM correction Standard PS Processing With ATM correction

median ? median ? median ? median ?

-11.7 4.6 -11.4 3.7 -0.1 1.9 0.2 1.6

-10.3 7.2 -11.3 5.7 -0.3 2.1 -0.4 1.5

After Fusion

6.

CONCLUSIONS

The vertical deformation map has identified a large continuous area of subsidence on the west of Faial Island, on the flank of Capelinhos eruption cone, with a maximum subsidence range of 10 mm/yr. It has also revealed the subsidence of the summit crater of Pico Island (9 mm/yr) and a large area of subsidence on the west of the island, corresponding mostly to creep movement. It was verified that the integration of GPS and ascending and descending INSAR data reduces the dispersion up to 50% and the correction of original interferogram by WRF-based estimates of tropospheric phase delay results in a further reduction of up to 30% in dispersion. 7. AKNOWLEDGMENTS

This research is part of the project KARMA, Kinematics And Rheological Modelling of the NubianEurasian plate boundary in the Azores, POCI/CTEGIN/57530/2004, supported by the Portuguese

[1] Catal?o, J; Miranda, JM; Lourenco, N. Deformation associated with the Faial (Capelinhos) 1957-1958 eruption: Inferences from 1937-1997 geodetic measurements. Journal of Volcanology and Geothermal Research, Vol 155, Issue: 3-4 Pages: 151-163, 2006 [2] Fernandes, RMS; Miranda, JM; Catal?o, J, et al.. Coseismic displacements of the M-W=6.1, July 9, 1998, Faial earthquake (Azores, North Atlantic). Geophysical Research Letters, Vol29, 16, An1774, 2002. [3] Catita, C; Feigl, KL; Catalao, J, et al.. InSAR time series analysis of the 9 July 1998 Azores earthquake. International Journal of Remote Sensing, Vol26, 13, 2715-2729,2005. [4] Gudmundsson, S; Sigmundsson, F; Carstensen, JM.. Three-dimensional surface motion maps estimated from combined interferometric synthetic aperture radar and GPS data. Journal of Geophysical Research-Solid Earth, 107, B10, Article Number: 2250, 2002. [5] Samsonov, S; Tiampo, K. Analytical optimization of a DInSAR and GPS dataset for derivation of threedimensional surface motion. IEEE Geoscience and Remote Sensing Letters, Vol 3, 1, 107-111, 2006. [6] Samsonov, S; Tiampo, K; Rundle, J, et al.. Application of DInSAR-GPS optimization for derivation of fine-scale surface motion maps of southern California. IEEE Transactions on Geoscience and Remote Sensing, Vol 45, 2, 512-521, 2007. [7] Hanssen, R. Radar Interferometry, data interpretation and error analysis, Kluwer Academic Publishers, 308 p., 2001 [8] Heiskanen, W.H. and H.Moritz. Physical Geodesy. W.H. Freeman and Co., San Francisco, 1967. [9] Catita, C., J. Catal?o e L.M. Victor. Kinematics of Faial-Pico Islands (Azores Archipelago) deduced from repeated GPS. In: Proceedings of VI Assembleia Luso Espanhola de Geodesia e Geofísica, Tomar, 11 a 14 Fevereiro de 2008. CD. Cambridge Univ Press ISB 848320-373-1, pp. 195-196, 2008. [10] King, R. and Y. Bock. Documentation for the GAMIT GPS software analysis, release 10.05. Mass. Inst. of Technol. Cambridge, Atmospheric, and Planetary Science. University of California at San Diego, 2001 [11] Feigl, K., D. Agnew, Y. Bock, D. Dong, A. Donnellan, B. Hager, T. Herring, D. Larsen, K. Larson, M. Murray, Z. Shen, and F. Webb (1993). Space geodetic measurement of the velocity field of central and southern California, 1984-1992, JGR, v. 98 (B12), 21677-21712. [12] Altamimi, Z., P. Sillard, and C. Boucher (2002). ITRF2000: A new release of the International Terrestrial Reference Frame for earth science

applications, J. Geophys. Res., 107(B10), 2214, doi:10.1029/2001JB000561 [13] Fang, P. and Y. Bock, "Scripps Orbit and Permanent Array Center Report to the IGS - 1995," 1994 Annual Report, International GPS Service for Geodynamics, J. F. Zumberge, R. Liu and R. E. Neilan, eds., IGS Central Bureau, Jet Propulsion Laboratory, Pasadena, 1995, (pp. 213-233). [14] Beutler, G., I.I. Mueller, and R.E. Neilan (1994). The International GPS Service for Geodynamics (IGS):

Development and start of official service on January 1, 1994, Bull. Geod., 68, 39-70, 1994. [15] Nico, G., R. Tomé, P. Benevides, J. Catal?o, and P. Miranda, 2009.Interferometric SAR analysis of atmospheric water vapor properties. Geophysical Research Abstracts, Vol. 11, EGU2009-7904, 2009. [16] Hooper, A; Segall, P; Zebker, H. Persistent scatterer interferometric synthetic aperture radar for crustal deformation analysis, with application to Volcan Alcedo, Galapagos. Journal Of Geophysical ResearchSolid Earth, 112, B7, B07407, JUL 10 2007.


相关文章:
P2 项目管理 乌龟图.pdf
P2 项目管理 乌龟图 - 项目管理乌龟图 项目起点 ※电脑及办公软件、投影仪;
室内P2参数.doc
室内P2参数 - 室内 P2 全彩 LED 技术参数 规格型号 像素间距(mm)
过程审核P2-P4.xls
过程审核P2-P4 - VDA6.3过 程审核提问表(A) 序号 星号 VDA6. 条款 3条款 P2 提问项目 项目管理 是否建立了项目组织机构(项目 管理),并且为项目管理负责人 ...
一张图读懂P2P的五大衍生模式(P2N、P2B、P2C、P2G都是啥?.doc
一张图读懂P2P的五大衍生模式(P2N、P2B、P2C、P2G都是啥? - 导语:你了解 P2P 的新模式吗?如果不知道,那么就看易 P2P 的分享 P2N、P2B、P2C、P2G 模式。...
ACCA P2 双语讲义.doc
ACCA P2 双语讲义 - Course Notes ACCA P2 Corp
51单片机的P1、P2、P3口的工作原理.doc
51单片机的P1、P2、P3口的工作原理 - 基于 proteus 的51单片机仿真实例五十、51单片机的 P1、P2、P 3口的工作原理 1、P1口某一位的内部电路结构如下图所示,...
P2实验室二级生物安全实验室主要技术指标的确定.doc
P2实验室二级生物安全实验室主要技术指标的确定 - P2 实验室即二级生物安全实
P2级实验室一般规范.ppt
P2级实验室一般规范 - P2 级实验室一般规范 z4c8e 实验室通风柜 ww
CATIA.P1,p2,p3等区别.doc
CATIA.P1,p2,p3等区别 - 由于P1、p2和p3应用平台都是在相同的
浅析P2摄像机的使用_图文.pdf
浅析P2摄像机的使用 - !‘ 钱 ?# ? ? 中央 电视 台 本文 介绍 %
如图,设P1,P2,…,P6为单位圆上逆时针均匀分布的六个点.....doc
如图,设P1,P2,…,P6为单位圆上逆时针均匀分布的六个点.现从这六个点中任选
P1P2P3储量定义.doc
P1P2P3储量定义 - 3.1 证实的 P1 储量的定义 3.1.1 证实储量
P2技术规格书.doc
P2技术规格书 - 高清P2显示屏、华星 高科 LEDTV... P2技术规格书_电子/电路_工程科技_专业资料。高清P2显示屏、华星 高科 LEDTV P2技术参数货物名称 室内全彩色...
P2背膘测定仪操作程序.doc
P2背膘测定仪操作程序 - P2 背膘测定仪操作程序 背膘测试仪用脉冲超声波来测
用松下P2 HD完成低成本电影.doc
用松下P2 HD完成低成本电影 - 用松下 P2 HD 完成低成本电影 谁是卡门
【ACCA考经分享】 一次性过P1 P2 P3的经验.pdf
泽稷网校-财务金融证书在线教育领导品牌 【ACCA 考经分享】不要轻易玩"3P",一次性过 P1P2P3 的 学员分享本人是江财 ACCA 专业的,12 月考试之前每天都昏天黑...
2018年9-12月p2题库(1)_图文.doc
2018年9-12月p2题库(1) - 1.成功经验 Describe a su
混合动力汽车P2系统介绍_图文.doc
混合动力汽车P2系统介绍 - 1. 首先解释一下 P2 的概念:P2 是一种并联
ACCA P2怎么学?需要准备多长时间.pdf
ACCA P2怎么学?需要准备多长时间_财会/金融考试_资格考试/认证_教育专区。ACCA考经分享 泽稷网校-财务金融证书在线教育领导品牌 ACCA P2 怎么学?需要准备多长时间 ...
美国注册管理会计师-P2-B.doc
美国注册管理会计师-P2-B - 题目 1: 2B2-LS41 债券受托人的法律