1 /***************************************************************************
2  * base.cpp is part of Math Graphic Library
3  * Copyright (C) 2007-2016 Alexey Balakin <mathgl.abalakin@gmail.ru>       *
4  *                                                                         *
5  *   This program is free software; you can redistribute it and/or modify  *
6  *   it under the terms of the GNU Lesser General Public License  as       *
7  *   published by the Free Software Foundation; either version 3 of the    *
8  *   License, or (at your option) any later version.                       *
9  *                                                                         *
10  *   This program is distributed in the hope that it will be useful,       *
11  *   but WITHOUT ANY WARRANTY; without even the implied warranty of        *
12  *   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the         *
13  *   GNU General Public License for more details.                          *
14  *                                                                         *
15  *   You should have received a copy of the GNU Lesser General Public     *
16  *   License along with this program; if not, write to the                 *
17  *   Free Software Foundation, Inc.,                                       *
18  *   59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.             *
19  ***************************************************************************/
20 #include "mgl2/font.h"
21 #include "mgl2/base.h"
22 #include "mgl2/eval.h"
23 
24 #if MGL_HAVE_FREETYPE
25 #include <ft2build.h>
26 #include FT_FREETYPE_H
27 #include FT_OUTLINE_H
28 #include FT_BBOX_H
29 #endif
30 
31 #if MGL_HAVE_OMP
32 #include <omp.h>
33 #endif
34 
35 //-----------------------------------------------------------------------------
36 static unsigned mgl_pb=0;
mgl_bsize(unsigned bsize)37 unsigned MGL_EXPORT mgl_bsize(unsigned bsize)
38 {
39 	if(!mgl_pb)	mgl_pb = (bsize>0 && bsize<100)?bsize:16;
40 	return mgl_pb;
41 }
mgl_bsize_(unsigned * bsize)42 unsigned MGL_EXPORT mgl_bsize_(unsigned *bsize)
43 {	return mgl_bsize(*bsize);	}
44 //-----------------------------------------------------------------------------
mgl_mutex_unlock(void * mutex)45 void MGL_EXPORT mgl_mutex_unlock(void *mutex)
46 {
47 #if MGL_HAVE_PTHREAD
48 	pthread_mutex_unlock((pthread_mutex_t *)mutex);
49 #elif MGL_HAVE_OMP
50 	omp_unset_lock((omp_lock_t *)mutex);
51 #endif
52 }
53 //-----------------------------------------------------------------------------
mgl_mutex_lock(void * mutex)54 void MGL_EXPORT mgl_mutex_lock(void *mutex)
55 {
56 #if MGL_HAVE_PTHREAD
57 	pthread_mutex_lock((pthread_mutex_t *)mutex);
58 #elif MGL_HAVE_OMP
59 	omp_set_lock((omp_lock_t *)mutex);
60 #endif
61 }
62 //-----------------------------------------------------------------------------
mgl_strdup(const char * s)63 char *mgl_strdup(const char *s)
64 {
65 	char *r = (char *)malloc((strlen(s)+1)*sizeof(char));
66 	if(r)	memcpy(r,s,(strlen(s)+1)*sizeof(char));
67 	return r;
68 }
69 //-----------------------------------------------------------------------------
mgl_create_cpp_font(HMGL gr,const wchar_t * how)70 void MGL_EXPORT mgl_create_cpp_font(HMGL gr, const wchar_t *how)
71 {
72 	unsigned long l=mgl_wcslen(how), i, n=0, m;
73 	wchar_t ch=*how;
74 	const mglFont *f = gr->GetFont();
75 	std::vector<wchar_t> s;	s.push_back(ch);
76 	for(i=1;i<l;i++)
77 		if(how[i]==',')	continue;
78 		else if(how[i]=='-' && i+1<l)
79 			for(ch++;ch<how[i+1];ch++)	s.push_back(ch);
80 		else	s.push_back(ch=how[i]);
81 	for(i=l=n=0;i<s.size();i++)
82 	{
83 		ch = f->Internal(s[i]);
84 		if(ch>=0)	{	l += 2*f->GetNl(0,ch);	n += 6*f->GetNt(0,ch);	}
85 	}
86 	printf("const unsigned long mgl_numg=%lu, mgl_cur=%lu;\n",(unsigned long)s.size(),l+n);
87 	printf("const float mgl_fact=%g;\n",f->GetFact(0)/mgl_fgen);
88 	printf("long mgl_gen_fnt[%lu][6] = {\n", (unsigned long)s.size());
89 	for(i=m=0;i<s.size();i++)	// first write symbols descriptions
90 	{
91 		ch = f->Internal(s[i]);
92 		if(ch<0)	continue;
93 		int m1 = f->GetNl(0,ch), m2 = f->GetNt(0,ch);
94 		printf("\t{0x%x,%d,%d,%lu,%d,%lu},\n",unsigned(s[i]),f->GetWidth(0,ch),m1,m,m2,m+2*m1);
95 		m += 2*m1+6*m2;
96 	}
97 	if(m!=l+n)	printf("#error \"%lu !=%lu + %lu\"",m,l,n);
98 	printf("};\nshort mgl_buf_fnt[%lu] = {\n",m);
99 	for(i=0;i<s.size();i++)		// now write data itself
100 	{
101 		ch = f->Internal(s[i]);
102 		if(ch<0)	continue;
103 		unsigned m1 = f->GetNl(0,ch), m2 = f->GetNt(0,ch);
104 		const short *ln = f->GetLn(0,ch), *tr = f->GetTr(0,ch);
105 		for(l=0;l<2*m1;l++)	printf("%d,",ln[l]);
106 		printf("\n");
107 		for(l=0;l<6*m2;l++)	printf("%d,",tr[l]);
108 		printf("\n");
109 	}
110 	printf("};\n");
111 }
112 //-----------------------------------------------------------------------------
mgl_strtrim(char * str)113 void MGL_EXPORT mgl_strtrim(char *str)
114 {
115 	if(!str || *str==0)	return;
116 	size_t n=strlen(str), k, i;
117 	for(k=0;k<n;k++)	if(str[k]>' ')	break;
118 	for(i=n;i>k;i--)	if(str[i-1]>' ')	break;
119 	memmove(str, str+k, (i-k));
120 	str[i-k]=0;
121 }
122 //-----------------------------------------------------------------------------
mgl_strlwr(char * str)123 void MGL_EXPORT mgl_strlwr(char *str)
124 {
125 	size_t l=strlen(str);
126 	for(size_t k=0;k<l;k++)
127 		str[k] = (str[k]>='A' && str[k]<='Z') ? str[k]+'a'-'A' : str[k];
128 }
129 //-----------------------------------------------------------------------------
mglBase()130 mglBase::mglBase()
131 {
132 	mgl_init();
133 	Flag=0;	saved=false;	PrmInd=NULL;
134 	limit_pm1 = false;
135 #if MGL_HAVE_PTHREAD
136 	pthread_mutex_init(&mutexPnt,0);
137 	pthread_mutex_init(&mutexTxt,0);
138 	pthread_mutex_init(&mutexSub,0);
139 	pthread_mutex_init(&mutexLeg,0);
140 	pthread_mutex_init(&mutexPrm,0);
141 	pthread_mutex_init(&mutexPtx,0);
142 	pthread_mutex_init(&mutexStk,0);
143 	pthread_mutex_init(&mutexGrp,0);
144 	pthread_mutex_init(&mutexGlf,0);
145 	pthread_mutex_init(&mutexAct,0);
146 	pthread_mutex_init(&mutexDrw,0);
147 	pthread_mutex_init(&mutexClf,0);
148 	Pnt.set_mutex(&mutexClf);
149 	Prm.set_mutex(&mutexClf);
150 //	Txt.set_mutex(&mutexClf);
151 #endif
152 #if MGL_HAVE_OMP
153 	lockClf = new omp_lock_t;
154 	omp_init_lock((omp_lock_t*)lockClf);
155 	Pnt.set_mutex(lockClf);
156 	Prm.set_mutex(lockClf);
157 //	Txt.set_mutex(&lockClf);
158 #else
159 	lockClf = NULL;
160 #endif
161 	fnt=0;	*FontDef=0;	fx=fy=fz=fa=fc=0;
162 	AMin.Set(0,0,0,0);	AMax.Set(1,1,1,1);
163 
164 	InUse = 1;	SetQuality();	FaceNum = 0;
165 	// Always create default palette txt[0] and default scheme txt[1]
166 	mglTexture t1(MGL_DEF_PAL,-1), t2(MGL_DEF_SCH,1);
167 	Txt.reserve(3);
168 	MGL_PUSH(Txt,t1,mutexTxt);
169 	MGL_PUSH(Txt,t2,mutexTxt);
170 
171 	strcpy(last_style,"__1 {dFFFF}k\0");
172 	MinS.Set(-1,-1,-1);	MaxS.Set(1,1,1);
173 	fnt = new mglFont;	fnt->gr = this;	PrevState=size_opt=NAN;
174 }
175 //-----------------------------------------------------------------------------
~mglBase()176 mglBase::~mglBase()
177 {
178 	ClearEq();	ClearPrmInd();	delete fnt;
179 	Pnt.set_mutex(0);	Prm.set_mutex(0);	//Txt.set_mutex(0);
180 #if MGL_HAVE_PTHREAD
181 	pthread_mutex_destroy(&mutexPnt);
182 	pthread_mutex_destroy(&mutexTxt);
183 	pthread_mutex_destroy(&mutexSub);
184 	pthread_mutex_destroy(&mutexLeg);
185 	pthread_mutex_destroy(&mutexPrm);
186 	pthread_mutex_destroy(&mutexPtx);
187 	pthread_mutex_destroy(&mutexStk);
188 	pthread_mutex_destroy(&mutexGrp);
189 	pthread_mutex_destroy(&mutexGlf);
190 	pthread_mutex_destroy(&mutexAct);
191 	pthread_mutex_destroy(&mutexDrw);
192 	pthread_mutex_destroy(&mutexClf);
193 #endif
194 #if MGL_HAVE_OMP
195 	omp_destroy_lock((omp_lock_t*)lockClf);
196 	delete ((omp_lock_t*)lockClf);
197 #endif
198 }
199 //-----------------------------------------------------------------------------
SetFontHscale(mreal val)200 void mglBase::SetFontHscale(mreal val)	{	fnt->HeightScale(val);	}
RestoreFont()201 void mglBase::RestoreFont()	{	fnt->Restore();	}
LoadFont(const char * name,const char * path)202 void mglBase::LoadFont(const char *name, const char *path)
203 {	if(name && *name)	fnt->Load(name,path);	else	fnt->Restore();	}
CopyFont(mglBase * gr)204 void mglBase::CopyFont(mglBase *gr)	{	fnt->Copy(gr->GetFont());	}
205 //-----------------------------------------------------------------------------
TextWidth(const char * text,const char * font,mreal size) const206 mreal mglBase::TextWidth(const char *text, const char *font, mreal size) const
207 {	return (size<0?-size*FontSize:size)*font_factor*fnt->Width(text,(font&&*font)?font:FontDef)/20.16;	}
TextWidth(const wchar_t * text,const char * font,mreal size) const208 mreal mglBase::TextWidth(const wchar_t *text, const char *font, mreal size) const
209 {	return (size<0?-size*FontSize:size)*font_factor*fnt->Width(text,(font&&*font)?font:FontDef)/20.16;	}
TextHeight(const char * text,const char * font,mreal size) const210 mreal mglBase::TextHeight(const char *text, const char *font, mreal size) const
211 {	float y1,y2;	fnt->Width(text,(font&&*font)?font:FontDef,&y1,&y2);
212 	return (size<0?-size*FontSize:size)*font_factor*(y2-y1)/20.16;	}
TextHeight(const wchar_t * text,const char * font,mreal size) const213 mreal mglBase::TextHeight(const wchar_t *text, const char *font, mreal size) const
214 {	float y1,y2;	fnt->Width(text,(font&&*font)?font:FontDef,&y1,&y2);
215 	return (size<0?-size*FontSize:size)*font_factor*(y2-y1)/20.16;	}
TextHeight(const char * font,mreal size) const216 mreal mglBase::TextHeight(const char *font, mreal size) const
217 {	return (size<0?-size*FontSize:size)*font_factor*fnt->Height(font?font:FontDef)/20.16; }
AddActive(long k,int n)218 void mglBase::AddActive(long k,int n)
219 {
220 	if(k<0 || (size_t)k>=Pnt.size())	return;
221 	mglActivePos p;
222 	const mglPnt &q=Pnt[k];
223 	int h=GetHeight();
224 	p.x = int(q.x);	p.y = h>1?h-1-int(q.y):int(q.y);
225 	p.id = ObjId;	p.n = n;
226 #pragma omp critical(act)
227 	MGL_PUSH(Act,p,mutexAct);
228 }
229 //-----------------------------------------------------------------------------
GetRatio() const230 mreal mglBase::GetRatio() const	{	return 1;	}
GetWidth() const231 int mglBase::GetWidth()  const	{	return 1;	}
GetHeight() const232 int mglBase::GetHeight() const	{	return 1;	}
233 //-----------------------------------------------------------------------------
StartGroup(const char * name,int id)234 void mglBase::StartGroup(const char *name, int id)
235 {
236 	LightScale(&B);
237 	char buf[128];
238 	snprintf(buf,128,"%s_%d",name,id);	buf[127]=0;
239 	StartAutoGroup(buf);
240 }
241 //-----------------------------------------------------------------------------
242 const char *mglWarn[mglWarnEnd] = {_("data dimension(s) is incompatible"),	//mglWarnDim
243 								_("data dimension(s) is too small"),		//mglWarnLow
244 								_("minimal data value is negative"),		//mglWarnNeg
245 								_("no file or wrong data dimensions"),		//mglWarnFile
246 								_("not enough memory"), 					//mglWarnMem
247 								_("data values are zero"),					//mglWarnZero
248 								_("no legend entries"),					//mglWarnLeg
249 								_("slice value is out of range"),			//mglWarnSlc
250 								_("number of contours is zero or negative"),//mglWarnCnt
251 								_("couldn't open file"),					//mglWarnOpen
252 								_("light: ID is out of range"),			//mglWarnLId
253 								_("size(s) is zero or negative"),			//mglWarnSize
254 								_("format is not supported for that build"),//mglWarnFmt
255 								_("axis ranges are incompatible"),			//mglWarnTern
256 								_("pointer is NULL"),						//mglWarnNull
257 								_("not enough space for plot"),			//mglWarnSpc
258 								_("There is wrong argument(s) in script"),	//mglScrArg
259 								_("There is wrong command(s) in script"),	//mglScrCmd
260 								_("There is too long string(s) in script"),	//mglScrLong
261 								_("There is unbalanced ' in script"),		//mglScrStr
262 								_("There is changing temporary data in script")};	//mglScrTemp
263 //-----------------------------------------------------------------------------
264 extern bool mglPrintWarn;
SetWarn(int code,const char * who)265 void mglBase::SetWarn(int code, const char *who)
266 {
267 	std::string warn;
268 	WarnCode = code>0 ? code:0;
269 	if(code>0 && code<mglWarnEnd)
270 	{
271 		if(who && *who)	warn = std::string(who)+": ";
272 		warn = warn+mglWarn[code-1];
273 	}
274 	else if(!code)	Mess="";
275 	else if(who && *who)	warn = who;
276 	if(mglPrintWarn && !warn.empty())	fprintf(stderr,_("MathGL message - %s\n"),warn.c_str());
277 	if(code && !warn.empty())	Mess = Mess+(code==-2?"":"\n")+warn;
278 	LoadState();
279 }
280 //-----------------------------------------------------------------------------
281 //		Add glyph to the buffer
282 //-----------------------------------------------------------------------------
283 #if MGL_HAVE_FREETYPE
mgl_glf_moveto(const FT_Vector * to,void * user)284 MGL_NO_EXPORT int mgl_glf_moveto(const FT_Vector *to, void *user)
285 {
286 	std::vector<double> *a = (std::vector<double> *)user;
287 	if(a->size()>1)	{	a->push_back(NAN);	a->push_back(NAN);	}
288 	a->push_back(to->x);	a->push_back(to->y);
289 	return 0;
290 }
mgl_glf_lineto(const FT_Vector * to,void * user)291 MGL_NO_EXPORT int mgl_glf_lineto(const FT_Vector *to, void *user)
292 {
293 	std::vector<double> *a = (std::vector<double> *)user;
294 	if(a->size()<2)	{	a->push_back(0);	a->push_back(0);	}
295 	a->push_back(to->x);	a->push_back(to->y);
296 	return 0;
297 }
mgl_glf_parabto(const FT_Vector * ctrl,const FT_Vector * to,void * user)298 MGL_NO_EXPORT int mgl_glf_parabto(const FT_Vector *ctrl, const FT_Vector *to, void *user)
299 {
300 	std::vector<double> *a = (std::vector<double> *)user;
301 	int x0=0,y0=0, x1=ctrl->x,y1=ctrl->y, x2=to->x,y2=to->y;
302 	size_t n=a->size();
303 	if(n<2)	{	a->push_back(0);	a->push_back(0);	}
304 	else	{	x0 = (*a)[n-2];	y0 = (*a)[n-1];	}
305 	for(int i=1;i<=10;i++)
306 	{
307 		double t = 0.1*i;
308 		a->push_back(x0*(1-t)*(1-t)+2*(1-t)*t*x1+t*t*x2);
309 		a->push_back(y0*(1-t)*(1-t)+2*(1-t)*t*y1+t*t*y2);
310 	}
311 	return 0;
312 }
mgl_glf_cubeto(const FT_Vector * ctrl1,const FT_Vector * ctrl2,const FT_Vector * to,void * user)313 MGL_NO_EXPORT int mgl_glf_cubeto(const FT_Vector *ctrl1, const FT_Vector *ctrl2, const FT_Vector *to, void *user)
314 {
315 	std::vector<double> *a = (std::vector<double> *)user;
316 	int x0=0,y0=0, x1=ctrl1->x,y1=ctrl1->y, x2=ctrl2->x,y2=ctrl2->y, x3=to->x,y3=to->y;
317 	size_t n=a->size();
318 	if(n<2)	{	a->push_back(0);	a->push_back(0);	}
319 	else	{	x0 = (*a)[n-2];	y0 = (*a)[n-1];	}
320 	for(int i=1;i<=30;i++)
321 	{
322 		double t = 0.1*i;
323 		a->push_back(x0*(1-t)*(1-t)*(1-t) + 3*(1-t)*t*(x1*(1-t)+t*x2) + t*t*t*x3);
324 		a->push_back(y0*(1-t)*(1-t)*(1-t) + 3*(1-t)*t*(y1*(1-t)+t*y2) + t*t*t*y3);
325 	}
326 	return 0;
327 }
328 #endif
Load(wchar_t id,const char * fname)329 void mglGlyph::Load(wchar_t id, const char *fname)
330 {
331 #if MGL_HAVE_FREETYPE
332 	FT_Library m_ftLibrary;
333 	if(FT_Init_FreeType(&m_ftLibrary))
334 		mgl_set_global_warn("Couldn't initialize the FreeType library.");
335 	else
336 	{
337 		FT_Face m_face;
338 		// For simplicity, always use the first face index.
339 		if(FT_New_Face(m_ftLibrary, fname, 0, &m_face))
340 			mgl_set_global_warn("Couldn't load the font file.");
341 		else
342 		{
343 			// For simplicity, use the charmap FreeType provides by default;
344 			// in most cases this means Unicode.
345 			FT_UInt index = FT_Get_Char_Index(m_face, id);
346 			if(FT_Load_Glyph(m_face, index, FT_LOAD_NO_SCALE|FT_LOAD_NO_BITMAP))
347 				mgl_set_global_warn("Couldn't load the glyph.");
348 			else
349 			{
350 				FT_GlyphSlot slot = m_face->glyph;
351 				FT_Outline &outline = slot->outline;
352 
353 				if(slot->format!=FT_GLYPH_FORMAT_OUTLINE || outline.n_contours <= 0 || outline.n_points <= 0 ||FT_Outline_Check(&outline))
354 					mgl_set_global_warn("Outline doesn't exist.");
355 				else
356 				{
357 					const FT_Fixed multiplier = 32768L;
358 					FT_Matrix matrix;
359 					matrix.xx = matrix.yy = multiplier;
360 					matrix.xy = matrix.yx = 0;
361 					FT_Outline_Transform(&outline, &matrix);
362 
363 					FT_Outline_Funcs callbacks;
364 					callbacks.move_to = mgl_glf_moveto;
365 					callbacks.line_to = mgl_glf_lineto;
366 					callbacks.conic_to = mgl_glf_parabto;
367 					callbacks.cubic_to = mgl_glf_cubeto;
368 					callbacks.shift = 0;
369 					callbacks.delta = 0;
370 					std::vector<double> xy_coor;
371 					if(FT_Outline_Decompose(&outline, &callbacks, &xy_coor))
372 						mgl_set_global_warn("Couldn't extract the outline.");
373 					nt = -id;	// TODO optimize and copy points. Q: actual width? Q: cmp with known.
374 					FT_BBox boundingBox;
375 					FT_Outline_Get_BBox(&outline, &boundingBox);
376 					FT_Pos xMin = boundingBox.xMin;
377 					FT_Pos yMin = boundingBox.yMin;
378 					FT_Pos xMax = boundingBox.xMax;
379 					FT_Pos yMax = boundingBox.yMax;
380 /*					m_xMin = xMin;
381 					m_yMin = yMin;
382 					m_width = xMax - xMin;
383 					m_height = yMax - yMin;*/
384 
385 				}
386 			}
387 			FT_Done_Face(m_face);
388 		}
389 	}
390 	FT_Done_FreeType(m_ftLibrary);
391 #endif
392 }
393 //-----------------------------------------------------------------------------
Create(long Nt,long Nl)394 void mglGlyph::Create(long Nt, long Nl)
395 {
396 //	if(Nt<0 || Nl<0)	return;
397 	nt=Nt;	nl=Nl;
398 #pragma omp critical(glf_create)
399 	{
400 		if(trig)	delete []trig;
401 		trig = nt>0?new short[6*nt]:0;
402 		if(line)	delete []line;
403 		line = nl>0?new short[2*nl]:0;
404 	}
405 }
406 //-----------------------------------------------------------------------------
operator ==(const mglGlyph & g) const407 bool mglGlyph::operator==(const mglGlyph &g) const
408 {
409 	if(nl!=g.nl || nt!=g.nt)	return false;
410 	if(trig && memcmp(trig,g.trig,6*nt*sizeof(short)))	return false;
411 	if(line && memcmp(line,g.line,2*nl*sizeof(short)))	return false;
412 	return true;
413 }
414 //-----------------------------------------------------------------------------
AddGlyph(int s,long j)415 long mglBase::AddGlyph(int s, long j)
416 {
417 	// first create glyph for current typeface
418 	s = s&3;
419 	mglGlyph g(fnt->GetNt(s,j), fnt->GetNl(s,j));
420 	memcpy(g.trig, fnt->GetTr(s,j), 6*g.nt*sizeof(short));
421 	memcpy(g.line, fnt->GetLn(s,j), 2*g.nl*sizeof(short));
422 	// now let find the similar glyph
423 	for(size_t i=0;i<Glf.size();i++)
424 		if(g!=Glf[i])	continue;	else	return i;
425 	// if no one then let add it
426 	long k;
427 #pragma omp critical(glf)
428 	{k=Glf.size();	MGL_PUSH(Glf,g,mutexGlf);}	return k;
429 }
430 //-----------------------------------------------------------------------------
AddGlyph(unsigned char id)431 long mglBase::AddGlyph(unsigned char id)
432 {
433 	size_t j=0;
434 	for(size_t i=0;i<UserGlf.size();i++)
435 		if(UserGlf[i].nt==-id)	j=i+1;
436 	if(j==0)	return -1;
437 	const mglGlyph &g=UserGlf[j-1];
438 	// let find the similar glyph
439 	for(size_t i=0;i<Glf.size();i++)
440 		if(g!=Glf[i])	continue;	else	return i;
441 	// if no one then let add it
442 	long k;
443 #pragma omp critical(glf)
444 	{k=Glf.size();	MGL_PUSH(Glf,g,mutexGlf);}	return k;
445 }
446 //-----------------------------------------------------------------------------
DefineGlyph(HCDT x,HCDT y,unsigned char id)447 void mglBase::DefineGlyph(HCDT x, HCDT y, unsigned char id)
448 {
449 	long n = x->GetNx();
450 	if(y->GetNx()!=n || n<2)	return;
451 	mglGlyph g(-id,n);
452 	mreal x1=1e10,x2=-1e10,y1=1e10,y2=-1e10;
453 	for(long i=0;i<n;i++)
454 	{
455 		mreal xx = x->v(i), yy = y->v(i);
456 		x1=x1>xx?xx:x1;	x2=x2<xx?xx:x2;
457 		y1=y1>yy?yy:y1;	y2=y2<yy?yy:y2;
458 	}
459 	mreal scale = 1;
460 	if(fabs(x1)<10 && fabs(x2)<10 && fabs(y1)<10 && fabs(y2)<10)
461 		scale=1000;
462 	for(long i=0;i<n;i++)
463 	{
464 		short sx = short(x->v(i)*scale), sy = short(y->v(i)*scale);
465 		g.line[2*i] = sx;	g.line[2*i+1] = sy;
466 	}
467 	UserGlf.push_back(g);
468 }
469 //-----------------------------------------------------------------------------
470 //		Add points to the buffer
471 //-----------------------------------------------------------------------------
PushPnts(size_t num,const mglPnt * qq)472 long mglBase::PushPnts(size_t num, const mglPnt *qq)
473 {
474 	long k;
475 #pragma omp critical(pnt)
476 	MGL_PUSHs({k=Pnt.size();Pnt.push_back(num,qq);},mutexPnt);
477 	return k;
478 }
479 //-----------------------------------------------------------------------------
AllocPnts(size_t num)480 long mglBase::AllocPnts(size_t num)
481 {
482 	long k;
483 #pragma omp critical(pnt)
484 	MGL_PUSHs({k=Pnt.allocate(num);},mutexPnt);
485 	return k;
486 }
487 //-----------------------------------------------------------------------------
mgl_put_inbox(mreal a1,mreal a2,mreal & a)488 void inline mgl_put_inbox(mreal a1, mreal a2, mreal &a)
489 {
490 	if(a1<a2)	{	if(a<a1)	a=a1;	if(a>a2)	a=a2;	}
491 	else		{	if(a<a2)	a=a2;	if(a>a1)	a=a1;	}
492 }
mgl_coor_box(HMGL gr,mglPoint & p)493 static void mgl_coor_box(HMGL gr, mglPoint &p)
494 {
495 	mgl_put_inbox(gr->Min.x, gr->Max.x, p.x);
496 	mgl_put_inbox(gr->Min.y, gr->Max.y, p.y);
497 	mgl_put_inbox(gr->Min.z, gr->Max.z, p.z);
498 }
AddPnt(const mglMatrix * mat,const mglPoint & p,mreal c,const mglPoint & n,mreal a,int scl)499 long mglBase::AddPnt(const mglMatrix *mat, const mglPoint &p, mreal c, const mglPoint &n, mreal a, int scl)
500 {
501 	mglPnt q;
502 	if(!AddPntQ(q,mat,p,c,n,a,scl))	return -1;
503 	long k;
504 #pragma omp critical(pnt)
505 	{k=Pnt.size();	MGL_PUSH(Pnt,q,mutexPnt);}	return k;
506 }
AddPntQ(mglPnt & q,const mglMatrix * mat,mglPoint p,mreal c,mglPoint n,mreal a,int scl)507 bool mglBase::AddPntQ(mglPnt &q, const mglMatrix *mat, mglPoint p, mreal c, mglPoint n, mreal a, int scl)
508 {
509 	// scl=0 -- no scaling
510 	// scl&1 -- usual scaling
511 	// scl&2 -- disable NAN at scaling
512 	// scl&4 -- disable NAN for normales if no light
513 	// scl&8 -- bypass palette for enabling alpha
514 	// scl&16 -- put points inside axis range
515 	if(mgl_isnan(c) || mgl_isnan(a))	{	q.x=NAN;	return false;	}
516 	bool norefr = mgl_isnan(n.x) && mgl_isnan(n.y) && !mgl_isnan(n.z), res=true;
517 	if(scl>0)
518 	{
519 		if(scl&16)	mgl_coor_box(this, p);
520 		res = ScalePoint(mat,p,n,!(scl&2));
521 	}
522 	if(mgl_isnan(p.x))	{	q.x=NAN;	return false;	}
523 	a = (a>=0 && a<=1) ? a : AlphaDef;
524 	c = (c>=0) ? c:CDef;
525 
526 	if(get(MGL_REDUCEACC))
527 	{
528 		q.x=q.xx=int(p.x*10)*0.1;	q.y=q.yy=int(p.y*10)*0.1;	q.z=q.zz=int(p.z*10)*0.1;
529 		q.c=int(c*100)*0.01;		q.ta=int(a*100)*0.01;
530 		q.u=mgl_isnum(n.x)?int(n.x*100)*0.01:NAN;
531 		q.v=mgl_isnum(n.y)?int(n.y*100)*0.01:NAN;
532 		q.w=mgl_isnum(n.z)?int(n.z*100)*0.01:NAN;
533 	}
534 	else
535 	{
536 		q.x=q.xx=p.x;	q.y=q.yy=p.y;	q.z=q.zz=p.z;
537 		q.c=c;	q.ta=a;	q.u=n.x;	q.v=n.y;	q.w=n.z;
538 	}
539 	long ci=long(c);
540 	if(ci<0 || ci>=(long)Txt.size())	ci=0;	// NOTE never should be here!!!
541 	const mglTexture &txt=Txt[ci];
542 	txt.GetC(c,a,q);	// RGBA color
543 	if(get(MGL_GRAY_MODE))
544 	{
545 		float h = 0.3*q.r + 0.59*q.g + 0.11*q.b;
546 		q.r = q.g = q.b = h;
547 	}
548 
549 	// add gap for texture coordinates for compatibility with OpenGL
550 	const mreal gap = 0./MGL_TEXTURE_COLOURS;
551 	q.c = ci+(q.c-ci)*(1-2*gap)+gap;
552 	q.ta = q.ta*(1-2*gap)+gap;
553 
554 	if(scl&8 && scl>0)	q.a=a;	// bypass palette for enabling alpha in Error()
555 	if(!get(MGL_ENABLE_ALPHA))	{	q.a=1;	if(txt.Smooth!=2)	q.ta=1-gap;	}
556 	if(norefr)	q.v=0;
557 	if(!get(MGL_ENABLE_LIGHT) && !(scl&4))	q.u=q.v=NAN;
558 	q.sub=mat->norot?-1*(int)Sub.size():Sub.size()-1;
559 	return (scl&16)?res:true;
560 }
561 //-----------------------------------------------------------------------------
CopyNtoC(long from,mreal c)562 long mglBase::CopyNtoC(long from, mreal c)
563 {
564 	mglPnt q;
565 	if(!CopyNtoC(q,from,c))	return -1;
566 	long k;
567 #pragma omp critical(pnt)
568 	{k=Pnt.size();	MGL_PUSH(Pnt,q,mutexPnt);}	return k;
569 }
570 //-----------------------------------------------------------------------------
CopyNtoC(mglPnt & q,long from,mreal c)571 bool mglBase::CopyNtoC(mglPnt &q, long from, mreal c)
572 {
573 	if(from<0)	return false;
574 	q = Pnt[from];
575 	if(mgl_isnum(c))	{	q.c=c;	q.ta=1;	Txt[long(c)].GetC(c,0,q);	q.a=1;	}
576 	else	q.x = NAN;
577 	return mgl_isnum(q.x);
578 }
579 //-----------------------------------------------------------------------------
CopyProj(long from,const mglPoint & p,const mglPoint & n,short sub)580 long mglBase::CopyProj(long from, const mglPoint &p, const mglPoint &n, short sub)
581 {
582 	mglPnt q;
583 	if(!CopyProj(q,from,p,n,sub))	return -1;
584 	long k;
585 #pragma omp critical(pnt)
586 	{k=Pnt.size();	MGL_PUSH(Pnt,q,mutexPnt);}	return k;
587 }
588 //-----------------------------------------------------------------------------
CopyProj(mglPnt & q,long from,const mglPoint & p,const mglPoint & n,short sub)589 bool mglBase::CopyProj(mglPnt &q, long from, const mglPoint &p, const mglPoint &n, short sub)
590 {
591 	if(from<0)	return false;
592 	q=Pnt[from];	q.sub = sub;
593 	q.x=q.xx=p.x;	q.y=q.yy=p.y;	q.z=q.zz=p.z;
594 	q.u = n.x;		q.v = n.y;		q.w = n.z;
595 	return mgl_isnum(q.x);
596 }
597 //-----------------------------------------------------------------------------
Reserve(long n)598 void mglBase::Reserve(long n)
599 {
600 	if(TernAxis&12)	n*=4;
601 #pragma omp critical(pnt)
602 	Pnt.reserve(n);
603 #pragma omp critical(prm)
604 	Prm.reserve(n);
605 }
606 //-----------------------------------------------------------------------------
607 //		Boundaries and scaling
608 //-----------------------------------------------------------------------------
RecalcCRange()609 bool mglBase::RecalcCRange()
610 {
611 	bool wrong=false;
612 	if(!fa)
613 	{	FMin.c = Min.c;	FMax.c = Max.c;	}
614 	else
615 	{
616 		FMin.c = INFINITY;	FMax.c = -INFINITY;
617 		int n=30;
618 		for(int i=0;i<=n;i++)
619 		{
620 			mreal a = fa->Calc(0,0,0,Min.c+i*(Max.c-Min.c)/n);
621 			if(mgl_isbad(a))	wrong=true;
622 			if(a<FMin.c)	FMin.c=a;
623 			if(a>FMax.c)	FMax.c=a;
624 		}
625 	}
626 	return wrong;
627 }
628 //-----------------------------------------------------------------------------
RecalcBorder()629 void mglBase::RecalcBorder()
630 {
631 	ZMin = 1.;
632 	bool wrong=false;
633 	if(!fx && !fy && !fz)
634 	{	FMin = Min;	FMax = Max;	}
635 	else
636 	{
637 		FMin.Set( INFINITY, INFINITY, INFINITY);
638 		FMax.Set(-INFINITY,-INFINITY,-INFINITY);
639 		int n=30;
640 		for(int i=0;i<=n;i++)	for(int j=0;j<=n;j++)	// x range
641 		{
642 			if(SetFBord(Min.x, Min.y+i*(Max.y-Min.y)/n, Min.z+j*(Max.z-Min.z)/n))	wrong=true;
643 			if(SetFBord(Max.x, Min.y+i*(Max.y-Min.y)/n, Min.z+j*(Max.z-Min.z)/n))	wrong=true;
644 		}
645 		for(int i=0;i<=n;i++)	for(int j=0;j<=n;j++)	// y range
646 		{
647 			if(SetFBord(Min.x+i*(Max.x-Min.x)/n, Min.y, Min.z+j*(Max.z-Min.z)/n))	wrong=true;
648 			if(SetFBord(Min.x+i*(Max.x-Min.x)/n, Max.y, Min.z+j*(Max.z-Min.z)/n))	wrong=true;
649 		}
650 		for(int i=0;i<=n;i++)	for(int j=0;j<=n;j++)	// x range
651 		{
652 			if(SetFBord(Min.x+i*(Max.x-Min.x)/n, Min.y+j*(Max.y-Min.y)/n, Min.z))	wrong=true;
653 			if(SetFBord(Min.x+i*(Max.x-Min.x)/n, Min.y+j*(Max.y-Min.y)/n, Max.z))	wrong=true;
654 		}
655 		if(!fx)	{	FMin.x = Min.x;	FMax.x = Max.x;	}
656 		else	{	mreal d=0.01*(FMax.x-FMin.x);	FMin.x-=d;	FMax.x+=d;	}
657 		if(!fy)	{	FMin.y = Min.y;	FMax.y = Max.y;	}
658 		else	{	mreal d=0.01*(FMax.y-FMin.y);	FMin.y-=d;	FMax.y+=d;	}
659 		if(!fz)	{	FMin.z = Min.z;	FMax.z = Max.z;	}
660 		else	{	mreal d=0.01*(FMax.z-FMin.z);	FMin.z-=d;	FMax.z+=d;	}
661 	}
662 	if(RecalcCRange())	wrong=true;
663 	if(wrong)	SetWarn(mglWarnTern, "Curved coordinates");
664 }
665 //-----------------------------------------------------------------------------
SetFBord(mreal x,mreal y,mreal z)666 bool mglBase::SetFBord(mreal x,mreal y,mreal z)
667 {
668 	bool wrong=false;
669 	if(fx)
670 	{
671 		mreal v = fx->Calc(x,y,z);
672 		if(mgl_isbad(v))	wrong = true;
673 		if(FMax.x < v)	FMax.x = v;
674 		if(FMin.x > v)	FMin.x = v;
675 	}
676 	if(fy)
677 	{
678 		mreal v = fy->Calc(x,y,z);
679 		if(mgl_isbad(v))	wrong = true;
680 		if(FMax.y < v)	FMax.y = v;
681 		if(FMin.y > v)	FMin.y = v;
682 	}
683 	if(fz)
684 	{
685 		mreal v = fz->Calc(x,y,z);
686 		if(mgl_isbad(v))	wrong = true;
687 		if(FMax.z < v)	FMax.z = v;
688 		if(FMin.z > v)	FMin.z = v;
689 	}
690 	return wrong;
691 }
692 //-----------------------------------------------------------------------------
ScalePoint(const mglMatrix *,mglPoint & p,mglPoint & n,bool use_nan) const693 bool mglBase::ScalePoint(const mglMatrix *, mglPoint &p, mglPoint &n, bool use_nan) const
694 {
695 	mreal &x=p.x, &y=p.y, &z=p.z;
696 	if(mgl_isnan(x) || mgl_isnan(y) || mgl_isnan(z))	{	x=NAN;	return false;	}
697 	mreal x1,y1,z1,x2,y2,z2;
698 	x1 = x>0?x*MGL_EPSILON:x/MGL_EPSILON;	x2 = x<0?x*MGL_EPSILON:x/MGL_EPSILON;
699 	y1 = y>0?y*MGL_EPSILON:y/MGL_EPSILON;	y2 = y<0?y*MGL_EPSILON:y/MGL_EPSILON;
700 	z1 = z>0?z*MGL_EPSILON:z/MGL_EPSILON;	z2 = z<0?z*MGL_EPSILON:z/MGL_EPSILON;
701 	bool res = true;
702 	if(x2>CutMin.x && x1<CutMax.x && y2>CutMin.y && y1<CutMax.y &&
703 		z2>CutMin.z && z1<CutMax.z)	res = false;
704 	if(fc && fc->Calc(x,y,z)!=0)	res = false;
705 
706 	if(get(MGL_ENABLE_CUT) || !use_nan)
707 	{
708 //		if(x1<Min.x || x2>Max.x || y1<Min.y || y2>Max.y || z1<Min.z || z2>Max.z)	res = false;
709 		if((x1-Min.x)*(x1-Max.x)>0 && (x2-Min.x)*(x2-Max.x)>0)	res = false;
710 		if((y1-Min.y)*(y1-Max.y)>0 && (y2-Min.y)*(y2-Max.y)>0)	res = false;
711 		if((z1-Min.z)*(z1-Max.z)>0 && (z2-Min.z)*(z2-Max.z)>0)	res = false;
712 	}
713 	else
714 	{
715 		if(Min.x<Max.x)
716 		{
717 			if(x1<Min.x)	{x=Min.x;	n.Set(1,0,0);}
718 			if(x2>Max.x)	{x=Max.x;	n.Set(1,0,0);}
719 		}
720 		else
721 		{
722 			if(x1<Max.x)	{x=Max.x;	n.Set(1,0,0);}
723 			if(x2>Min.x)	{x=Min.x;	n.Set(1,0,0);}
724 		}
725 		if(Min.y<Max.y)
726 		{
727 			if(y1<Min.y)	{y=Min.y;	n.Set(0,1,0);}
728 			if(y2>Max.y)	{y=Max.y;	n.Set(0,1,0);}
729 		}
730 		else
731 		{
732 			if(y1<Max.y)	{y=Max.y;	n.Set(0,1,0);}
733 			if(y2>Min.y)	{y=Min.y;	n.Set(0,1,0);}
734 		}
735 		if(Min.z<Max.z)
736 		{
737 			if(z1<Min.z)	{z=Min.z;	n.Set(0,0,1);}
738 			if(z2>Max.z)	{z=Max.z;	n.Set(0,0,1);}
739 		}
740 		else
741 		{
742 			if(z1<Max.z)	{z=Max.z;	n.Set(0,0,1);}
743 			if(z2>Min.z)	{z=Min.z;	n.Set(0,0,1);}
744 		}
745 	}
746 
747 	x1=x;	y1=y;	z1=z;
748 	mreal xx=1,xy=0,xz=0,yx=0,yy=1,yz=0,zx=0,zy=0,zz=1;
749 	if(fx)	{	x1 = fx->Calc(x,y,z);	xx = fx->CalcD('x',x,y,z);	xy = fx->CalcD('y',x,y,z);	xz = fx->CalcD('z',x,y,z);	}
750 	if(fy)	{	y1 = fy->Calc(x,y,z);	yx = fy->CalcD('x',x,y,z);	yy = fy->CalcD('y',x,y,z);	yz = fy->CalcD('z',x,y,z);	}
751 	if(fz)	{	z1 = fz->Calc(x,y,z);	zx = fz->CalcD('x',x,y,z);	zy = fz->CalcD('y',x,y,z);	zz = fz->CalcD('z',x,y,z);	}
752 	if(mgl_isnan(x1) || mgl_isnan(y1) || mgl_isnan(z1))	{	x=NAN;	return false;	}
753 
754 	mreal d;
755 	d = 1/(FMax.x - FMin.x);	x = (2*x1 - FMin.x - FMax.x)*d;	xx /= d;	xy /= d;	xz /= d;
756 	d = 1/(FMax.y - FMin.y);	y = (2*y1 - FMin.y - FMax.y)*d;	yx /= d;	yy /= d;	yz /= d;
757 	d = 1/(FMax.z - FMin.z);	z = (2*z1 - FMin.z - FMax.z)*d;	zx /= d;	zy /= d;	zz /= d;
758 	mreal nx=n.x, ny=n.y, nz=n.z;
759 	n.x = nx*xx+ny*xy+nz*xz;
760 	n.y = nx*yx+ny*yy+nz*yz;
761 	n.z = nx*zx+ny*zy+nz*zz;
762 	if((TernAxis&3)==1)	// usual ternary axis
763 	{
764 		if(x+y>0)
765 		{
766 			if(get(MGL_ENABLE_CUT))	res = false;
767 			else	y = -x;
768 		}
769 		x += (y+1)/2;	n.x += n.y/2;
770 	}
771 	else if((TernAxis&3)==2)	// quaternary axis
772 	{
773 		if(x+y+z>-1)
774 		{
775 			if(get(MGL_ENABLE_CUT))	res = false;
776 			else	z = -1-y-x;
777 		}
778 		x += 1+(y+z)/2;		y += (z+1)/3;
779 		n.x += (n.y+n.z)/2;	n.y += n.z/3;
780 	}
781 	if(fabs(x)>MGL_FEPSILON || fabs(y)>MGL_FEPSILON || fabs(z)>MGL_FEPSILON)	res = false;
782 
783 	if(!res && use_nan)	x = NAN;	// extra sign that point shouldn't be plotted
784 	else if(limit_pm1)
785 	{
786 		x = x>1?1:(x<-1?-1:x);
787 		y = y>1?1:(y<-1?-1:y);
788 		z = z>1?1:(z<-1?-1:z);
789 	}
790 	return res;
791 }
792 //-----------------------------------------------------------------------------
793 //		Ranges
794 //-----------------------------------------------------------------------------
mglScaleAxis(mreal & v1,mreal & v2,mreal & v0,mreal x1,mreal x2)795 void mglScaleAxis(mreal &v1, mreal &v2, mreal &v0, mreal x1, mreal x2)
796 {
797 	if(!mgl_isrange(x1,x2) || !mgl_isrange(v1,v2))	return;
798 	mreal dv,d0;	x2-=1;
799 	if(v1*v2>0 && (v2/v1>=100 || v2/v1<=0.01))	// log scale
800 	{
801 		dv=log(v2/v1);	d0 = log(v0/v1)/log(v2/v1);
802 		v1*=exp(dv*x1);	v2*=exp(dv*x2);	v0=v1*exp(d0*log(v2/v1));
803 	}
804 	else
805 	{
806 		dv=v2-v1;	d0=(v0-v1)/(v2-v1);
807 		v1+=dv*x1;	v2+=dv*x2;	v0=v1+d0*(v2-v1);
808 	}
809 }
810 //-----------------------------------------------------------------------------
SetOrigin(mreal x0,mreal y0,mreal z0,mreal c0)811 void mglBase::SetOrigin(mreal x0, mreal y0, mreal z0, mreal c0)
812 {
813 	Org.Set(x0,y0,z0,c0);
814 	if((TernAxis&3)==0)
815 	{
816 		Min = OMin;	Max = OMax;
817 		mglScaleAxis(Min.x, Max.x, Org.x, AMin.x, AMax.x);
818 		mglScaleAxis(Min.y, Max.y, Org.y, AMin.y, AMax.y);
819 		mglScaleAxis(Min.z, Max.z, Org.z, AMin.z, AMax.z);
820 		mglScaleAxis(Min.c, Max.c, Org.c, AMin.c, AMax.c);
821 	}
822 }
823 //-----------------------------------------------------------------------------
SetRanges(const mglPoint & m1,const mglPoint & m2)824 void mglBase::SetRanges(const mglPoint &m1, const mglPoint &m2)
825 {
826 	if(mgl_isrange(m1.x, m2.x))	{	Min.x=m1.x;	Max.x=m2.x;	}
827 	if(mgl_isrange(m1.y, m2.y))	{	Min.y=m1.y;	Max.y=m2.y;	}
828 	if(mgl_isrange(m1.z, m2.z))	{	Min.z=m1.z;	Max.z=m2.z;	}
829 	if(mgl_isrange(m1.c, m2.c))	{	Min.c=m1.c;	Max.c=m2.c;	}
830 	else	{	Min.c=Min.z;Max.c=Max.z;}
831 
832 	if(Org.x<Min.x && mgl_isnum(Org.x))	Org.x = Min.x;
833 	if(Org.x>Max.x && mgl_isnum(Org.x))	Org.x = Max.x;
834 	if(Org.y<Min.y && mgl_isnum(Org.y))	Org.y = Min.y;
835 	if(Org.y>Max.y && mgl_isnum(Org.y))	Org.y = Max.y;
836 	if(Org.z<Min.z && mgl_isnum(Org.z))	Org.z = Min.z;
837 	if(Org.z>Max.z && mgl_isnum(Org.z))	Org.z = Max.z;
838 
839 	if((TernAxis&3)==0)
840 	{
841 		OMax = Max;	OMin = Min;
842 		mglScaleAxis(Min.x, Max.x, Org.x, AMin.x, AMax.x);
843 		mglScaleAxis(Min.y, Max.y, Org.y, AMin.y, AMax.y);
844 		mglScaleAxis(Min.z, Max.z, Org.z, AMin.z, AMax.z);
845 		mglScaleAxis(Min.c, Max.c, Org.c, AMin.c, AMax.c);
846 	}
847 
848 	CutMin.Set(0,0,0);	CutMax.Set(0,0,0);
849 	RecalcBorder();
850 }
851 //-----------------------------------------------------------------------------
CRange(HCDT a,bool add,mreal fact)852 void mglBase::CRange(HCDT a,bool add, mreal fact)
853 {
854 	mreal v1=a->Minimal(), v2=a->Maximal(), dv;
855 	dv=(v2-v1)*fact;	v1 -= dv;	v2 += dv;
856 	CRange(v1,v2,add);
857 }
CRange(mreal v1,mreal v2,bool add)858 void mglBase::CRange(mreal v1,mreal v2,bool add)
859 {
860 	if(!mgl_isrange(v1,v2) && !add)	return;
861 	if(!add)
862 	{
863 		if(mgl_isnum(v1))	Min.c = v1;
864 		if(mgl_isnum(v2))	Max.c = v2;
865 	}
866 	else if(Min.c<Max.c)
867 	{
868 		if(Min.c>v1)	Min.c=v1;
869 		if(Max.c<v2)	Max.c=v2;
870 	}
871 	else
872 	{
873 		mreal dv = Min.c;
874 		Min.c = v1<Max.c ? v1:Max.c;
875 		Max.c = v2>dv ? v2:dv;
876 	}
877 	if(Org.c<Min.c && mgl_isnum(Org.c))	Org.c = Min.c;
878 	if(Org.c>Max.c && mgl_isnum(Org.c))	Org.c = Max.c;
879 	if((TernAxis&3)==0)
880 	{
881 		OMax.c = Max.c;	OMin.c = Min.c;
882 		mglScaleAxis(Min.c, Max.c, Org.c, AMin.c, AMax.c);
883 	}
884 	RecalcCRange();
885 }
886 //-----------------------------------------------------------------------------
XRange(HCDT a,bool add,mreal fact)887 void mglBase::XRange(HCDT a,bool add,mreal fact)
888 {
889 	mreal v1=a->Minimal(), v2=a->Maximal(), dv;
890 	dv=(v2-v1)*fact;	v1 -= dv;	v2 += dv;
891 	XRange(v1,v2,add);
892 }
XRange(mreal v1,mreal v2,bool add)893 void mglBase::XRange(mreal v1,mreal v2,bool add)
894 {
895 	if(!mgl_isrange(v1,v2) && !add)	return;
896 	if(!add)
897 	{
898 		if(mgl_isnum(v1))	Min.x = v1;
899 		if(mgl_isnum(v2))	Max.x = v2;
900 	}
901 	else if(Min.x<Max.x)
902 	{
903 		if(Min.x>v1)	Min.x=v1;
904 		if(Max.x<v2)	Max.x=v2;
905 	}
906 	else
907 	{
908 		mreal dv = Min.x;
909 		Min.x = v1<Max.x ? v1:Max.x;
910 		Max.x = v2>dv ? v2:dv;
911 	}
912 	if(Org.x<Min.x && mgl_isnum(Org.x))	Org.x = Min.x;
913 	if(Org.x>Max.x && mgl_isnum(Org.x))	Org.x = Max.x;
914 	if((TernAxis&3)==0)
915 	{
916 		OMax.x = Max.x;	OMin.x = Min.x;
917 		mglScaleAxis(Min.x, Max.x, Org.x, AMin.x, AMax.x);
918 	}
919 	RecalcBorder();
920 }
921 //-----------------------------------------------------------------------------
YRange(HCDT a,bool add,mreal fact)922 void mglBase::YRange(HCDT a,bool add,mreal fact)
923 {
924 	mreal v1=a->Minimal(), v2=a->Maximal(), dv;
925 	dv=(v2-v1)*fact;	v1 -= dv;	v2 += dv;
926 	YRange(v1,v2,add);
927 }
YRange(mreal v1,mreal v2,bool add)928 void mglBase::YRange(mreal v1,mreal v2,bool add)
929 {
930 	if(!mgl_isrange(v1,v2) && !add)	return;
931 	if(!add)
932 	{
933 		if(mgl_isnum(v1))	Min.y = v1;
934 		if(mgl_isnum(v2))	Max.y = v2;
935 	}
936 	else if(Min.y<Max.y)
937 	{
938 		if(Min.y>v1)	Min.y=v1;
939 		if(Max.y<v2)	Max.y=v2;
940 	}
941 	else
942 	{
943 		mreal dv = Min.y;
944 		Min.y = v1<Max.y ? v1:Max.y;
945 		Max.y = v2>dv ? v2:dv;
946 	}
947 	if(Org.y<Min.y && mgl_isnum(Org.y))	Org.y = Min.y;
948 	if(Org.y>Max.y && mgl_isnum(Org.y))	Org.y = Max.y;
949 	if((TernAxis&3)==0)
950 	{
951 		OMax.y = Max.y;	OMin.y = Min.y;
952 		mglScaleAxis(Min.y, Max.y, Org.y, AMin.y, AMax.y);
953 
954 	}
955 	RecalcBorder();
956 }
957 //-----------------------------------------------------------------------------
ZRange(HCDT a,bool add,mreal fact)958 void mglBase::ZRange(HCDT a,bool add,mreal fact)
959 {
960 	mreal v1=a->Minimal(), v2=a->Maximal(), dv;
961 	dv=(v2-v1)*fact;	v1 -= dv;	v2 += dv;
962 	ZRange(v1,v2,add);
963 }
ZRange(mreal v1,mreal v2,bool add)964 void mglBase::ZRange(mreal v1,mreal v2,bool add)
965 {
966 	if(!mgl_isrange(v1,v2) && !add)	return;
967 	if(!add)
968 	{
969 		if(mgl_isnum(v1))	Min.z = v1;
970 		if(mgl_isnum(v2))	Max.z = v2;
971 	}
972 	else if(Min.z<Max.z)
973 	{
974 		if(Min.z>v1)	Min.z=v1;
975 		if(Max.z<v2)	Max.z=v2;
976 	}
977 	else
978 	{
979 		mreal dv = Min.z;
980 		Min.z = v1<Max.z ? v1:Max.z;
981 		Max.z = v2>dv ? v2:dv;
982 	}
983 	if(Org.z<Min.z && mgl_isnum(Org.z))	Org.z = Min.z;
984 	if(Org.z>Max.z && mgl_isnum(Org.z))	Org.z = Max.z;
985 	if((TernAxis&3)==0)
986 	{
987 		OMax.z = Max.z;	OMin.z = Min.z;
988 		mglScaleAxis(Min.z, Max.z, Org.z, AMin.z, AMax.z);
989 	}
990 	RecalcBorder();
991 }
992 //-----------------------------------------------------------------------------
SetAutoRanges(mreal x1,mreal x2,mreal y1,mreal y2,mreal z1,mreal z2,mreal c1,mreal c2)993 void mglBase::SetAutoRanges(mreal x1, mreal x2, mreal y1, mreal y2, mreal z1, mreal z2, mreal c1, mreal c2)
994 {
995 	if(mgl_isrange(x1,x2))	{	Min.x = x1;	Max.x = x2;	}
996 	if(mgl_isrange(y1,y2))	{	Min.y = y1;	Max.y = y2;	}
997 	if(mgl_isrange(z1,z2))	{	Min.z = z1;	Max.z = z2;	}
998 	if(mgl_isrange(c1,c2))	{	Min.c = c1;	Max.c = c2;	}
999 }
1000 //-----------------------------------------------------------------------------
Ternary(int t)1001 void mglBase::Ternary(int t)
1002 {
1003 	static mglPoint x1(-1,-1,-1),x2(1,1,1),o(NAN,NAN,NAN);
1004 	static bool c = true;
1005 	TernAxis = t;
1006 	if(t&3)
1007 	{
1008 		if(c)	{	x1 = Min;	x2 = Max;	o = Org;	}
1009 		SetRanges(mglPoint(0,0,0),mglPoint(1,1,(t&3)==1?0:1));
1010 		Org.Set(0,0,(t&3)==1?NAN:0);	c = false;
1011 	}
1012 	else if(!c)	{	SetRanges(x1,x2);	Org=o;	c=true;	}
1013 }
1014 //-----------------------------------------------------------------------------
1015 //		Transformation functions
1016 //-----------------------------------------------------------------------------
SetFunc(const char * EqX,const char * EqY,const char * EqZ,const char * EqA)1017 void mglBase::SetFunc(const char *EqX,const char *EqY,const char *EqZ,const char *EqA)
1018 {
1019 	if(fa)	delete fa;
1020 	if(fx)	delete fx;
1021 	if(fy)	delete fy;
1022 	if(fz)	delete fz;
1023 	if(EqX && *EqX && (EqX[0]!='x' || EqX[1]!=0))
1024 		fx = new mglFormula(EqX);
1025 	else	fx = 0;
1026 	if(EqY && *EqY && (EqY[0]!='y' || EqY[1]!=0))
1027 		fy = new mglFormula(EqY);
1028 	else	fy = 0;
1029 	if(EqZ && *EqZ && (EqZ[0]!='z' || EqZ[1]!=0))
1030 		fz = new mglFormula(EqZ);
1031 	else	fz = 0;
1032 	if(EqA && *EqA && ((EqA[0]!='c' && EqA[0]!='a') || EqA[1]!=0))
1033 		fa = new mglFormula(EqA);
1034 	else	fa = 0;
1035 	RecalcBorder();
1036 }
1037 //-----------------------------------------------------------------------------
CutOff(const char * EqC)1038 void mglBase::CutOff(const char *EqC)
1039 {
1040 #pragma omp critical(eq)
1041 	{
1042 		if(fc)	delete fc;
1043 		fc = (EqC && EqC[0])?new mglFormula(EqC):0;
1044 	}
1045 }
1046 //-----------------------------------------------------------------------------
SetCoor(int how)1047 void mglBase::SetCoor(int how)
1048 {
1049 	switch(how)
1050 	{
1051 	case mglCartesian:	SetFunc(0,0);	break;
1052 	case mglPolar:
1053 		SetFunc("x*cos(y)","x*sin(y)");	break;
1054 	case mglSpherical:
1055 		SetFunc("x*sin(y)*cos(z)","x*sin(y)*sin(z)","x*cos(y)");	break;
1056 	case mglParabolic:
1057 		SetFunc("x*y","(x*x-y*y)/2");	break;
1058 	case mglParaboloidal:
1059 		SetFunc("(x*x-y*y)*cos(z)/2","(x*x-y*y)*sin(z)/2","x*y");	break;
1060 	case mglOblate:
1061 		SetFunc("cosh(x)*cos(y)*cos(z)","cosh(x)*cos(y)*sin(z)","sinh(x)*sin(y)");	break;
1062 //		SetFunc("x*y*cos(z)","x*y*sin(z)","(x*x-1)*(1-y*y)");	break;
1063 	case mglProlate:
1064 		SetFunc("sinh(x)*sin(y)*cos(z)","sinh(x)*sin(y)*sin(z)","cosh(x)*cos(y)");	break;
1065 	case mglElliptic:
1066 		SetFunc("cosh(x)*cos(y)","sinh(x)*sin(y)");	break;
1067 	case mglToroidal:
1068 		SetFunc("sinh(x)*cos(z)/(cosh(x)-cos(y))","sinh(x)*sin(z)/(cosh(x)-cos(y))",
1069 			"sin(y)/(cosh(x)-cos(y))");	break;
1070 	case mglBispherical:
1071 		SetFunc("sin(y)*cos(z)/(cosh(x)-cos(y))","sin(y)*sin(z)/(cosh(x)-cos(y))",
1072 			"sinh(x)/(cosh(x)-cos(y))");	break;
1073 	case mglBipolar:
1074 		SetFunc("sinh(x)/(cosh(x)-cos(y))","sin(y)/(cosh(x)-cos(y))");	break;
1075 	case mglLogLog:	SetFunc("lg(x)","lg(y)");	break;
1076 	case mglLogX:	SetFunc("lg(x)","");	break;
1077 	case mglLogY:	SetFunc("","lg(y)");	break;
1078 	default:	SetFunc(0,0);	break;
1079 	}
1080 }
1081 //-----------------------------------------------------------------------------
ClearEq()1082 void mglBase::ClearEq()
1083 {
1084 #pragma omp critical(eq)
1085 	{
1086 		if(fx)	delete fx;
1087 		if(fy)	delete fy;
1088 		if(fz)	delete fz;
1089 		if(fa)	delete fa;
1090 		if(fc)	delete fc;
1091 		fx = fy = fz = fc = fa = 0;
1092 	}
1093 	RecalcBorder();
1094 }
1095 //-----------------------------------------------------------------------------
1096 //		Colors ids
1097 //-----------------------------------------------------------------------------
1098 MGL_EXPORT mglColorID mglColorIds[31] = {{'k', mglColor(0,0,0)},
1099 	{'r', mglColor(1,0,0)},		{'R', mglColor(0.5,0,0)},
1100 	{'g', mglColor(0,1,0)},		{'G', mglColor(0,0.5,0)},
1101 	{'b', mglColor(0,0,1)},		{'B', mglColor(0,0,0.5)},
1102 	{'w', mglColor(1,1,1)},		{'W', mglColor(0.7,0.7,0.7)},
1103 	{'c', mglColor(0,1,1)},		{'C', mglColor(0,0.5,0.5)},
1104 	{'m', mglColor(1,0,1)},		{'M', mglColor(0.5,0,0.5)},
1105 	{'y', mglColor(1,1,0)},		{'Y', mglColor(0.5,0.5,0)},
1106 	{'h', mglColor(0.5,0.5,0.5)},	{'H', mglColor(0.3,0.3,0.3)},
1107 	{'l', mglColor(0,1,0.5)},	{'L', mglColor(0,0.5,0.25)},
1108 	{'e', mglColor(0.5,1,0)},	{'E', mglColor(0.25,0.5,0)},
1109 	{'n', mglColor(0,0.5,1)},	{'N', mglColor(0,0.25,0.5)},
1110 	{'u', mglColor(0.5,0,1)},	{'U', mglColor(0.25,0,0.5)},
1111 	{'q', mglColor(1,0.5,0)},	{'Q', mglColor(0.5,0.25,0)},
1112 	{'p', mglColor(1,0,0.5)},	{'P', mglColor(0.5,0,0.25)},
1113 	{' ', mglColor(-1,-1,-1)},	{0, mglColor(-1,-1,-1)}	// the last one MUST have id=0
1114 };
1115 //-----------------------------------------------------------------------------
mgl_chrrgb(char p,float c[3])1116 void MGL_EXPORT mgl_chrrgb(char p, float c[3])
1117 {
1118 	c[0]=c[1]=c[2]=-1;
1119 	for(long i=0; mglColorIds[i].id; i++)
1120 		if(mglColorIds[i].id==p)
1121 		{
1122 			c[0]=mglColorIds[i].col.r;
1123 			c[1]=mglColorIds[i].col.g;
1124 			c[2]=mglColorIds[i].col.b;
1125 			break;
1126 		}
1127 }
1128 //-----------------------------------------------------------------------------
mgl_get_num_color(const char * s,int smooth)1129 size_t MGL_EXPORT_PURE mgl_get_num_color(const char *s, int smooth)
1130 {
1131 	if(!s || !s[0])	return 0;
1132 	size_t l=strlen(s), n=0;	long j=0;
1133 	for(size_t i=0;i<l;i++)		// find number of colors
1134 	{
1135 		if(smooth>=0 && s[i]==':' && j<1)	break;
1136 		if(s[i]=='{' && strchr(MGL_COLORS"x",s[i+1]) && j<1)	n++;
1137 		if(s[i]=='[' || s[i]=='{')	j++;
1138 		if(s[i]==']' || s[i]=='}')	j--;
1139 		if(strchr(MGL_COLORS,s[i]) && j<1)	n++;
1140 //		if(smooth && s[i]==':')	break;	// NOTE: should use []
1141 	}
1142 	return n;
1143 }
1144 //-----------------------------------------------------------------------------
Set(const char * s,int smooth,mreal alpha)1145 void mglTexture::Set(const char *s, int smooth, mreal alpha)
1146 {
1147 	// NOTE: New syntax -- colors are CCCCC or {CNCNCCCN}; options inside []
1148 	if(!s || !s[0])	return;
1149 	mgl_strncpy(Sch,s,259);	Smooth=smooth;	Alpha=alpha;	Clear();
1150 
1151 	long l=strlen(s);
1152 	bool map = smooth==2 || mglchr(s,'%'), sm = smooth>=0 && !strchr(s,'|');	// Use mapping, smoothed colors
1153 	n = mgl_get_num_color(s,smooth);
1154 	if(!n)
1155 	{
1156 		if(strchr(s,'|') && !smooth)	// sharp colors
1157 		{	n=l=6;	s=MGL_DEF_SCH;	sm = false;	}
1158 		else if(smooth==0)		// none colors but color scheme
1159 		{	n=l=6;	s=MGL_DEF_SCH;	}
1160 	}
1161 	if(n<=0)	return;
1162 	bool man=sm;
1163 	c0 = new mglColor[2*n];	// Colors itself
1164 	val = new float[n];
1165 	for(long i=0, m=0, j=n=0;i<l;i++)	// fill colors
1166 	{
1167 		if(smooth>=0 && s[i]==':' && j<1)	break;
1168 		if(s[i]=='[')	j++;
1169 		if(s[i]==']')	j--;
1170 		if(s[i]=='{')	m++;
1171 		if(s[i]=='}')	m--;
1172 		if(strchr(MGL_COLORS,s[i]) && j<1 && (m==0 || s[i-1]=='{'))	// {CN,val} format, where val in [0,1]
1173 		{
1174 			if(m>0 && s[i+1]>'0' && s[i+1]<='9')// ext color
1175 			{	c0[2*n].Set(s[i],(s[i+1]-'0')/5.f);	i++;	}
1176 			else	c0[2*n].Set(s[i]);	// usual color
1177 			val[n]=-1;	c0[2*n].a = -1;	n++;
1178 		}
1179 		if(s[i]=='x' && i>0 && s[i-1]=='{' && j<1)	// {xRRGGBB,val} format, where val in [0,1]
1180 		{
1181 			uint32_t id = strtoul(s+1+i,0,16);
1182 			if(memchr(s+i+1,'}',8) || memchr(s+i+1,',',8))	c0[2*n].a = -1;
1183 			else	{	c0[2*n].a = (id%256)/255.;	id /= 256;	}
1184 			c0[2*n].b = (id%256)/255.;	id /= 256;
1185 			c0[2*n].g = (id%256)/255.;	id /= 256;
1186 			c0[2*n].r = (id%256)/255.;
1187 			while(strchr("0123456789abcdefABCDEFx",s[i]))	i++;
1188 			val[n]=-1;	n++;	i--;
1189 		}
1190 		if(s[i]==',' && m>0 && j<1 && n>0)
1191 			val[n-1] = atof(s+i+1);
1192 		// NOTE: User can change alpha if it placed like {AN}
1193 		if(s[i]=='A' && j<1 && m>0 && s[i+1]>'0' && s[i+1]<='9')
1194 		{	man=false;	alpha = 0.1*(s[i+1]-'0');	i++;	}
1195 	}
1196 	for(long i=0;i<n;i++)	// default texture
1197 	{
1198 		if(c0[2*i].a<0)	c0[2*i].a=alpha;
1199 		c0[2*i+1]=c0[2*i];
1200 		if(man)	c0[2*i].a=0;
1201 	}
1202 	if(map && sm && n>1)		// map texture
1203 	{
1204 		if(n==2)
1205 		{	c0[1]=c0[2];	c0[2]=c0[0];	c0[0]=BC;	c0[3]=c0[1]+c0[2];	}
1206 		else if(n==3)
1207 		{	c0[1]=c0[2];	c0[2]=c0[0];	c0[0]=BC;	c0[3]=c0[4];	n=2;}
1208 		else
1209 		{	c0[1]=c0[4];	c0[3]=c0[6];	n=2;	}
1210 		c0[0].a = c0[1].a = c0[2].a = c0[3].a = alpha;
1211 		val[0]=val[1]=-1;
1212 	}
1213 	// TODO if(!sm && n==1)	then try to find color in palette ???
1214 
1215 	// fill missed values  of val[]
1216 	float  v1=0,v2=1;
1217 	std::vector <long>  def;
1218 	val[0]=0;	val[n-1]=1;	// boundary have to be [0,1]
1219 	for(long i=0;i<n;i++) if(val[i]>0 && val[i]<1) 	def.push_back(i);
1220 	def.push_back(n-1);
1221 	long i1=0;
1222 	for(size_t j=0;j<def.size();j++)	for(long i=i1+1;i<def[j];i++)
1223 	{
1224 		i1 = j>0?def[j-1]:0;	long i2 = def[j];
1225 		v1 = val[i1];	v2 = val[i2];
1226 		v2 = i2-i1>1?(v2-v1)/(i2-i1):0;
1227 		val[i]=v1+v2*(i-i1);
1228 	}
1229 	// fill texture itself
1230 	mreal v=sm?(n-1)/255.:n/256.;
1231 	if(!sm)	for(long i=0;i<256;i++)
1232 	{
1233 		long j = 2*long(v*i);	//u-=j;
1234 		col[2*i] = c0[j];	col[2*i+1] = c0[j+1];
1235 	}
1236 	else	for(long i=i1=0;i<256;i++)
1237 	{
1238 		mreal u = v*i;	long j = long(u);	//u-=j;
1239 		if(j<n-1)	// advanced scheme using val
1240 		{
1241 			for(;i1<n-1 && i>=255*val[i1];i1++);
1242 			v2 = i1<n?1/(val[i1]-val[i1-1]):0;
1243 			j=i1-1;	u=(i/255.-val[j])*v2;
1244 			col[2*i] = c0[2*j]*(1-u)+c0[2*j+2]*u;
1245 			col[2*i+1]=c0[2*j+1]*(1-u)+c0[2*j+3]*u;
1246 		}
1247 		else
1248 		{	col[2*i] = c0[2*n-2];col[2*i+1] = c0[2*n-1];	}
1249 	}
1250 }
1251 //-----------------------------------------------------------------------------
GetC(mreal u,mreal v) const1252 mglColor mglTexture::GetC(mreal u,mreal v) const
1253 {
1254 	u -= long(u);
1255 	long i=long(255*u);	u = u*255-i;
1256 	const mglColor *s=col+2*i;
1257 	return (s[0]*(1-u)+s[2]*u)*(1-v) + (s[1]*(1-u)+s[3]*u)*v;
1258 }
1259 //-----------------------------------------------------------------------------
GetC(mreal u,mreal v,mglPnt & p) const1260 void mglTexture::GetC(mreal u,mreal v,mglPnt &p) const
1261 {
1262 	u -= long(u);
1263 	long i=long(255*u);	u = u*255-i;
1264 	const mglColor &s0=col[2*i], &s1=col[2*i+1], &s2=col[2*i+2], &s3=col[2*i+3];
1265 	p.r = (s0.r*(1-u)+s2.r*u)*(1-v) + (s1.r*(1-u)+s3.r*u)*v;
1266 	p.g = (s0.g*(1-u)+s2.g*u)*(1-v) + (s1.g*(1-u)+s3.g*u)*v;
1267 	p.b = (s0.b*(1-u)+s2.b*u)*(1-v) + (s1.b*(1-u)+s3.b*u)*v;
1268 	p.a = (s0.a*(1-u)+s2.a*u)*(1-v) + (s1.a*(1-u)+s3.a*u)*v;
1269 //	p.a = (s0.a*(1-u)+s2.a*u)*v + (s1.a*(1-u)+s3.a*u)*(1-v);	// for alpha use inverted
1270 }
1271 //-----------------------------------------------------------------------------
AddTexture(const char * cols,int smooth)1272 long mglBase::AddTexture(const char *cols, int smooth)
1273 {
1274 	if(smooth>=0)	SetMask(cols);
1275 	mglTexture t(cols,smooth,smooth==2?AlphaDef:1);
1276 	if(t.n==0)	return smooth<0 ? 0:1;
1277 	if(smooth<0)	CurrPal=0;
1278 	// check if already exist
1279 	for(size_t i=0;i<Txt.size();i++)
1280 		if(!t.IsSame(Txt[i]))	continue;	else	return i;
1281 	// create new one
1282 	long k;
1283 #pragma omp critical(txt)
1284 	{k=Txt.size();	MGL_PUSH(Txt,t,mutexTxt);}	return k;
1285 }
1286 //-----------------------------------------------------------------------------
AddTexture(mglColor c)1287 mreal mglBase::AddTexture(mglColor c)
1288 {
1289 	if(!c.Valid())	return -1;
1290 	// first lets try an existed one
1291 	for(size_t i=0;i<Txt.size();i++)	for(int j=0;j<255;j++)
1292 		if(c==Txt[i].col[2*j])
1293 			return i+j/255.;
1294 	// add new texture
1295 	mglTexture t;
1296 	for(long i=0;i<MGL_TEXTURE_COLOURS;i++)	t.col[i]=c;
1297 	long k;
1298 #pragma omp critical(txt)
1299 	{k=Txt.size();	MGL_PUSH(Txt,t,mutexTxt);}	return k;
1300 }
1301 //-----------------------------------------------------------------------------
1302 //		Coloring and palette
1303 //-----------------------------------------------------------------------------
NextColor(long & id)1304 mreal mglBase::NextColor(long &id)
1305 {
1306 	long i=labs(id)/256, n=Txt[i].n, p=labs(id)&0xff;
1307 	if(id>=0)	{	p=(p+1)%n;	id = 256*i+p;	}
1308 	CDef = i + (n>0 ? (p+0.5)/n : 0);	CurrPal++;
1309 	sprintf(last_style+11,"{&%g}",CDef);
1310 	if(!leg_str.empty())
1311 	{	AddLegend(leg_str.c_str(),last_style);	leg_str.clear();	}
1312 	return CDef;
1313 }
1314 //-----------------------------------------------------------------------------
NextColor(long id,long sh)1315 mreal mglBase::NextColor(long id, long sh)
1316 {
1317 	long i=labs(id)/256, n=Txt[i].n, p=labs(id)&0xff;
1318 	if(id>=0)	p=(p+sh)%n;
1319 	mreal cc = i + (n>0 ? (p+0.5)/n : 0);
1320 	sprintf(last_style+11,"{&%g}",cc);
1321 	return cc;
1322 }
1323 //-----------------------------------------------------------------------------
mglchrs(const char * str,const char * chr)1324 MGL_EXPORT_PURE const char *mglchrs(const char *str, const char *chr)
1325 {
1326 	if(!str || !str[0] || !chr || !chr[0])	return NULL;
1327 	size_t l=strlen(chr);
1328 	for(size_t i=0;i<l;i++)
1329 	{
1330 		const char *res = mglchr(str,chr[i]);
1331 		if(res)	return res;
1332 	}
1333 	return NULL;
1334 }
1335 //-----------------------------------------------------------------------------
mglchr(const char * str,char ch)1336 MGL_EXPORT_PURE const char *mglchr(const char *str, char ch)
1337 {
1338 	if(!str || !str[0])	return NULL;
1339 	size_t l=strlen(str),k=0;
1340 	for(size_t i=0;i<l;i++)
1341 	{
1342 		char c = str[i];
1343 		if(c=='{')	k++;
1344 		if(c=='}')	k--;
1345 		if(c==ch && k==0)	return str+i;
1346 	}
1347 	return NULL;
1348 }
1349 //-----------------------------------------------------------------------------
SetPenPal(const char * p,long * Id,bool pal)1350 char mglBase::SetPenPal(const char *p, long *Id, bool pal)
1351 {
1352 	char mk=0;
1353 	PDef = 0xffff;	// reset to solid line
1354 	strcpy(last_style,"__1 {dFFFF}k\0");
1355 
1356 	const char *s;
1357 	Arrow1 = Arrow2 = 0;	PenWidth = 1;
1358 	if(p && *p)
1359 	{
1360 //		const char *col = "wkrgbcymhRGBCYMHWlenuqpLENUQP";
1361 		unsigned val[8] = {0x0000, 0xffff, 0x00ff, 0x0f0f, 0x1111, 0x087f, 0x2727, 0x3333};
1362 		const char *stl = " -|;:ji=";
1363 		const char *mrk = "*o+xsd.^v<>";
1364 		const char *MRK = "YOPXSDCTVLR";
1365 		const char *wdh = "123456789";
1366 		const char *arr = "AKDTVISOX_";
1367 		long m=0;
1368 		size_t l=strlen(p);
1369 		for(size_t i=0;i<l;i++)
1370 		{
1371 			if(p[i]=='{')	m++;
1372 			if(p[i]=='}')	m--;
1373 			if(m>0 && p[i]=='d')	PDef = strtol(p+i+1,0,16);
1374 			if(m>0)	continue;
1375 			s = mglchr(stl,p[i]);
1376 			if(s)
1377 			{	PDef = val[s-stl];	sprintf(last_style+6,"%04x",PDef);	last_style[10]='}';	}
1378 			else if(mglchr(mrk,p[i]))
1379 			{	mk = p[i];	last_style[3] = mk;	}
1380 			else if(mglchr(wdh,p[i]))
1381 			{	PenWidth = p[i]-'0';	last_style[2] = p[i];	}
1382 			else if(mglchr(arr,p[i]))
1383 			{
1384 				if(!Arrow2)	Arrow2 = p[i];
1385 				else	Arrow1 = p[i];
1386 			}
1387 		}
1388 		if(!Arrow1)	Arrow1='_';
1389 		if(!Arrow2)	Arrow2='_';
1390 		if(mglchr(p,'#'))
1391 		{
1392 			s = mglchr(mrk,mk);
1393 			if(s)
1394 			{	mk = MRK[s-mrk];	last_style[3] = mk;	}
1395 		}
1396 		if((s=strstr(p,"{&"))!=0)
1397 		{	mk = last_style[3] = p[3];	strcpy(last_style+11,s);	}
1398 		else if(mk && mglchr(p,'&'))
1399 		{	mk += 128;	last_style[3] = mk;	}
1400 		last_style[0] = Arrow1;	last_style[1] = Arrow2;
1401 	}
1402 	if(pal)
1403 	{
1404 		if(p && (s=strstr(p,"{&"))!=0)
1405 		{
1406 			CDef = atof(s+2);
1407 //			if(Id)	*Id=long(tt)*256+(n+CurrPal-1)%n;
1408 		}
1409 		else
1410 		{
1411 			long tt, n;
1412 			tt = AddTexture(p?p:MGL_DEF_PAL,-1);	n=Txt[tt].n;
1413 			CDef = tt+((n+CurrPal-1)%n+0.5)/n;
1414 			if(Id)	*Id=long(tt)*256+(n+CurrPal-1)%n;
1415 			sprintf(last_style+11,"{&%g}",CDef);
1416 		}
1417 	}
1418 	if(Arrow1=='_')	Arrow1=0;
1419 	if(Arrow2=='_')	Arrow2=0;
1420 	return mk;
1421 }
1422 //-----------------------------------------------------------------------------
1423 // keep this for restore default mask
1424 MGL_EXPORT uint64_t mgl_mask_def[16]={
1425 	0x000000FF00000000,	0x080808FF08080808,	0x0000FF00FF000000,	0x0000000F00000000,
1426 	0x0000182424180000,	0x0000183C3C180000,	0x00003C24243C0000,	0x00003C3C3C3C0000,
1427 	0x0000060990600000,	0x0060584658600000,	0x00061A621A060000,	0x0000002700000000,
1428 	0x0008083E08080000,	0x0139010010931000,	0x0000001818000000,	0x101010FF010101FF};
1429 MGL_EXPORT uint64_t mgl_mask_val[16]={
1430 	0x000000FF00000000,	0x080808FF08080808,	0x0000FF00FF000000,	0x0000000F00000000,
1431 	0x0000182424180000,	0x0000183C3C180000,	0x00003C24243C0000,	0x00003C3C3C3C0000,
1432 	0x0000060990600000,	0x0060584658600000,	0x00061A621A060000,	0x0000002700000000,
1433 	0x0008083E08080000,	0x0139010010931000,	0x0000001818000000,	0x101010FF010101FF};
1434 // 	0x000000FF00000000,	0x080808FF08080808,	0x0000FF00FF000000,	0x0000007700000000,
1435 // 	0x0000182424180000,	0x0000183C3C180000,	0x00003C24243C0000,	0x00003C3C3C3C0000,
1436 // 	0x0000060990600000,	0x0060584658600000,	0x00061A621A060000,	0x0000005F00000000,
1437 //	0x0008142214080000,	0x00081C3E1C080000,	0x8142241818244281,	0x0000001824420000};
SetMask(const char * p)1438 void mglBase::SetMask(const char *p)
1439 {
1440 	mask = MGL_SOLID_MASK;	// reset to solid face
1441 	PenWidth = 1;	MaskAn=DefMaskAn;
1442 	if(p && *p)
1443 	{
1444 		const char *msk = MGL_MASK_ID, *s;
1445 		const char *wdh = "123456789";
1446 		long m=0, l=strlen(p);
1447 		for(long i=0;i<l;i++)
1448 		{
1449 			if(p[i]=='{')	m++;
1450 			if(p[i]=='}')	m--;
1451 			if(m>0 && p[i]=='s')	mask = strtoull(p+i+1,0,16);
1452 			if(m>0)	continue;
1453 			if(p[i]==':')	break;
1454 			s = mglchr(msk, p[i]);
1455 			if(s)	mask = mgl_mask_val[s-msk];
1456 			else if(mglchr(wdh,p[i]))	PenWidth = p[i]-'0';
1457 			else if(p[i]=='I')	MaskAn=90;
1458 			else if(p[i]=='/')	MaskAn=315;	// =360-45
1459 			else if(p[i]=='\\')	MaskAn=45;
1460 		}
1461 		// use line if rotation only specified
1462 		if(mask==MGL_SOLID_MASK && MaskAn!=0)	mask = mgl_mask_val[0];
1463 	}
1464 }
1465 //-----------------------------------------------------------------------------
GetA(mreal a) const1466 mreal mglBase::GetA(mreal a) const
1467 {
1468 	if(fa)	a = fa->Calc(0,0,0,a);
1469 	a = (a-FMin.c)/(FMax.c-FMin.c);
1470 	a = (a>1?1:(a<0?0:a))/MGL_FEPSILON;	// for texture a must be <1 always!!!
1471 //	a = (a<1?(a>0?a:0):1)/MGL_FEPSILON;	// for texture a must be <1 always!!!
1472 	return a;
1473 }
1474 //-----------------------------------------------------------------------------
GetX(HCDT x,int i,int j,int k)1475 mglPoint GetX(HCDT x, int i, int j, int k)
1476 {
1477 	k = k<x->GetNz() ? k : 0;
1478 	if(x->GetNy()>1)
1479 		return mglPoint(x->v(i,j,k),x->dvx(i,j,k),x->dvy(i,j,k));
1480 	else
1481 		return mglPoint(x->v(i),x->dvx(i),0);
1482 }
1483 //-----------------------------------------------------------------------------
GetY(HCDT y,int i,int j,int k)1484 mglPoint GetY(HCDT y, int i, int j, int k)
1485 {
1486 	k = k<y->GetNz() ? k : 0;
1487 	if(y->GetNy()>1)
1488 		return mglPoint(y->v(i,j,k),y->dvx(i,j,k),y->dvy(i,j,k));
1489 	else
1490 		return mglPoint(y->v(j),0,y->dvx(j));
1491 }
1492 //-----------------------------------------------------------------------------
GetZ(HCDT z,int i,int j,int k)1493 mglPoint GetZ(HCDT z, int i, int j, int k)
1494 {
1495 	if(z->GetNy()>1)
1496 		return mglPoint(z->v(i,j,k),z->dvx(i,j,k),z->dvy(i,j,k));
1497 	else
1498 		return mglPoint(z->v(k),0,0);
1499 }
1500 //-----------------------------------------------------------------------------
vect_plot(long p1,long p2,mreal s)1501 void mglBase::vect_plot(long p1, long p2, mreal s)
1502 {
1503 	if(p1<0 || p2<0)	return;
1504 	const mglPnt &q1=Pnt[p1], &q2=Pnt[p2];
1505 	mglPnt s1=q2,s2=q2;
1506 	s = s<=0 ? 0.1 : s*0.1;
1507 	s1.x=s1.xx = q2.x - 3*s*(q2.x-q1.x) + s*(q2.y-q1.y);
1508 	s2.x=s2.xx = q2.x - 3*s*(q2.x-q1.x) - s*(q2.y-q1.y);
1509 	s1.y=s1.yy = q2.y - 3*s*(q2.y-q1.y) - s*(q2.x-q1.x);
1510 	s2.y=s2.yy = q2.y - 3*s*(q2.y-q1.y) + s*(q2.x-q1.x);
1511 	s1.z=s1.zz=s2.z=s2.zz = q2.z - 3*s*(q2.z-q1.z);
1512 	long n1,n2;
1513 #pragma omp critical(pnt)
1514 	{
1515 		n1=Pnt.size();	MGL_PUSH(Pnt,s1,mutexPnt);
1516 		n2=Pnt.size();	MGL_PUSH(Pnt,s2,mutexPnt);
1517 	}
1518 	line_plot(p1,p2);	line_plot(n1,p2);	line_plot(p2,n2);
1519 }
1520 //-----------------------------------------------------------------------------
mglFindArg(const char * str)1521 int MGL_LOCAL_PURE mglFindArg(const char *str)
1522 {
1523 	long l=0,k=0,len=strlen(str);
1524 	for(long i=0;i<len;i++)
1525 	{
1526 		if(str[i]=='\'') l++;
1527 		if(str[i]=='{') k++;
1528 		if(str[i]=='}') k--;
1529 		if(l%2==0 && k==0)
1530 		{
1531 			if(str[i]=='#' || str[i]==';')	return -i;
1532 			if(str[i]<=' ')	return i;
1533 		}
1534 	}
1535 	return 0;
1536 }
1537 //-----------------------------------------------------------------------------
SetAmbient(mreal bright)1538 void mglBase::SetAmbient(mreal bright)
1539 {	AmbBr=bright;	size_t n=Sub.size();	if(n>0)	Sub[n-1].AmbBr=bright;	}
SetDiffuse(mreal bright)1540 void mglBase::SetDiffuse(mreal bright)
1541 {	DifBr=bright;	size_t n=Sub.size();	if(n>0)	Sub[n-1].DifBr=bright;	}
1542 //-----------------------------------------------------------------------------
SaveState(const char * opt)1543 mreal mglBase::SaveState(const char *opt)
1544 {
1545 	if(saved)	return PrevState;
1546 	if(!opt || !opt[0])	return NAN;
1547 	MSS=MarkSize;	ASS=ArrowSize;
1548 	FSS=FontSize;	ADS=AlphaDef;
1549 	MNS=MeshNum;	CSS=Flag;	LSS=AmbBr;
1550 	MinS=Min;		MaxS=Max;	saved=true;
1551 	// parse option
1552 	char *qi=mgl_strdup(opt),*q=qi, *s;
1553 	mgl_strtrim(q);
1554 	// NOTE: not consider '#' inside legend entry !!!
1555 	s=strchr(q,'#');	if(s)	*s=0;
1556 	mreal res=NAN;
1557 	while(q && *q)
1558 	{
1559 		s=q;	q=strchr(s,';');
1560 		if(q)	{	*q=0;	q++;	}
1561 		mgl_strtrim(s);		char *a=s;
1562 		long n=mglFindArg(s);	if(n>0)	{	s[n]=0;	s=s+n+1;	}
1563 		mgl_strtrim(a);		char *b=s;
1564 		n=mglFindArg(s);	if(n>0)	{	s[n]=0;		s=s+n+1;	}
1565 		mgl_strtrim(b);
1566 
1567 		mreal ff=atof(b),ss;
1568 		size_opt = NAN;
1569 		if(!strcmp(b,"on"))	ff=1;
1570 		if(!strcmp(a+1,"range"))
1571 		{
1572 			n=mglFindArg(s);	char *c=s;
1573 			if(n>0)	{	s[n]=0;	s=s+n+1;	}
1574 			mgl_strtrim(c);		ss = atof(c);
1575 			if(a[0]=='x')		{	Min.x=ff;	Max.x=ss;	}
1576 			else if(a[0]=='y')	{	Min.y=ff;	Max.y=ss;	}
1577 			else if(a[0]=='z')	{	Min.z=ff;	Max.z=ss;	}
1578 // NOTE!	else if(a[0]=='c')	{	Min.c=ff;	Max.c=ss;	}	// Bad idea since there is formula for coloring
1579 		}
1580 		else if(!strcmp(a,"cut"))		SetCut(ff!=0);
1581 		else if(!strcmp(a,"meshnum"))	SetMeshNum(ff);
1582 		else if(!strcmp(a,"alpha"))		{Alpha(true);	SetAlphaDef(ff);}
1583 		else if(!strcmp(a,"light"))		Light(ff!=0);
1584 		else if(!strcmp(a,"ambient"))	SetAmbient(ff);
1585 		else if(!strcmp(a,"diffuse"))	SetDifLight(ff);
1586 		else if(!strcmp(a,"size"))
1587 		{	SetMarkSize(ff);	SetFontSize(ff);	SetArrowSize(ff);	size_opt=ff;	}
1588 		else if(!strcmp(a,"num") || !strcmp(a,"number") || !strcmp(a,"value"))	res=ff;
1589 		else if(!strcmp(a,"legend"))
1590 		{	if(*b=='\'')	{	b++;	b[strlen(b)-1]=0;	}	leg_str = b;	}
1591 	}
1592 	free(qi);	PrevState=res;	return res;
1593 }
1594 //-----------------------------------------------------------------------------
LoadState()1595 void mglBase::LoadState()
1596 {
1597 	if(!saved)	return;
1598 	MarkSize=MSS;	ArrowSize=ASS;
1599 	FontSize=FSS;	AlphaDef=ADS;
1600 	MeshNum=MNS;	Flag=CSS;	AmbBr=LSS;
1601 	Min=MinS;		Max=MaxS;	saved=false;
1602 }
1603 //-----------------------------------------------------------------------------
AddLegend(const wchar_t * text,const char * style)1604 void mglBase::AddLegend(const wchar_t *text,const char *style)
1605 {
1606 	if(text)
1607 #pragma omp critical(leg)
1608 		MGL_PUSH(Leg,mglText(text,style),mutexLeg);
1609 }
1610 //-----------------------------------------------------------------------------
AddLegend(const char * str,const char * style)1611 void mglBase::AddLegend(const char *str,const char *style)
1612 {
1613 	MGL_TO_WCS(str,AddLegend(wcs, style));
1614 }
1615 //-----------------------------------------------------------------------------
mgl_check_dim2(HMGL gr,HCDT x,HCDT y,HCDT z,HCDT a,const char * name,bool less)1616 bool MGL_EXPORT mgl_check_dim2(HMGL gr, HCDT x, HCDT y, HCDT z, HCDT a, const char *name, bool less)
1617 {
1618 //	if(!gr || !x || !y || !z)	return true;		// if data is absent then should be segfault!!!
1619 	long n=z->GetNx(),m=z->GetNy();
1620 	if(n<2 || m<2)	{	gr->SetWarn(mglWarnLow,name);	return true;	}
1621 	if(a && z->GetNN()!=a->GetNN())
1622 	{	gr->SetWarn(mglWarnDim,name);	return true;	}
1623 	if(less)
1624 	{
1625 		if(x->GetNx()<n)
1626 		{	gr->SetWarn(mglWarnDim,name);	return true;	}
1627 		if(y->GetNx()<m && (x->GetNy()<m || y->GetNx()<n || y->GetNy()<m))
1628 		{	gr->SetWarn(mglWarnDim,name);	return true;	}
1629 	}
1630 	else
1631 	{
1632 		if(x->GetNx()!=n)
1633 		{	gr->SetWarn(mglWarnDim,name);	return true;	}
1634 		if(y->GetNx()!=m && (x->GetNy()!=m || y->GetNx()!=n || y->GetNy()!=m))
1635 		{	gr->SetWarn(mglWarnDim,name);	return true;	}
1636 	}
1637 	return false;
1638 }
1639 //-----------------------------------------------------------------------------
mgl_check_dim0(HMGL gr,HCDT x,HCDT y,HCDT z,HCDT r,const char * name,bool less)1640 bool MGL_EXPORT mgl_check_dim0(HMGL gr, HCDT x, HCDT y, HCDT z, HCDT r, const char *name, bool less)
1641 {
1642 //	if(!gr || !x || !y)	return true;		// if data is absent then should be segfault!!!
1643 	long n=y->GetNx();
1644 	if(less)
1645 	{
1646 		if(x->GetNx()<n)		{	gr->SetWarn(mglWarnDim,name);	return true;	}
1647 		if(z && z->GetNx()<n)	{	gr->SetWarn(mglWarnDim,name);	return true;	}
1648 		if(r && r->GetNx()<n)	{	gr->SetWarn(mglWarnDim,name);	return true;	}
1649 	}
1650 	else
1651 	{
1652 		if(x->GetNx()!=n)		{	gr->SetWarn(mglWarnDim,name);	return true;	}
1653 		if(z && z->GetNx()!=n)	{	gr->SetWarn(mglWarnDim,name);	return true;	}
1654 		if(r && r->GetNx()!=n)	{	gr->SetWarn(mglWarnDim,name);	return true;	}
1655 	}
1656 	return false;
1657 }
1658 //-----------------------------------------------------------------------------
mgl_check_dim1(HMGL gr,HCDT x,HCDT y,HCDT z,HCDT r,const char * name,bool less)1659 bool MGL_EXPORT mgl_check_dim1(HMGL gr, HCDT x, HCDT y, HCDT z, HCDT r, const char *name, bool less)
1660 {
1661 //	if(!gr || !x || !y)	return true;		// if data is absent then should be segfault!!!
1662 	long n=y->GetNx();
1663 	if(n<2)	{	gr->SetWarn(mglWarnLow,name);	return true;	}
1664 	if(less)
1665 	{
1666 		if(x->GetNx()<n)		{	gr->SetWarn(mglWarnDim,name);	return true;	}
1667 		if(z && z->GetNx()<n)	{	gr->SetWarn(mglWarnDim,name);	return true;	}
1668 		if(r && r->GetNx()<n)	{	gr->SetWarn(mglWarnDim,name);	return true;	}
1669 	}
1670 	else
1671 	{
1672 		if(x->GetNx()!=n)		{	gr->SetWarn(mglWarnDim,name);	return true;	}
1673 		if(z && z->GetNx()!=n)	{	gr->SetWarn(mglWarnDim,name);	return true;	}
1674 		if(r && r->GetNx()!=n)	{	gr->SetWarn(mglWarnDim,name);	return true;	}
1675 	}
1676 	return false;
1677 }
1678 //-----------------------------------------------------------------------------
mgl_check_dim3(HMGL gr,bool both,HCDT x,HCDT y,HCDT z,HCDT a,HCDT b,const char * name)1679 bool MGL_EXPORT mgl_check_dim3(HMGL gr, bool both, HCDT x, HCDT y, HCDT z, HCDT a, HCDT b, const char *name)
1680 {
1681 // 	if(!gr || !x || !y || !z || !a)	return true;		// if data is absent then should be segfault!!!
1682 	long n=a->GetNx(),m=a->GetNy(),l=a->GetNz();
1683 	if(n<2 || m<2 || l<2)
1684 	{	gr->SetWarn(mglWarnLow,name);	return true;	}
1685 	if(!both && (x->GetNx()!=n || y->GetNx()!=m || z->GetNx()!=l))
1686 	{	gr->SetWarn(mglWarnDim,name);	return true;	}
1687 	if(b && b->GetNN()!=n*m*l)
1688 	{	gr->SetWarn(mglWarnDim,name);	return true;	}
1689 	return false;
1690 }
1691 //-----------------------------------------------------------------------------
mgl_check_trig(HMGL gr,HCDT nums,HCDT x,HCDT y,HCDT z,HCDT a,const char * name,int d)1692 bool MGL_EXPORT mgl_check_trig(HMGL gr, HCDT nums, HCDT x, HCDT y, HCDT z, HCDT a, const char *name, int d)
1693 {
1694 // 	if(!gr || !x || !y || !z || !a || !nums)	return true;		// if data is absent then should be segfault!!!
1695 	long n = x->GetNN(), m = nums->GetNy();
1696 	if(nums->GetNx()<d)	{	gr->SetWarn(mglWarnLow,name);	return true;	}
1697 	if(y->GetNN()!=n || z->GetNN()!=n)	{	gr->SetWarn(mglWarnDim,name);	return true;	}
1698 	if(a->GetNN()!=m && a->GetNN()!=n)	{	gr->SetWarn(mglWarnDim,name);	return true;	}
1699 	return false;
1700 }
1701 //-----------------------------------------------------------------------------
mgl_isnboth(HCDT x,HCDT y,HCDT z,HCDT a)1702 bool MGL_EXPORT mgl_isnboth(HCDT x, HCDT y, HCDT z, HCDT a)
1703 {
1704 	long n=a->GetNN();
1705 	return x->GetNN()!=n || y->GetNN()!=n || z->GetNN()!=n;
1706 }
1707 //-----------------------------------------------------------------------------
mgl_isboth(HCDT x,HCDT y,HCDT z,HCDT a)1708 bool MGL_EXPORT mgl_isboth(HCDT x, HCDT y, HCDT z, HCDT a)
1709 {
1710 	long n=a->GetNN();
1711 	return x->GetNN()==n && y->GetNN()==n && z->GetNN()==n;
1712 }
1713 //-----------------------------------------------------------------------------
mgl_check_vec3(HMGL gr,HCDT x,HCDT y,HCDT z,HCDT ax,HCDT ay,HCDT az,const char * name)1714 bool MGL_EXPORT mgl_check_vec3(HMGL gr, HCDT x, HCDT y, HCDT z, HCDT ax, HCDT ay, HCDT az, const char *name)
1715 {
1716 // 	if(!gr || !x || !y || !z || !ax || !ay || !az)	return true;		// if data is absent then should be segfault!!!
1717 	long n=ax->GetNx(),m=ax->GetNy(),l=ax->GetNz(), nn=n*m*l;
1718 	if(nn!=ay->GetNN() || nn!=az->GetNN())
1719 	{	gr->SetWarn(mglWarnDim,name);	return true;	}
1720 	if(n<2 || m<2 || l<2)	{	gr->SetWarn(mglWarnLow,name);	return true;	}
1721 	bool both = x->GetNN()==nn && y->GetNN()==nn && z->GetNN()==nn;
1722 	if(!(both || (x->GetNx()==n && y->GetNx()==m && z->GetNx()==l)))
1723 	{	gr->SetWarn(mglWarnDim,name);	return true;	}
1724 	return false;
1725 }
1726 //-----------------------------------------------------------------------------
ClearUnused()1727 void mglBase::ClearUnused()
1728 {
1729 #if MGL_HAVE_PTHREAD
1730 	pthread_mutex_lock(&mutexPnt);	pthread_mutex_lock(&mutexPrm);
1731 #endif
1732 #pragma omp critical
1733 	{
1734 		size_t l=Prm.size();
1735 		// find points which are actually used
1736 		long *used = new long[Pnt.size()];	memset(used,0,Pnt.size()*sizeof(long));
1737 		for(size_t i=0;i<l;i++)
1738 		{
1739 			const mglPrim &p=Prm[i];
1740 			if(p.n1<0)	continue;
1741 			used[p.n1] = 1;
1742 			switch(p.type)
1743 			{
1744 			case 1:	case 4:	if(p.n2>=0)	used[p.n2] = 1;
1745 				break;
1746 			case 2:	if(p.n2>=0 && p.n3>=0)	used[p.n2] = used[p.n3] = 1;
1747 				break;
1748 			case 3:	if(p.n2>=0 && p.n3>=0 && p.n4>=0)
1749 					used[p.n2] = used[p.n3] = used[p.n4] = 1;
1750 				break;
1751 			}
1752 		}
1753 		// now add proper indexes
1754 		l=Pnt.size();
1755 		mglStack<mglPnt> pnt;	pnt.reserve(l);
1756 		for(size_t i=0;i<l;i++)	if(used[i])
1757 		{	pnt.push_back(Pnt[i]);	used[i]=pnt.size()-1;	}
1758 		Pnt = pnt;	pnt.clear();
1759 		// now replace point id
1760 		l=Prm.size();
1761 		for(size_t i=0;i<l;i++)
1762 		{
1763 			mglPrim &p=Prm[i];	p.n1=used[p.n1];
1764 			if(p.type==1 || p.type==4)	p.n2=used[p.n2];
1765 			if(p.type==2)	{	p.n2=used[p.n2];	p.n3=used[p.n3];	}
1766 			if(p.type==3)	{	p.n2=used[p.n2];	p.n3=used[p.n3];	p.n4=used[p.n4];	}
1767 		}
1768 		delete []used;
1769 	}
1770 #if MGL_HAVE_PTHREAD
1771 	pthread_mutex_unlock(&mutexPnt);	pthread_mutex_unlock(&mutexPrm);
1772 #endif
1773 }
1774 //-----------------------------------------------------------------------------
ClearPrmInd()1775 void mglBase::ClearPrmInd()
1776 {
1777 #pragma omp critical(prmind)
1778 	{	if(PrmInd)	delete []PrmInd;	PrmInd=NULL;	}
1779 }
1780 //-----------------------------------------------------------------------------
curve_plot(size_t num,size_t k0,size_t step)1781 void mglBase::curve_plot(size_t num, size_t k0, size_t step)
1782 {
1783 	// exclude straight-line parts
1784 	if(get(MGL_FULL_CURV))	for(size_t i=0;i+1<num;i++)
1785 		line_plot(k0+i*step,k0+(i+1)*step);
1786 	else	for(size_t i=0;i+1<num;i++)
1787 	{
1788 		const mglPoint p1(GetPntP(k0+i*step)), ps(GetPntP(k0+(i+1)*step));
1789 		if(mgl_isnan(p1.x) || mgl_isnan(ps.x))	continue;
1790 		const mglColor c1(GetPntC(k0+i*step));
1791 		// remove duplicates
1792 		for(;i+1<num;i++)
1793 		{
1794 			size_t ii = k0+(i+1)*step;
1795 			const mglPoint pp(GetPntP(ii));
1796 			if(p1!=pp || mgl_isnan(pp.x))	break;
1797 		}
1798 		if(i+1>=num)	break;
1799 
1800 		float t1=-100, t2=100;		// XY angle boundary
1801 		float rg1=-100, rg2=100;	// RG angle boundary
1802 		float gb1=-100, gb2=100;	// GB angle boundary
1803 		size_t k;
1804 		for(k=i+1;k<num;k++)
1805 		{
1806 			const mglPoint p2(GetPntP(k0+k*step)-p1);
1807 			if(mgl_isnan(p2.x))	break;
1808 			float dd=p2.x*p2.x+p2.y*p2.y+p2.z*p2.z;
1809 			if(dd<=0)	continue;	// the same point (micro-loop? :) )
1810 			float t = atan2(p2.y,p2.x), d = atan(0.03/dd);
1811 			if(t1 > t+d || t2 < t-d)	break;		// too curved
1812 			const mglColor c2(GetPntC(k0+(k-1)*step)-c1);	dd = c2.NormS();
1813 			if(dd>0)	// color are different
1814 			{
1815 				float rg = atan2(c2.r,c2.g), gb = atan2(c2.g,c2.b);	d = atan(1e-4/dd);
1816 				if(rg1 > rg+d || rg2 < rg-d || gb1 > gb+d || gb2 < gb-d)	break;		// too curved
1817 				rg1 = rg1<rg-d?rg-d:rg1;	rg2 = rg2>rg+d?rg+d:rg2;	// new RG range
1818 				gb1 = gb1<gb-d?gb-d:gb1;	gb2 = gb2>gb+d?gb+d:gb2;	// new GB range
1819 			}
1820 			t1 = t1<t-d?t-d:t1;	t2 = t2>t+d?t+d:t2;	// new range
1821 		}
1822 		k--;	line_plot(k0+i*step,k0+k*step);	i = k-1;
1823 	}
1824 }
1825 //-----------------------------------------------------------------------------
1826