util.c 4.7 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365
  1. #include <u.h>
  2. #include <libc.h>
  3. #include <bio.h>
  4. #include "sky.h"
  5. double PI_180 = 0.0174532925199432957692369;
  6. double TWOPI = 6.2831853071795864769252867665590057683943387987502;
  7. double LN2 = 0.69314718055994530941723212145817656807550013436025;
  8. static double angledangle=(180./PI)*MILLIARCSEC;
  9. int
  10. rint(char *p, int n)
  11. {
  12. int i=0;
  13. while(*p==' ' && n)
  14. p++, --n;
  15. while(n--)
  16. i=i*10+*p++-'0';
  17. return i;
  18. }
  19. DAngle
  20. dangle(Angle angle)
  21. {
  22. return angle*angledangle;
  23. }
  24. Angle
  25. angle(DAngle dangle)
  26. {
  27. return dangle/angledangle;
  28. }
  29. double
  30. rfloat(char *p, int n)
  31. {
  32. double i, d=0;
  33. while(*p==' ' && n)
  34. p++, --n;
  35. if(*p == '+')
  36. return rfloat(p+1, n-1);
  37. if(*p == '-')
  38. return -rfloat(p+1, n-1);
  39. while(*p == ' ' && n)
  40. p++, --n;
  41. if(n == 0)
  42. return 0.0;
  43. while(n-- && *p!='.')
  44. d = d*10+*p++-'0';
  45. if(n <= 0)
  46. return d;
  47. p++;
  48. i = 1;
  49. while(n--)
  50. d+=(*p++-'0')/(i*=10.);
  51. return d;
  52. }
  53. int
  54. sign(int c)
  55. {
  56. if(c=='-')
  57. return -1;
  58. return 1;
  59. }
  60. char*
  61. hms(Angle a)
  62. {
  63. static char buf[20];
  64. double x;
  65. int h, m, s, ts;
  66. x=DEG(a)/15;
  67. x += 0.5/36000.; /* round up half of 0.1 sec */
  68. h = floor(x);
  69. x -= h;
  70. x *= 60;
  71. m = floor(x);
  72. x -= m;
  73. x *= 60;
  74. s = floor(x);
  75. x -= s;
  76. ts = 10*x;
  77. sprint(buf, "%dh%.2dm%.2d.%ds", h, m, s, ts);
  78. return buf;
  79. }
  80. char*
  81. dms(Angle a)
  82. {
  83. static char buf[20];
  84. double x;
  85. int sign, d, m, s, ts;
  86. x = DEG(a);
  87. sign='+';
  88. if(a<0){
  89. sign='-';
  90. x=-x;
  91. }
  92. x += 0.5/36000.; /* round up half of 0.1 arcsecond */
  93. d = floor(x);
  94. x -= d;
  95. x *= 60;
  96. m = floor(x);
  97. x -= m;
  98. x *= 60;
  99. s = floor(x);
  100. x -= s;
  101. ts = floor(10*x);
  102. sprint(buf, "%c%d°%.2d'%.2d.%d\"", sign, d, m, s, ts);
  103. return buf;
  104. }
  105. char*
  106. ms(Angle a)
  107. {
  108. static char buf[20];
  109. double x;
  110. int d, m, s, ts;
  111. x = DEG(a);
  112. x += 0.5/36000.; /* round up half of 0.1 arcsecond */
  113. d = floor(x);
  114. x -= d;
  115. x *= 60;
  116. m = floor(x);
  117. x -= m;
  118. x *= 60;
  119. s = floor(x);
  120. x -= s;
  121. ts = floor(10*x);
  122. if(d != 0)
  123. sprint(buf, "%d°%.2d'%.2d.%d\"", d, m, s, ts);
  124. else
  125. sprint(buf, "%.2d'%.2d.%d\"", m, s, ts);
  126. return buf;
  127. }
  128. char*
  129. hm(Angle a)
  130. {
  131. static char buf[20];
  132. double x;
  133. int h, m, n;
  134. x = DEG(a)/15;
  135. x += 0.5/600.; /* round up half of tenth of minute */
  136. h = floor(x);
  137. x -= h;
  138. x *= 60;
  139. m = floor(x);
  140. x -= m;
  141. x *= 10;
  142. n = floor(x);
  143. sprint(buf, "%dh%.2d.%1dm", h, m, n);
  144. return buf;
  145. }
  146. char*
  147. hm5(Angle a)
  148. {
  149. static char buf[20];
  150. double x;
  151. int h, m;
  152. x = DEG(a)/15;
  153. x += 2.5/60.; /* round up 2.5m */
  154. h = floor(x);
  155. x -= h;
  156. x *= 60;
  157. m = floor(x);
  158. m -= m % 5;
  159. sprint(buf, "%dh%.2dm", h, m);
  160. return buf;
  161. }
  162. char*
  163. dm(Angle a)
  164. {
  165. static char buf[20];
  166. double x;
  167. int sign, d, m, n;
  168. x = DEG(a);
  169. sign='+';
  170. if(a<0){
  171. sign='-';
  172. x=-x;
  173. }
  174. x += 0.5/600.; /* round up half of tenth of arcminute */
  175. d = floor(x);
  176. x -= d;
  177. x *= 60;
  178. m = floor(x);
  179. x -= m;
  180. x *= 10;
  181. n = floor(x);
  182. sprint(buf, "%c%d°%.2d.%.1d'", sign, d, m, n);
  183. return buf;
  184. }
  185. char*
  186. deg(Angle a)
  187. {
  188. static char buf[20];
  189. double x;
  190. int sign, d;
  191. x = DEG(a);
  192. sign='+';
  193. if(a<0){
  194. sign='-';
  195. x=-x;
  196. }
  197. x += 0.5; /* round up half degree */
  198. d = floor(x);
  199. sprint(buf, "%c%d°", sign, d);
  200. return buf;
  201. }
  202. char*
  203. getword(char *ou, char *in)
  204. {
  205. int c;
  206. for(;;) {
  207. c = *in++;
  208. if(c == ' ' || c == '\t')
  209. continue;
  210. if(c == 0)
  211. return 0;
  212. break;
  213. }
  214. if(c == '\'')
  215. for(;;) {
  216. if(c >= 'A' && c <= 'Z')
  217. c += 'a' - 'A';
  218. *ou++ = c;
  219. c = *in++;
  220. if(c == 0)
  221. return 0;
  222. if(c == '\'') {
  223. *ou = 0;
  224. return in-1;
  225. }
  226. }
  227. for(;;) {
  228. if(c >= 'A' && c <= 'Z')
  229. c += 'a' - 'A';
  230. *ou++ = c;
  231. c = *in++;
  232. if(c == ' ' || c == '\t' || c == 0) {
  233. *ou = 0;
  234. return in-1;
  235. }
  236. }
  237. }
  238. /*
  239. * Read formatted angle. Must contain no embedded blanks
  240. */
  241. Angle
  242. getra(char *p)
  243. {
  244. Rune r;
  245. char *q;
  246. Angle f, d;
  247. int neg;
  248. neg = 0;
  249. d = 0;
  250. while(*p == ' ')
  251. p++;
  252. for(;;) {
  253. if(*p == ' ' || *p=='\0')
  254. goto Return;
  255. if(*p == '-') {
  256. neg = 1;
  257. p++;
  258. }
  259. if(*p == '+') {
  260. neg = 0;
  261. p++;
  262. }
  263. q = p;
  264. f = strtod(p, &q);
  265. if(q > p) {
  266. p = q;
  267. }
  268. p += chartorune(&r, p);
  269. switch(r) {
  270. default:
  271. Return:
  272. if(neg)
  273. d = -d;
  274. return RAD(d);
  275. case 'h':
  276. d += f*15;
  277. break;
  278. case 'm':
  279. d += f/4;
  280. break;
  281. case 's':
  282. d += f/240;
  283. break;
  284. case L'°':
  285. d += f;
  286. break;
  287. case '\'':
  288. d += f/60;
  289. break;
  290. case '\"':
  291. d += f/3600;
  292. break;
  293. }
  294. }
  295. }
  296. double
  297. xsqrt(double a)
  298. {
  299. if(a < 0)
  300. return 0;
  301. return sqrt(a);
  302. }
  303. Angle
  304. dist(Angle ra1, Angle dec1, Angle ra2, Angle dec2)
  305. {
  306. double a;
  307. a = sin(dec1) * sin(dec2) +
  308. cos(dec1) * cos(dec2) *
  309. cos(ra1 - ra2);
  310. a = atan2(xsqrt(1 - a*a), a);
  311. if(a < 0)
  312. a = -a;
  313. return a;
  314. }
  315. int
  316. dogamma(Pix c)
  317. {
  318. float f;
  319. f = c - gam.min;
  320. if(f < 1)
  321. f = 1;
  322. if(gam.absgamma == 1)
  323. c = f * gam.mult2;
  324. else
  325. c = exp(log(f*gam.mult1) * gam.absgamma) * 255;
  326. if(c > 255)
  327. c = 255;
  328. if(gam.neg)
  329. c = 255-c;
  330. return c;
  331. }