Re: Problem with plotting polygons using PSXY and also problem with GMTSPATIAL
Joaquim Luis <[email protected]> Fri, 1 Aug 2014 22:37:03 -0700
| Newsgroups | gmane.comp.gis.gmt.user |
|---|---|
| Message-ID | <[email protected]> |
> I can also confirm all of Kara=E2=80=99s observations. If this is a proble= m that > is unlikely to be fixed in the near term, is it relatively straightforward= > to extract GMT 4.5.8 (or similar) from the GMT SVN? Otherwise I=E2=80=99ll= have to > try dig up an old compiled folder=E2=80=A6 >=20 This is a known issue. For those of you that need it badly, the simplest sol= ution is to have installed both GMT4 and GMT5. See the cookbook for how to h= ave both running side by side.=20 Joaquim > Thanks for the help and time in advance! >=20 > Cheers, > Sabin >=20 > PS: I attach a sample of my own problem, the grey polygons (LIPs) that > cross the east map boundary cause the artefact. However, as Kara > demonstrated, often the west boundary is also affected. >=20 >=20 > -- >=20 > MR SABIN ZAHIROVIC | PhD Candidate > School of Geosciences | Faculty of Science >=20 >=20 >=20 >=20 > THE UNIVERSITY OF SYDNEY > Rm 414, Madsen Building F09 | The University of Sydney | NSW | 2006 > M +61 416 775 589 =20 > E [email protected] | W http://www.earthbyte.org > <http://www.earthbyte.org/> >=20 > CRICOS 00026A > This email plus any attachments to it are confidential. Any unauthorised > use is strictly prohibited. If you receive this email in error, please > delete it and any attachments. >=20 >=20 >=20 >=20 >> On 1/08/2014 2:52 pm, "Kara Matthews" <[email protected]> wrote= : >>=20 >> Hi Christian, and hi everyone, >>=20 >> Thanks very much for your reply. The script I previously sent also >> included my work-around in it so I could show what steps I had to take to= >> overcome the problem. Sorry, I should have said which lines needed >> commenting out. Sabin downloaded the developer version 5.1.2 (r13381, 64 >> bit) and it didn=C2=B9t work for him. >>=20 >> I did some more testing and this is what I have concluded: >> Using v5.1.1 (r12968, 64 bit) and 5.1.2 (r13381, 64 bit) I cannot use psx= y >> to make a mercator map that contains 180deg longitude and contains >> polygons that are truncated by both the left and right map margins - I ge= t >> banding across the map. Depending on whether you specify your region in >> the form 0-360 or =C2=AD180/=C2=AD180 you can get a map with polygons tha= t are >> truncated by one of the margins without banding, but not both. >>=20 >> In the attached scripts I plot 8 rectangles (10 degrees wide) in each map= . >> 4 are cut in half by the right margin and 4 are cut in half by the left >> margin. The 4 rectangles at each margin differ as follows: >> - 1 is digitised clockwise and coordinates are =C2=AD180/180 >> - 1 is digitised counterclockwise and coordinates are =C2=AD180/180 >> - 1 is digitised clockwise and coordinates are 0/360 >> - 1 is digitised counterclockwise and coordinates are 0/360 >>=20 >> I ran tests in 3 locations (three directories within >> polygon_fill_problem.zip): eastern US, Pacific, Australia >>=20 >> For each location, the top two maps have the same =C2=ADR region but are >> specified either in the form =C2=AD180/180 or 0/360, and the bottom two m= aps >> have the same location but the =C2=ADR=C2=B9s are specified differently. I= have >> annotated this. >>=20 >> You will see that as long as the map doesn=C2=B9t include 180 degrees >> longitude, see attached gmt_poly_plotting_test_E-US.pdf (from test1_E-US >> example), then it is possible to make the map without banding, you just >> have to choose the right =C2=ADR format. However, if the map includes 180= >> degrees longitude, see gmt_poly_plotting_test_AUS_2.pdf and >> gmt_poly_plotting_test_PAC_2.pdf (from test2_PAC top two maps, and >> test3_AUS bottom two maps), then it is not possible to have polygons >> truncated on each side of the map without banding. >>=20 >> In the attached directories contained within polygon_fill_problem.zip all= >> the files are in place and you just have to run the scripts as they are. I= >> am using 5.1.1 (r12972) 64-bit. >>=20 >> Interestingly, I did some further testing and this problem does not occur= >> using v4.5.8 (64-bit)!! >>=20 >> I=C2=B9m sorry to be a bother, but I had a bunch of plotting scripts that= I >> made late last year using an older developer version of GMT5 (I can=C2=B9= t >> remember which version sorry =C2=AD but I can find out tonight if needed)= that >> work and now they don=C2=B9t. Plus I know a couple of others in the group= >> having this problem. >>=20 >>=20 >> Does anyone have any recommendations? To the developers, is there a way t= o >> change the GMT thing responsible for this back to the way it was in >> v4.5.8? >>=20 >> Thank you very very much for your help Christian and to anyone else who >> can offer any suggestions. >>=20 >> Cheers, >> Kara >>=20 >>=20 >>=20 >>=20 >>=20 >>=20 >>=20 >>> On 24/07/2014 7:40 am, "Christian Heine" <[email protected]> wrot= e: >>>=20 >>> Hi Kara, >>>=20 >>>> I am using GMT 5.1.1 (r12972) on OSX Mavericks and am having a problem >>>> with plotting polygons (I=C2=B9m using the Mercator projection). Specif= ically >>>> I am getting banding across my map in instances when a polygon is being= >>>> cropped by the map border. It seems that the polygon is wrapping the >>>> wrong way around the world. Compare =C2=B3test_0-360.pdf=C2=B2 (artefac= t in bottom >>>> right) with =C2=B3test_required_plotting.pdf=C2=B2 (this is what I want= the map to >>>> look like). A couple of my other colleagues at work are having the same= >>>> problem. >>>=20 >>> I just ran your test.sh script using Version 5.1.2_r13178 [64-bit] on OS= >>> 10.9.4 and get the attached result (test_0-360_chhei_crop.pdf). It is >>> quite similar to the look of your 'desired' map. The only difference >>> between your desired map and the one I generated using your script is th= e >>> central part where there are white/unfilled patches (missing polygons?)= >>> as opposed to filled space in your required map. >>>=20 >>> Plotting COBreconstructed_${age}.00Ma.xy (2.pdf) and >>> WARSreconstructed_${age}.00Ma.xy (1.pdf) independently results in the >>> attached maps.=20 >>>=20 >>> For what it's worth it might be helpful to dump the reconstructed data a= s >>> *.shp file and throw this into QGIS and do basic geometry checking: In >>> QGIS there's the menu Vector -> Check Geometry Validity (btw this should= >>> also work on OGR GMT files). But I suspect that this is more something >>> related to dateline crossings. >>>=20 >>> Also, I ran your polygons through gmtspatial -Q+ for a handedness check -= >>> some few of your polygons are CW whereas the others are CCW handed. Not >>> sure whether this causes GMT to misbehave (wouldn't assume so). >>>=20 >>> When exporting from GPlates, I can highly recommend to use the GMT OGR >>> format (*.gmt) as this also exports the GPlates Feature ID which lets yo= u >>> search and identify a single rogue feature much much easier. >>>=20 >>> Possible workaround: ogr2ogr from the GDAL library has a "[-clipsrc [xmi= n >>> ymin xmax ymax]" option which you can also use to clip the data. >>>=20 >>> Sorry for not being able to help with specific debugging info with your >>> setup. Maybe an update to the recent SVN versions will fix some of the >>> problems?=20 >>>=20 >>> Greetings from NL, >>> Christian >>>=20 >>>=20 >>>=20 >>> Mailing list for GMT discussions of all kinds. If you are not sure you >>> have found a bug, discuss it here first. >>> To formally report bugs or request features, please register and add New= >>> Issue on gmt.soest.hawaii.edu >>> To unsubscribe, send the message "signoff gmt-help" to >>> [email protected] >>> Note: gmt-help will become obsolete on Sept 1, 2014 - please use forum o= n >>> gmt.soest.hawaii.edu instead. >>>=20 >>>=20 >>>>=20 >>>> To try and get around this issue I used gmtspatial and the =C2=ADC opti= on to >>>> close the polygons at the map border using my -R. This sometimes works,= >>>> but not all the time. In some cases it seems to be inventing some >>>> strange polygons that I can=C2=B9t account for. See >>>> "test_0-360_gmtspatial.pdf" >>>>=20 >>>> So when gmtspatial wasn=C2=B9t working I used gmtconvert to separate my= >>>> polygons out into separate files to then isolate the problem polygons. >>>> It seems that gmtspatial was making up some polygons with strange >>>> coordinates, coordinates that weren=C2=B9t in the original polygon file= . See >>>> prob_poly.txt that coincidentally has the same longitudes as my map >>>> region. >>>>=20 >>>> I was able to get around my problem by working out which polygons >>>> weren=C2=B9t working (a bit of trial and error with awk and bash if >>>> statements) and then not plotting these=C5=A0 but my clunky work-around= is >>>> not ideal for making lots of figures, or plotting lots of polygon files= . >>>>=20 >>>> I have also attached a test script and the input files in >>>> problem_with_psxy.zip. I have currently uncommented the set of lines I >>>> needed to make my ideal figure. You should be able to run this script >>>> from the folder if needed. >>>>=20 >>>> Thanks kindly for any help you can offer!! >>>>=20 >>>> Kind regards, >>>> Kara >>>>=20 >>>>=20 >>>> --- >>>>=20 >>>> Dr Kara J. Matthews | Postdoctoral Researcher >>>> School of Geosciences | Faculty of Science >>>>=20 >>>> THE UNIVERSITY OF SYDNEY >>>> Rm 403 | Madsen Building (F09) | The University of Sydney | NSW | 2006 >>>> E: [email protected] >>>> T: +61 2 9351 3625 >>>> W: www.earthbyte.org >>>> Mailing list for GMT discussions of all kinds. If you are not sure you >>>> have found a bug, discuss it here first. To formally report bugs or >>>> request features, please register and add New Issue on >>>> gmt.soest.hawaii.edu To unsubscribe, send the message "signoff gmt-help= " >>>> to [email protected] Note: gmt-help will become obsolete on Sep= t >>>> 1, 2014 - please use forum on gmt.soest.hawaii.edu instead. >>>> <test_0-360.pdf><test_0-360_gmtspatial.pdf><test_required_plotting.pdf>= < >>>> p >>>> rob_poly.txt><problem_with_psxy.zip> >>>=20 >>>=20 >>> Mailing list for GMT discussions of all kinds. If you are not sure you >>> have found a bug, discuss it here first. >>> To formally report bugs or request features, please register and add New= >>> Issue on gmt.soest.hawaii.edu >>> To unsubscribe, send the message "signoff gmt-help" to >>> [email protected] >>> Note: gmt-help will become obsolete on Sept 1, 2014 - please use forum o= n >>> gmt.soest.hawaii.edu instead. >>=20 >>=20 >> Mailing list for GMT discussions of all kinds. If you are not sure you >> have found a bug, discuss it here first. >> To formally report bugs or request features, please register and add New >> Issue on gmt.soest.hawaii.edu >> To unsubscribe, send the message "signoff gmt-help" to >> [email protected] >> Note: gmt-help will become obsolete on Sept 1, 2014 - please use forum on= >> gmt.soest.hawaii.edu instead. >=20 >=20 > Mailing list for GMT discussions of all kinds. If you are not sure you ha= ve found a bug, discuss it here first. > To formally report bugs or request features, please register and add New I= ssue on gmt.soest.hawaii.edu > To unsubscribe, send the message "signoff gmt-help" to [email protected]= ii.edu > Note: gmt-help will become obsolete on Sept 1, 2014 - please use forum on g= mt.soest.hawaii.edu instead. > <new_gmt_0025.png> Mailing list for GMT discussions of all kinds. If you are not sure you have found a bug, discuss it here first. To formally report bugs or request features, please register and add New Issue on gmt.soest.hawaii.edu To unsubscribe, send the message "signoff gmt-help" to [email protected] Note: gmt-help will become obsolete on Sept 1, 2014 - please use forum on gmt.soest.hawaii.edu instead.