Re: GRAVFFT query/check

Hi Joaquim,

Finally got all working as expected with gravfft (v5)

gmtset PROJ_ELLIPSOID Sphere
---- Create land and sea masks ----
grdlandmask -R20/50/35/50 -Di -A10000/0/0 -I1m -N0/1 -GBLACK1_land_mask.grd -V -r
grdlandmask -R20/50/35/50 -Di -A10000/0/0 -I1m -N1/0 -GBLACK1_sea_mask.grd -V -r
--- Compute the Bouguer correction for land (2670) and sea (2670-1030 - 1640)
gravfft BLACK1-topo.grd -GBLACK1_corr_land.grd -D2670 -E3 -V -N3600/1800
gravfft BLACK1-topo.grd -GBLACK1_corr_sea.grd -D1640 -E3 -V -N3600/1800
---- Mask the land and sea Bouguer correction grids -----
grdmath BLACK1_corr_land.grd BLACK1_land_mask.grd MUL -V = BLACK1-BAA_corr_land-clip.grd
grdmath BLACK1_corr_sea.grd BLACK1_sea_mask.grd MUL -V = BLACK1-BAA_corr_sea-clip.grd
grdmath BLACK1-BAA_corr_land-clip.grd BLACK1-BAA_corr_sea-clip.grd ADD -V = BLACK1-Land-Sea_Corr.grd
---- Upward continue the full correction to smooth the high-frequency ----
grdfft BLACK1-Land-Sea_Corr.grd -GBLACK1-Land-Sea_Corr_UC4km.grd -fg -C4000 -V
---- Bouguer anomaly = Free-air - Bouguer correction
grdmath BLACK1-FAA.grd BLACK1-Land-Sea_Corr_UC4km.grd ADD -V = BLACK1-Bouguer_anomaly.grd
---- Some post-filtering (if noisy data) ----
grdfilter BLACK1-Bouguer_anomaly.grd -GBLACK1-Bouguer_anomaly_LP20km.grd -Fc25 -D1 -V

The results are exactly the same as working out the corrections via grdmath. Yes the syntax is slightly different in GMT 5 but the basic elements stay the same; not a big challenge to learn.

Did notice that if one sets the Parker expansion order too high (eg E6) then the correction has large errors. Up to E3 seems to be ok.

Thanks
Lester

Date: Sun, 2 Feb 2014 15:43:39 +0000
From: [email protected]
Subject: Re: [GMT-HELP] GRAVFFT query/check
To: [email protected]


  
    
  
  
     Lester,

    

    There is no difference between GMT4 & GMT5 procedures (except
    that some option syntax did actually change since there never was a
    gravfft GMT4)

    And yes, you need to blank (set to zero) the sea-side of the land
    grid and vice-versa using a mask and add the results at the end.

    

    Joaquim

    

    
      
      Hello

        

        I just wanted to check the use of gravfft for the Bouguer
        correction (land and Sea). Windows 7 operation.

        

        Under v4 I would split the computation for land/sea and so this
        would become (under v5):

        

        gravfft test_topo.grd -Gtopo_corr_land.grd -D2670 -E6 -V (use
        land mask -N0/1) 

        gravfft test_bathy.grd -Gtopo_corr_sea.grd -D1640 -E6 -V (use
        sea mask -N1/0)

        

        grdmath free-air.grd topo_corr_land.grd SUB -V =
        bouguer_land.grd

        grdmath free-air.grd topo_corr_sea.grd SUB -V = bouguer_sea.grd

        grdmath bouguer_land.grd bouguer_sea.grd ADD -V =
        Bouguer_land_sea.grd

        

        What I noticed was that, using the land-only grid, the values in
        the "sea" were say around a near constant -53, when the land
        mask has the sea set to 0. This would clearly give a wrong value
        when added to the sea-only data. Clearly, one cannot have NaN
        values in the grid as the FFT would fail, in fact when tested it
        caused gravfft to crash.

        

        If it is all land or all sea, there is no problem, but how can
        this be done for a combined region properly? Do I compute the
        land and sea density grids over the whole area, and then mask
        the results - just trying to see the logical path here.

        

        Under GMT 4 I would simply use grdmath and the mask grids to get
        the separate files to sum, no issues there. 

        

        I am sure it is a simple answer (hopefully)

        

        Cheers

        Lester

        

        

        

        

      
      To unsubscribe, send the message "signoff gmt-help" to
      [email protected]
    
    

  

To unsubscribe, send the message "signoff gmt-help" to [email protected] 		 	   		  
To unsubscribe, send the message "signoff gmt-help" to [email protected]