fix a bug in coupe ( -A+wwidth is useless when using -N )

[CovMat]

gmt/src/seis/pscoupe.c

Lines 1331 to 1333 in 58850df

	if (!pscoupe_dans_coupe (in[GMT_X], in[GMT_Y], depth, Ctrl->A.xlonref, Ctrl->A.ylatref, Ctrl->A.fuseau, Ctrl->A.PREF.str,
	Ctrl->A.PREF.dip, Ctrl->A.p_length, Ctrl->A.p_width, &distance, &n_dep) && !Ctrl->N.active)
	continue;

if -N is used, !Ctrl->N.active is always 0, so the whole condition is always 0 too.

However, if user wants to use -A+wwidth to choose data points only in the range of the

width in km of the cross-section on each side of a vertical plane or above and under an oblique plane
https://docs.generic-mapping-tools.org/6.6/supplements/seis/coupe.html#a

, although the data points outside the range let pscoupe_dans_coupe returns 0 and !pscoupe_dans_coupe equal to 1:

gmt/src/seis/pscoupe.c

Lines 353 to 371 in 58850df

	GMT_LOCAL int pscoupe_dans_coupe (double lon, double lat, double depth, double xlonref, double ylatref, int fuseau, double str, double dip, double p_length, double p_width, double *distance, double *n_dep) {
	/* if fuseau < 0, cartesian coordinates */
	double xlon, ylat, largeur, sd, cd, ss, cs;
	if (fuseau >= 0)
	pscoupe_gutm (lon, lat, &xlon, &ylat, fuseau);
	else {
	xlon = lon;
	ylat = lat;
	}
	sincosd (dip, &sd, &cd);
	sincosd (str, &ss, &cs);
	largeur = (xlon - xlonref) * cs - (ylat - ylatref) * ss;
	*n_dep = depth * sd + largeur * cosd (dip);
	largeur = depth * cosd (dip) - largeur * sd;
	*distance = (ylat - ylatref) * cs + (xlon - xlonref) * ss;
	return (*distance >= 0. && *distance <= p_length && fabs (largeur) <= p_width);
	}

these data points will not be skipped because !Ctrl->N.active is always 0.

Actually, in docs: https://docs.generic-mapping-tools.org/6.6/supplements/seis/coupe.html#n , -N means

Does not skip symbols that fall outside map border [Default plots points inside border only].

it is well done in lines:

gmt/src/seis/pscoupe.c

Line 1259 in 58850df

gmt_set_basemap_orders (GMT, Ctrl->N.active ? GMT_BASEMAP_FRAME_BEFORE : GMT_BASEMAP_FRAME_AFTER, GMT_BASEMAP_GRID_BEFORE, GMT_BASEMAP_ANNOT_BEFORE);

gmt/src/seis/pscoupe.c

Line 1265 in 58850df

if (!Ctrl->N.active) gmt_map_clip_on (GMT, GMT->session.no_rgb, 3);

gmt/src/seis/pscoupe.c

Lines 1338 to 1341 in 58850df

	if (!Ctrl->N.active) {
	gmt_map_outside (GMT, xy[GMT_X], xy[GMT_Y]);
	if (abs (GMT->current.map.this_x_status) > 1 || abs (GMT->current.map.this_y_status) > 1) continue;
	}

gmt/src/seis/pscoupe.c

Line 1684 in 58850df

if (!Ctrl->N.active) gmt_map_clip_off (GMT);

It makes no sense to use !Ctrl->N.active in

if (!pscoupe_dans_coupe (in[GMT_X], in[GMT_Y], depth, Ctrl->A.xlonref, Ctrl->A.ylatref, Ctrl->A.fuseau, Ctrl->A.PREF.str,
     Ctrl->A.PREF.dip, Ctrl->A.p_length, Ctrl->A.p_width, &distance, &n_dep) && !Ctrl->N.active)
     continue;

[CovMat]

If we choose to keep && !Ctrl->N.active in

gmt/src/seis/pscoupe.c

Lines 1331 to 1333 in 58850df

	if (!pscoupe_dans_coupe (in[GMT_X], in[GMT_Y], depth, Ctrl->A.xlonref, Ctrl->A.ylatref, Ctrl->A.fuseau, Ctrl->A.PREF.str,
	Ctrl->A.PREF.dip, Ctrl->A.p_length, Ctrl->A.p_width, &distance, &n_dep) && !Ctrl->N.active)
	continue;

The description of -N in https://docs.generic-mapping-tools.org/6.6/supplements/seis/coupe.html#n should be modified as:

Ignore the limitations of map border and the cross-section, and plot all points [Default plots points inside border and inside the cross-section]. Note that the limitation of cross-section is 3D when using -A+wwidth

[CovMat]

Here is a example explaining this bug

data: nichong.txt

I want to plot data points inside the width=4km on each side of a vertical plane (+w4). The range of map in depth direction is 5 km ~ 20 km, but there are some data points above 5km depth. I want to plot all +w4 data points, including those above depth 5km. I use following script:

gmt begin figure2 png
gmt basemap -R0/30/5/20 -JX15c/-5c -BWSne -Byaf+l"Focal depth (km)" -Bxaf+l"Distance (km)"
gmt coupe nichong.txt -Sa0.9c -Aa78.43/41.26/78.70/41.05+d90+w4 -W0.4p -GLIGHTSEAGREEN -Q -N
gmt end

and get following figure:

This figure is wrong, because the data points outside the width=4km are also plotted. -A+w4 is useless when using -N.

If we use following script without -N, we are able to find the actual number of data points inside the width=4km (+w4)

gmt begin figure1 png
gmt basemap -R0/30/5/20 -JX15c/-5c -BWSen -Byaf+l"Focal depth (km)" -Bxaf+l"Distance (km)"
gmt coupe nichong.txt -Sa0.9c -Aa78.43/41.26/78.70/41.05+d90+w4 -W0.4p -GLIGHTSEAGREEN -Q
gmt end

So this bug makes me can not plot data points inside the width=4km (+w4) and above depth 5km, in a map with depth range 5km ~ 20km

[Esteban82]

Hi @CovMat. I'm not quite sure what behavior you consider to be a bug. Could you tell me specifically what problem you want to solve?

[CovMat]

Hi @CovMat. I'm not quite sure what behavior you consider to be a bug. Could you tell me specifically what problem you want to solve?

@Esteban82, sorry I updated the example in above comment to explain this bug. Thank you for reading it

[joa-quim]

It took me a little while to fully understand this issue. We certainly do not want to change the -N meaning/description, which is transversal to several GMT modules and means simply - clip by map borders.