image.c (3197B)
1 #include <u.h> 2 #include <libc.h> 3 #include <bio.h> 4 #include "sky.h" 5 6 char rad28[] = "0123456789abcdefghijklmnopqr"; 7 8 Picture* 9 image(Angle ra, Angle dec, Angle wid, Angle hig) 10 { 11 Pix *p; 12 uchar *b, *up; 13 int i, j, sx, sy, x, y; 14 char file[50]; 15 Picture *pic; 16 Img* ip; 17 int lowx, lowy, higx, higy; 18 int slowx, slowy, shigx, shigy; 19 Header *h; 20 Angle d, bd; 21 Plate *pp, *bp; 22 23 if(gam.gamma == 0) 24 gam.gamma = -1; 25 if(gam.max == gam.min) { 26 gam.max = 17600; 27 gam.min = 2500; 28 } 29 gam.absgamma = gam.gamma; 30 gam.neg = 0; 31 if(gam.absgamma < 0) { 32 gam.absgamma = -gam.absgamma; 33 gam.neg = 1; 34 } 35 gam.mult1 = 1. / (gam.max - gam.min); 36 gam.mult2 = 255. * gam.mult1; 37 38 if(nplate == 0) 39 getplates(); 40 41 bp = 0; 42 bd = 0; 43 for(i=0; i<nplate; i++) { 44 pp = &plate[i]; 45 d = dist(ra, dec, pp->ra, pp->dec); 46 if(bp == 0 || d < bd) { 47 bp = pp; 48 bd = d; 49 } 50 } 51 52 if(debug) 53 Bprint(&bout, "best plate: %s %s disk %d %s\n", 54 hms(bp->ra), dms(bp->dec), 55 bp->disk, bp->rgn); 56 57 h = getheader(bp->rgn); 58 xypos(h, ra, dec, 0, 0); 59 if(wid <= 0 || hig <= 0) { 60 lowx = h->x; 61 lowy = h->y; 62 lowx = (lowx/500) * 500; 63 lowy = (lowy/500) * 500; 64 higx = lowx + 500; 65 higy = lowy + 500; 66 } else { 67 lowx = h->x - wid*ARCSECONDS_PER_RADIAN*1000 / 68 (h->param[Pxpixelsz]*h->param[Ppltscale]*2); 69 lowy = h->y - hig*ARCSECONDS_PER_RADIAN*1000 / 70 (h->param[Pypixelsz]*h->param[Ppltscale]*2); 71 higx = h->x + wid*ARCSECONDS_PER_RADIAN*1000 / 72 (h->param[Pxpixelsz]*h->param[Ppltscale]*2); 73 higy = h->y + hig*ARCSECONDS_PER_RADIAN*1000 / 74 (h->param[Pypixelsz]*h->param[Ppltscale]*2); 75 } 76 free(h); 77 78 if(lowx < 0) lowx = 0; 79 if(higx < 0) higx = 0; 80 if(lowy < 0) lowy = 0; 81 if(higy < 0) higy = 0; 82 if(lowx > 14000) lowx = 14000; 83 if(higx > 14000) higx = 14000; 84 if(lowy > 14000) lowy = 14000; 85 if(higy > 14000) higy = 14000; 86 87 if(debug) 88 Bprint(&bout, "xy on plate: %d,%d %d,%d\n", 89 lowx,lowy, higx, higy); 90 91 if(lowx >= higx || lowy >=higy) { 92 Bprint(&bout, "no image found\n"); 93 return 0; 94 } 95 96 b = malloc((higx-lowx)*(higy-lowy)*sizeof(*b)); 97 if(b == 0) { 98 emalloc: 99 fprint(2, "malloc error\n"); 100 return 0; 101 } 102 memset(b, 0, (higx-lowx)*(higy-lowy)*sizeof(*b)); 103 104 slowx = lowx/500; 105 shigx = (higx-1)/500; 106 slowy = lowy/500; 107 shigy = (higy-1)/500; 108 109 for(sx=slowx; sx<=shigx; sx++) 110 for(sy=slowy; sy<=shigy; sy++) { 111 if(sx < 0 || sx >= nelem(rad28) || sy < 0 || sy >= nelem(rad28)) { 112 fprint(2, "bad subplate %d %d\n", sy, sx); 113 free(b); 114 return 0; 115 } 116 sprint(file, "%s/%s/%s.%c%c", 117 dssmount(bp->disk), 118 bp->rgn, bp->rgn, 119 rad28[sy], 120 rad28[sx]); 121 122 ip = dssread(file); 123 if(ip == 0) { 124 fprint(2, "can't read %s: %r\n", file); 125 free(b); 126 return 0; 127 } 128 129 x = sx*500; 130 y = sy*500; 131 for(j=0; j<ip->ny; j++) { 132 if(y+j < lowy || y+j >= higy) 133 continue; 134 p = &ip->a[j*ip->ny]; 135 up = b + (higy - (y+j+1))*(higx-lowx) + (x - lowx); 136 for(i=0; i<ip->nx; i++) { 137 if(x+i >= lowx && x+i < higx) 138 *up = dogamma(*p); 139 up++; 140 p += 1; 141 } 142 } 143 free(ip); 144 } 145 146 pic = malloc(sizeof(Picture)); 147 if(pic == 0){ 148 free(b); 149 goto emalloc; 150 } 151 pic->minx = lowx; 152 pic->miny = lowy; 153 pic->maxx = higx; 154 pic->maxy = higy; 155 pic->data = b; 156 strcpy(pic->name, bp->rgn); 157 return pic; 158 }