/* test routine for C access to TRACE housekeeping data */ #include #include #include #define ABS(x) ((x)>=0?(x):-(x)) #define MIN(a,b) (((a)<(b))?(a):(b)) #define MAX(a,b) (((a)>(b))?(a):(b)) #define X_INT(x,n) (((int) x[n+1]) & 0xff) + \ ((((int) x[n])<<8 )& 0xff00) + \ ((((int) x[n+3])<<16 )& 0xff0000) + \ ((((int) x[n+2])<<24 )& 0xff000000) #define X_INT2(x,n) (((int) x[n+3]) & 0xff) + \ ((((int) x[n+2])<<8 )& 0xff00) + \ ((((int) x[n+1])<<16 )& 0xff0000) + \ ((((int) x[n])<<24 )& 0xff000000) #define X_SHORT(x,n) (((int) x[n+1]) & 0xff) + \ ((((int) x[n])<<8 )& 0xff00) /* these are times of offset changes */ double trace_wave_offsets_times[] = {685757277.0, 725760479.0}; /* offsets for time bands for the 5 UV bands and a common one for the EUV bands */ float trace_wave_offsets_euv [2][2] = { { 4.3,-7.2}, {2.0, -6.4} }; float trace_wave_offsets_1216 [2][2] = { { -6.2, 2.5}, {-4.2, 2.3} }; float trace_wave_offsets_1600 [2][2] = { {-1.2, 1.2 }, {-1.1, 1.5} }; float trace_wave_offsets_1700 [2][2] = { {-1.2, 1.8 }, {-1.3, 2.4} }; float trace_wave_offsets_1550 [2][2] = { {-0.5, 2.5 }, {-0.6, 1.3} }; float trace_wave_offsets_wl [2][2] = { {0,0 }, {0,0} }; float trace_plate_scale = 0.504; /*--------------------------------------------------------------------------*/ /* wedge2solar This function returns the x and z solar angles associated with a specified set of wedge positions (w1,w2). Call - wedge2solar (w1,w2,&x,&z) w1 and w2 are int wedge positions x and z are pointers to floats for the angles Method - Cribbed from Ted Tarbell's IDL proc. Bugs - Doesn't allow for xsun,zsun, rmax to change. Will add some environment variables later */ void wedge2solar (int w1,int w2, float *x, float *z) { int m1,m2; double th1,th2; double xsun,zsun,rmax; double radeg = 57.2957795; // first flight values (from R. McMullen) xsun = -13.; zsun = -82.; rmax = 1800.0 * 1.040; m1 = w1; m2 = 180 - w2; th1 = 2*m1 / radeg; th2 = -2*m2/radeg; /* changed both signs (x and z) 4/22/98 */ *x = -rmax*(sin(th1)+sin(th2)) + xsun; *z = -rmax*(cos(th1)+cos(th2)) + zsun; //printf("wedge2solar result = %f, %f\n", *x, *z); } /*--------------------------------------------------------------------------*/ int trace_ccdpositions(dp_head,k,ix,iy,ix2,iy2,nx,ny,bm,xsum,amplifier) /* returns CCD position information extracted from header array (in dp_head) for image k, the remaining arguments are returned values so must be passed by reference */ char *dp_head; int k, *ix, *iy, *ix2, *iy2, *nx, *ny, *bm, *xsum, *amplifier; { int sai, ysum, amp, btyo, koff, xq, sa, tmp, nc, nr, rl, ixs, iys; char *head, *parea; /* the sai is the source area index, 0-2 in a 16 bit word, just need the low byte (big endian style) */ koff = k*1010; /* be careful, no range checking on image # */ head = dp_head + koff; sai = head[413*2+1]; switch (sai) { case 0: parea = head + 350*2; break; case 1: parea = head + 366*2; break; case 2: parea = head + 382*2; break; default: { printf("illegal source area index in trace_ccdpositions = %d\n", sai); return 0; } } *xsum = head[102*2+1]; ysum = head[103*2]; if (*xsum != ysum) { printf("problem with camera summing factors, %d, %d\n", *xsum, ysum); return 0; } amp = (head[102*2]>>4 ) & 0x1; //printf("xsum, ysum, amp, %d, %d %d\n", *xsum, ysum,amp); btyo = (int) ((head[312*2]<<8) & 0xff00) + (int) (head[312*2+1]& 0xff); nc = X_SHORT(parea, 10); nr = X_SHORT(parea, 12); xq = X_SHORT(parea, 16); xq = MAX(xq, 1); //printf("nc, nr, xq = %d %d %d\n", nc, nr, xq); *bm = xq; /* nc and nr are image size before binning in the DHC but after summing */ /* image size, after binning in nx and ny */ *nx = nc/xq; *ny = nr/xq; tmp = X_INT(parea, 6); sa = (tmp & 0x000fffff) - 1300; /* starting address */ /* row length in the DHC, varies with camera sum */ rl = (int) ((parea[14]<<8) & 0xff00) + (int) (parea[15]& 0xff); if (rl <= 0) { printf("illegal row length = %d, probably a bad header\n",rl); return 0; } /* location in the summed image (but not binned) */ ixs = sa%rl; iys = sa/rl; /* add the first row (important for partial reads) */ iys = iys + btyo; /* here we get ix, iy for corner of subarea */ *ix = ixs * *xsum; *iy = iys * ysum; /* if this was A amp, ix, iy are the starting pixels in the full CCD the A amp case has amp = 0 */ if (amp) { /* the B amp case */ *ix = 1024 - (ixs + nc)* *xsum; /* note that the equivalent expression is 1023 - ((ixs + nc -1)*xsum + xsum - 1) */ *iy = 1024 - (iys + nr)*ysum; } /* then the other edge of the cutout is */ *ix2 = *ix + nc* *xsum - 1; *iy2 = *iy + nr*ysum - 1; //printf(" ix, ix2, iy, iy2 = %d %d %d %d\n", *ix, *ix2, *iy, *iy2); *amplifier = amp; return 1; } /*------------------------------------------------------------------------- */ int trace_wave_offset(wave, tai, xoff, yoff) int wave; double tai; float *xoff, *yoff; { float offset; int it; /* returns a time dependent offset for the given wavelength number */ /* only 2 sets available right now */ if (tai > trace_wave_offsets_times[0]) it = 1; else it = 0; //printf("wave = %d, it = %d\n", wave, it); if (wave < 12) { /* the EUV case */ *xoff = trace_wave_offsets_euv[it][0]; *yoff = trace_wave_offsets_euv[it][1]; return 1; } switch (wave) { case 13: /* 1216 */ *xoff = trace_wave_offsets_1216[it][0]; *yoff = trace_wave_offsets_1216[it][1]; break; case 15: /* wl */ *xoff = trace_wave_offsets_wl[it][0]; *yoff = trace_wave_offsets_wl[it][1]; break; case 16: /* 1600 */ *xoff = trace_wave_offsets_1600[it][0]; *yoff = trace_wave_offsets_1600[it][1]; break; case 18: /* 1700 */ *xoff = trace_wave_offsets_1700[it][0]; *yoff = trace_wave_offsets_1700[it][1]; break; case 21: /* 1550 */ *xoff = trace_wave_offsets_1550[it][0]; *yoff = trace_wave_offsets_1550[it][1]; break; } return 1; } /*------------------------------------------------------------------------- */ int trace_mechs(dp_head,k, quad_encoder, quad, f1, f2, w1, w2, focus) /* returns some mechanism information extracted from header array (in dp_head) for image k, the remaining arguments are returned values so must be passed by reference */ char *dp_head; int k, *quad, *f1, *f2, *w1, *w2, *focus, *quad_encoder; { char *head, *parea; int koff; koff = k*1010; /* be careful, no range checking on image # */ head = dp_head + koff; *quad_encoder = head[529] & 0x3f; *quad = head[529] >> 6; *f1 = head[525]; *f2 = head[526]; *w1 = head[527]; *w2 = 180 - head[528]; *focus = head[524]; return 1; } /*------------------------------------------------------------------------- */ int show_header(dp_head,k,ny) /* show_header */ char *dp_head; int k, ny; /* dp_head must point to a dp header array, 1010 bytes/image as read from a raw TRACE hourly file k is the first image in the header to show, ny the number, there is no check on the size of ny */ /* this is primarily a test code */ { int nx, i, iq, sl; int start, j, temp, ik, n; int seqid, fn_inseq, frm_num, trgid, frmid, mseqid, xsum, amp; int ix, iy, ix2, iy2, nx_ccd, ny_ccd, bm, stat; int quad_encoder, quad, f1, f2, w1, w2, focus, wave; float xsun,ysun,xoff,yoff,trace_xc,trace_yc,trace_pixel_size,trace_exp; double dq, tai; short mask, *pvar_short; char *dp_ptr, *phead; char *ppq, *sq, *q, *qb; nx = 1010; /* number of bytes per image header */ dp_ptr = dp_head + nx*k; printf(" seqid frm_num fn_inseq trgid frmid mseqid TAI\n"); phead = dp_ptr; while (ny--) { /* loop over images */ /* sequence and various other program ID's */ seqid = X_INT(phead,198*2); frm_num = X_INT2(phead,202*2); fn_inseq = X_SHORT(phead,201*2); trgid = X_SHORT(phead,205*2); frmid = X_SHORT(phead,294*2); mseqid = X_SHORT(phead,244*2); /* note - corrected 4/26/2001, was 228 (wrong)*/ printf("%10d%10d%10d%10d%10d%10d",seqid,frm_num,fn_inseq,trgid,frmid,mseqid); /* image position on detector */ stat = trace_ccdpositions(phead,0, &ix, &iy, &ix2, &iy2, &nx_ccd, &ny_ccd, &bm, &xsum, &); if (stat == 0) { return -1; } /* some mechanism positions */ stat = trace_mechs(phead,0, &quad_encoder, &quad, &f1, &f2, &w1, &w2, &focus); //printf("mechanisms - quad_encoder, quad, f1, f2, w1, w2, focus\n%10d%10d%10d%10d%10d%10d%10d\n", // quad_encoder, quad, f1, f2, w1, w2, focus); if (stat == 0) { return -1; } /* spacecraft pointing */ wedge2solar(w1, w2, &xsun, &ysun); //printf("xsun, ysun = %f %f\n", xsun, ysun); /* (xsun, ysun) is the solar position (approx.) for the white light image, have to apply offsets for the others */ /* get the wavelength index */ wave = (X_SHORT(phead,296*2) >> 7) & 0x1f; //printf("wl = %d\n", wave); /* time, TAI from spacecraft time */ iq = X_INT2(phead,3*2); dq = (double) iq; iq = X_SHORT(phead,3*2+4); dq = dq + ((double) iq)/65536. - 271641594.0 + 25.0; tai = dq; printf(" %.2f\n", tai); /* position offsets, function of wavelength and time */ stat = trace_wave_offset(wave, tai, &xoff, &yoff); //printf("xoff, yoff = %f %f\n", xoff, yoff); //printf("xsun, ysun = %f %f\n", xsun, ysun); /* the final position depends on telescope pointing, position on detector, and time varying offsets for each band */ trace_xc = xsun - trace_plate_scale*(.5*(iy+iy2) - 511.5 - xoff); trace_yc = ysun - trace_plate_scale*(.5*(ix+ix2) - 511.5 - yoff); //printf("trace_xc, trace_yc = %g %g\n",trace_xc, trace_yc); trace_pixel_size = trace_plate_scale*bm*xsum; /* exposure time */ { int dark_flag, jq, tmp2, tmp3, tmp4, n, sweeps; float nq, t, dt; dark_flag = 0; if (phead[296*2] & 0x80) dark_flag = 1; //printf("dark flag = %d\n", dark_flag); jq = X_SHORT(phead, 265*2); jq = jq & 0xffff; tmp2 = X_SHORT(phead, 209*2); tmp3 = X_SHORT(phead, 210*2); sweeps = phead[268*2+1]; if (dark_flag) { trace_exp = jq * 0.004; } else { /* the open/close times, lsb is 1/256 ms */ dt = 0.00000390625 * (double) (tmp2 - tmp3); /* might be negative */ /* use the coarse exposure to resolve wrap arounds */ nq = (jq - 250.*dt)/64.; n = (int) (nq+.5); /* get wrap count */ t = .256*n + dt; /* each wrap is 256 ms, add dt for total */ /* if we get zero, it was a short exposure and we need sweep count */ if (t == 0.0) t = 0.0016 * (1 + sweeps); //printf("dt, n, sweeps = %g %d %d\n", dt, n, sweeps); trace_exp = t; } } printf("trace_exp = %g\n", trace_exp); phead += nx; } return 1; } /*------------------------------------------------------------------------- */