*! graph commands for fig2, stb hidlin insert,  Guy van Melle
*!
prog def fig2 
*
	version 5.0
*------ need some data, convention: 0 means outside grid 
	tempname A
	#delim ;
	mat `A'=(
		 52,46,37,26,30,38,43,45,47,50
		\ 0, 0,37,41,50,55,56,56,55, 0 
		\ 0,46,47,46,44,46,52,61, 0, 0 )
	; #delim cr
	loc nc= colsof(`A')
	loc nr= rowsof(`A')
	loc L =`nr'*`nc'
	loc nzer=0
	while `A'[2,`nzer'+1]==0 { loc nzer=`nzer'+1 }
*------	demo is for crossing any place on LAST row (nr) !
	loc cros 0
	loc ok 0
	loc i 1
	while `i'<`nc'-1 & !`ok' {
		loc i= `i'+1
		if `A'[`nr'-1,`i']>0 & `A'[`nr',`i'] < `A'[`nr'-1,`i'] {
			loc ok 1
			loc cros=`i'
		}
	}

	drop _all
qui{

*------ want a circle 
	set obs 101
	g a =_n*_pi/50
	g v = sin(a)
	g u = cos(a)

*------ data long layout 
	g k = 1+int((_n-1)/`nc')
	g x = 1+mod((_n-1),`nc')
	g y = .

*------ transfer matrix data 
	loc i 0
	while `i' < `L' {
		loc i = `i'+1
		replace y = `A'[k,x] in `i'
	}
	replace y=. if y==0

*--- copy the 4 points to be zoomed 
	loc L1= (`nr'-2)*`nc'+`cros'
	loc L2= (`nr'-1)*`nc'+`cros'	

	replace y = y[`L1'-1]	if _n==`L'+1
	replace y = y[`L1'-0]	if _n==`L'+2
	replace y = y[`L2'-1]	if _n==`L'+3
	replace y = y[`L2'-0]	if _n==`L'+4

	su y if _n> `L'
	replace y = 2*(y-_result(5)) 	if _n>`L' 
	replace x = 2*(x-2*(x>2))+`nc' 	if _n>`L'
	replace k = k+1			if _n> `L'+2
	replace x = . 			if _n> `L'+4
	replace k = . 			if _n> `L'+4 | y==.
	su x
	loc xrg=_result(6)-_result(5)
	su y
	loc yrg=_result(6)-_result(5)

	replace u= 1.05*(u-.5)+x[`cros']
	replace v= 2.10*v*`yrg'/`xrg'+(y[`L1']+y[`L1'-1]+y[`L2']+y[`L2'-1])/4

	set gr off
	gr y x
	set gr on

	loc ay=_result(5)
	loc by=_result(6)
	loc ax=_result(7)
	loc bx=_result(8)

	replace y= .95*(y- 4)  if _n<=`L'
	replace v= .95*(v- 4) 
	replace y= y+12	 if _n> `L'
	replace x= x-1.5 if _n> `L'
	replace y= int(`ay'* y   +`by')
	replace v= int(`ay'* v   +`by')
	replace x= int(`ax'*(x-1)+`bx')
	replace u= int(`ax'*(u-1)+`bx')

*------ find intersection for hidden part 
	loc i 0
	while `i'<4 {
		loc i= `i'+1
		loc y`i'= y[`L'+`i']
		loc x`i'= x[`L'+`i']
	}
	loc alf= (`y3'-`y1') / (`y3'-`y1'-`y4'+`y2')
	loc xP = `x1' + `alf'* (`x2'-`x1')
	loc yP = `y1' + `alf'* (`y2'-`y1')
}
*-------
gph open, saving(3dfig2.gph,replace)
	gph font 850 380
	gph text 1500 800 0 -1 Hidlin: example showing 3 lines	
	gph font 650 280
	loc lox   400 
	loc hix=32000-`lox'
	loc loy   300
	loc hiy=22000-`loy'

*------ draw a box and the circle 
	gph  line `loy' `lox' `loy' `hix'
	gph  line `loy' `hix' `hiy' `hix'
	gph  line `hiy' `hix' `hiy' `lox'
	gph  line `hiy' `lox' `loy' `lox'
	gph vline  v u
	qui {
		su u
		loc uP=_result(3)-300
		su v
		loc vP=_result(3)-700
	}
	gph text `vP' `uP' 0 0 P	

*------ draw the data, join within rows, label rows as 'lines' 
	loc j 0
	loc p 1
	while `j'< `nr'{
		loc j= `j'+1
		loc i =`nc'*(`j'-1)+`nzer'+3-`j'
		loc yy= y[`i']
		loc xx= x[`i']-3500
		loc p= `p'+1
		gph pen `p'	
		gph vpoint y x if k*`j'==1
		gph vpoint y x if `j'==`nr'-1 & (_n==`L1'-1 | _n==`L1')
		gph vpoint y x if `j'==`nr'   & (_n==`L2'-1 | _n==`L2')
		gph vline  y x if k==`j'
		gph text `yy' `xx' 0 -1 curve `j'
		
		loc yy= `yy'+800
		if `j'< `nr' {
			gph text `yy' `xx' 0 -1 visible
		}
		else     {
			gph text `yy' `xx' 0 -1 partly
			loc yy= `yy'+800
			gph text `yy' `xx' 0 -1 hidden
		}		
	}
	su v
	loc yy = _result(5)
	loc xx = x[`cros'-1]
	gph text `yy' `xx' 0 1 disappears
	loc yy = y[`L'+1-`nr']+500
	loc xx = x[`L'+1-`nr']+200
	gph text `yy' `xx' 0 -1 may reappear

*------ show vertical alignment 
	tempname B
	mat `B'=(0,1,4,5,4,6,5,6,3,0)
	gph pen 1
	loc d1 330
	loc d2 450
	loc k 0
	while `k'<10 {
		loc k=`k'+1
		loc xx=x[`k']-50
		loc t1=y[`k']+`d1'
		loc t2= `t1' +`d2'
		gph line `t1' `xx' `t2' `xx'
		loc t2=y[`k']
		loc j 0
		while `j'<`B'[1,`k'] {
			loc j=`j'+1
			loc t1= `t2' -`d1'
			loc t2= `t1' -`d2'
			gph line `t1' `xx' `t2' `xx'
		}
	}
	
 
*------ redraw a zoom of the circled area 

	loc j= `nr'+1
	loc p= `j'-1
	gph pen `p'	
	gph vpoint y x if k==`j'
	gph vline  y x if k==`j'

	su y if _n>`L'
	loc yy= _result(5)-1200
	su x if _n>`L'
	loc xx= (_result(5)+_result(6))/2-400
	gph pen 1
	gph text  `yy' `xx' 0 0 INTERPOLATION  
	gph point `yP' `xP' 400 4
	loc yy=`yP'-1000
	gph font 850 380
	gph text  `yy' `xP' 0 0 P
	gph font 650 280
	loc j= `nr'+2
	loc p= `j'-1
	gph pen `p'	
	gph vpoint y x if k==`j'
	gph  line  `y3' `x3' `yP' `xP'

*------ label points 
	loc i 0
	loc s -1
	while `i'<4 {
		loc i = `i'+1
		loc s = - `s'
		loc t = (-1)^`i'
		loc yy= `y`i''-`s'*100*`t'
		loc xx= `x`i''-`s'*200-200	
		loc p = `nr'+(`i'>2)
		gph pen `p'
		gph text `yy' `xx' 0 `s' `i'
		if mod(`i',2)==0 { 
			loc cu=1+`i'/2
			loc cu curve `cu'
			loc xx=`xx'+1500
			gph text `yy' `xx' 0 -1 `cu'
		}
	}
	
*------ show dotted line to indicate what's hidden 

	loc p 1
	while `alf'< 1 {                  
		loc alf= min(`alf'+.03*(1+1.5*`p'),1)
		loc p= 1-`p'
		gph pen `p'

		loc xQ = `x3'+`alf'*(`x4'-`x3')
		loc yQ = `y3'+`alf'*(`y4'-`y3')		
		gph  line  `yP' `xP' `yQ' `xQ'
		loc yP= `yQ'
		loc xP= `xQ'
	} 		
	gph pen 1
	loc xx   19600
	gph text 18000 `xx' 0 -1 xP = x1 + alfa * ( x2-x1 )
	gph text 19200 `xx' 0 -1 yP = y1 + alfa * ( y2-y1 )
	loc xx= `xx'-2500
	gph text 20600 `xx' 0 -1 where   alfa= ( y3-y1 ) / ( y3-y1 + y2-y4 )

gph close
*--------
end
*