ezprojections.pas

很管用的GIS控件

Unit EzProjections;

{***********************************************************}
{     EzGIS/CAD Components                                  }
{   (c) 2003 EzSoft Engineering                             }
{         All Rights Reserved                               }
{***********************************************************}

{$I EZ_FLAG.PAS}
Interface

Uses
  SysUtils, Classes, Windows, ezbase, ezlib ;

Const
  EN_SIZE = 5;
  APA_SIZE = 3;
  HUGE_VAL = 1.0E+300;
  HALFPI = PI / 2;
  FORTPI = PI / 4;
  TWOPI = 2 * PI;

  IS_ANAL_XL_YL = 01; //* derivatives of lon analytic */
  IS_ANAL_XP_YP = 02; //* derivatives of lat analytic */
  IS_ANAL_HK = 04; //* h and k analytic */
  IS_ANAL_CONV = 010; //* convergence analytic */

Type

  // commented projections will be implemented in near future
  TEzProjectionCode = (
    aea,
    //aeqd     ,
    airy,
    //aitoff   ,
    //alsk     ,
    //apian    ,
    //august   ,
    //bacon    ,
    //bipc     ,
    //boggs    ,
    bonne,
    cass,
    cc,
    cea,
    //chamb    ,
    //collg    ,
    //crast    ,
    //denoy    ,
    //eck1     ,
    //eck2     ,
    //eck3     ,
    eck4,
    eck5,
    eck6,
    //eqc      ,
    //eqdc     ,
    //euler    ,
    //fahey    ,
    //fouc     ,
    //fouc_s   ,
    //gall     ,
    //gins8    ,
    gn_sinu,
    //gnom     ,
    //goode    ,
    //gs48     ,
    //gs50     ,
    //hammer   ,
    //hatano   ,
    imw_p,
    //kav5     ,
    //kav7     ,
    //labrd    ,
    laea,
    //lagrng   ,
    //larr     ,
    //lask     ,
    lcc,
    leac,
    //lee_os   ,
    //loxim    ,
    //lsat     ,
    //mbt_s    ,
    //mbt_fps  ,
    //mbtfpp   ,
    //mbtfpq   ,
    mbtfps,
    merc,
    //mil_os   ,
    mill,
    //mpoly    ,
    moll,
    //murd1    ,
    //murd2    ,
    //murd3    ,
    //nell     ,
    //nell_h   ,
    //nicol    ,
    //nsper    ,
    //nzmg     ,
    //ob_tran  ,
    //ocea     ,
    //oea      ,
    omerc,
    //ortel    ,
    ortho,
    //pconic   ,
    poly,
    //putp1    ,
    //putp2    ,
    //putp3    ,
    //putp3p   ,
    //putp4p   ,
    //putp5    ,
    //putp5p   ,
    //putp6    ,
    //putp6p   ,
    //qua_aut  ,
    //robin    ,
    //rpoly    ,
    sinu,
    //somerc   ,
    stere,
    //tcc      ,
    tcea,
    //tissot   ,
    tmerc,
    tpeqd,
    //tpers    ,
    ups,
    //urm5     ,
    //urmfps   ,
    utm,
    //vandg    ,
    //vandg2   ,
    //vandg3   ,
    vandg4,
    //vitk1    ,
    //wag1     ,
    //wag2     ,
    //wag3     ,
    wag4,
    wag5,
    //wag6     ,
    wag7
    //weren    ,
    //wink1    ,
    //wink2    ,
    //wintri
    );

  TEzEllipsoidCode = (
    ecMERIT,
    ecSGS85,
    ecGRS80,
    ecIAU76,
    ecairy,
    ecAPL4_9,
    ecNWL9D,
    ecmod_airy,
    ecandrae,
    ecaust_SA,
    ecGRS67,
    ecbessel,
    ecbess_nam,
    ecclrk66,
    ecclrk80,
    ecCPM,
    ecdelmbr,
    ecengelis,
    ecevrst30,
    ecevrst48,
    ecevrst56,
    ecevrst69,
    ecevrstSS,
    ecfschr60,
    ecfschr60m,
    ecfschr68,
    echelmert,
    echough,
    ecintl,
    eckrass,
    eckaula,
    eclerch,
    ecmprts,
    ecnew_intl,
    ecplessis,
    ecSEasia,
    ecwalbeck,
    ecWGS60,
    ecWGS66,
    ecWGS72,
    ecWGS84,
    ecITRFMEX
    );

  { future use
  TUV = record
     u, v: Double;
  end;

  TCOMPLEX = record
     r, i: Double;
  end; }

  TXY = Record
    x, y: double;
  End;

  TLP = Record
    lam, phi: double;
  End;

  TPJ_ELLPS = Record
    //id: String;
    major: String;
    ell: String;
    name: String;
  End;

  TDERIVS = Record
    x_l, x_p: double;
    y_l, y_p: double;
  End;

  TFACTORS = Record
    der: TDERIVS;
    h, k: double;
    omega, thetap: double;
    conv: double;
    s: double;
    //a,
    b: double;
    code: integer;
  End;

  TEzGeoConvert = Class;

  TForward = Function( lp: TLP; P: TEzGeoConvert ): TXY;
  TInverse = Function( xy: TXY; P: TEzGeoConvert ): TLP;
  TProj = Function( P: TEzGeoConvert; init: boolean ): TEzGeoConvert;
  TSpecial = Procedure( lp: TLP; P: TEzGeoConvert; Var fac: TFACTORS );

  { the list of projections }
  TPJ_LIST = Record
    //ID: String;
    Descr: String;
    Proj: TProj; //projection entry point
  End;

  { projection parameters handling }
  TEzProjectParam = Class
  Private
    fGeoConvert: TEzGeoConvert;
  Public
    Constructor Create( GeoConvert: TEzGeoConvert );
    Function Defined( Const opt: String ): boolean;
    Function AsString( Const opt: String ): String;
    Function AsInteger( Const opt: String ): Integer;
    Function AsFloat( Const opt: String ): double;
    Function AsRadians( Const opt: String ): double;
    Function AsBoolean( Const opt: String ): boolean;
  End;

  //--------- the low level conversion object -------------------------//
  TEzGeoConvert = Class
  Private
    fParaList: TStringList;
    fpj_param: TEzProjectParam;
    Function pj_inv( Var xy: TXY ): TLP;
    Function pj_fwd( Var lp: TLP ): TXY;
    Function pj_ell_set( Var a, es: double ): integer;
  Public
    pj_errno: integer;

    { projection variables }
    fwd: TFarProc; // procedure pointer for converting from lat/long to projection
    inv: TFarProc; // procedure pointer for converting from projection to lat/lon
    spc: TFarProc; // procedure pointer for special processing

    over: boolean;
    geoc: boolean;
    a: double; { major axis or radius if es=0 }
    e: double; { eccentricity }
    es: double; { e ^ 2 }
    ra: double; { 1/A }
    one_es: double; { 1 - e^2 }
    rone_es: double; { 1/one_es }
    lam0, phi0: double; { central longitude, latitude }
    x0, y0: double; { false easting and northing }
    k0: double; { general scaling factor }
    to_meter, fr_meter: double; { cartesian scaling }
    units: String[8]; { for fast access}

    { additional info for every projection }
    ec: double;
    n: double;
    c: double;
    dd: double;
    n2: double;
    rho0: double;
    rho: double;
    phi1: double;
    phi2: double;
    ellips: boolean; // ellipsoid or spheroid

    sinph0: double;
    cosph0: double;
    M1: double;
    N1: double;
    Mp: double;
    He: double;
    G: double;
    mode: Integer;

    p_halfpi: double;
    Cb: double;
    no_cut: boolean; { do not cut at hemisphere limit }

    cosphi1: double;

    bacn: Integer;
    ortl: Integer;

    noskew: Integer;

    cphi1: double;
    am1: double;

    m0: double;
    t: double;
    a1: double;
    r: double;
    d2: double;
    a2: double;
    tn: double;

    ap: double;

    esp: double;
    ml0: double;
    en: Array[0..EN_SIZE] Of double;
    apa: Array[0..APA_SIZE] Of double;
    { oblique mercator data }
    bl: double;
    singam: double;
    al: double;
    el: double;
    cosgam: double;
    u_0: double;
    rot: boolean;
    cosrot: double;
    sinrot: double;
    Gamma: double;
    lamc: double;
    alpha: double;
    lam1: double;
    lam2: double;
    { mollweide}
    C_x: double;
    C_y: double;
    C_p: double;
    sinb1: double;
    cosb1: double;
    xmf: double;
    ymf: double;
    mmf: double;
    rq: double;
    phits: double;
    sinX1: double;
    cosX1: double;
    akm1: double;
    P, Pp, Q, Qp,
      R_1, R_2,
      sphi_1, sphi_2, C2: double;
    phi_1, phi_2,
      lam_1, m,
      sp1, cp1, sp2,
      dlam2, cp2,
      z02, r2z0, ccs, cs, sc, hz0,
      thz0, rhshz0, lp, sa, ca, rk0: double;

    { methods }
    Constructor Create;
    Destructor Destroy; Override;

    Procedure Geo_CoordSysInit( Params: TStringList );
    Procedure Geo_CoordSysFromLatLong( Const Long, Lat: Double; Var x, y: Double );
    Procedure Geo_CoordSysToLatLong( Const x, y: Double; Var Long, Lat: Double );
    Function Geo_Distance( Const Long1, Lat1, Long2, Lat2: Double ): Double;

    Property ParaList: TStringList Read fParaList;
    Property pj_param: TEzProjectParam Read fpj_param;
  End;

  {-----------------------------------------------------------------------------}
  {    TEzProjector - the component used to calculate for projections           }
  {-----------------------------------------------------------------------------}