/******************************************************************************* * E.S.O. - VLT project * * "@(#) $Id: AP_list.C,v 2.43 2004/06/17 23:04:27 vltsccm Exp $" * * who when whatSexte * -------- -------- ---------------------------------------------- * hummel 02/03/98 created * aaguayo 04/12/02 removed warning (APR2003) * aaguayo 28/10/03 PPRS10127 - Added SIMPLE option methods */ /************************************************************************ * NAME * TARG_set - compares target positions with a mask * * SYNOPSIS * TARG_set(); * * PARENT CLASS * AP_mask * * DESCRIPTION * This class maintains a pointer array of targets. * the targets are defined in the parent class AP_target * TARG_set is able to read the target list from an external file. * The purpose if this class is to compare the target list * with an instrumental mask (here the VLT-UT1-FORS1-MOS mask). * and find a optimum position in that way, that the number of * targets matching one of the 19 MOS slitlets is maximized. * E.g. having 30 targets on the field of view at any * distribution; find a telescope pointing position and an * appropriate flansh rotator angle that at least every of the 19 * MOS slits can be associated with one of the objects. * * * PUBLIC METHODS * void add(char* tn, float ta, float td, float tx, * float ty, float tm, int tp); * add one target to the pointer list * * void print(); * print the whole pointer list * * void delete_all(); * delete the whole pointer list * * void set_by_index(int i, char* tn="ds", float ta=0.0, * float td=0.0, float tx=0.0, * float ty=0.0, float tm=0.0, int tp=9); * edit individual target * * void delete_by_index(int indax); * delete individual target from list * * particular MOS applications * * * void load_from_file (char *filename); * open file, read control parameters, read list, * create pointer list * * void SaveBest(char *filename, int slits[]); * open file (p)mos(1/2).(counter).autopos.fims * use ESO-VLTSW format for target aquisition files * saves the best mask position * saves the best slit-let positions * saves the target id for each MOS-slit * this output file can be loaded by FIMS * * void fit_list(); * selects one of the fitting methods * dependent on private variable kind * * PUBLIC DATA MEMBERS * * PROTECTED METHODS * * PROTECTED DATA MEMBERS * * PRIVATE METHODS * float wcs_quality(float mask_ra, float mask_dec); * returns the number of matches for a * non-rotated mask (SHIFT-option); * input is the mask position * The problem is solved in WCS * * float pixel_mos_quality(float mask_x, float mask_y, float mask_rot); * returns the number of MOS matches for a rotated mask * The problem is solved in pixel coordinates * * float pixel_mxu_quality(float mask_x, float mask_y, float mask_rot, int f); * returns the number of created MXU slits for a rotated mask * The problem is solved in pixel coordinates * * void Static_Method(); * just calls pixel_mos_quality with given mask coordinates * * void Static_Simple(); * * void Rotat_Method(); * simple position angle loop and * just calls pixel_mos_quality with given mask coordinates * * void Shift_Method(); * simple fit method for option SHIFT * Puts the mask on the target coords of all * members of the catalog * just calls pixel_mos_quality with given mask coordinates * * void Shift_Method_ATF(); * simple fit method for option SHIFT * Puts the mask on the target coords of all * members of the catalog * ATF= All Targets Fit * just calls wcs_quality with given mask coordinates * * void NoNameFitMethod(); * dummy for the next fitting method * * float a_mean(); * returns the mean (average) RA value of the list * * float d_mean(); * returns the mean (average) DEC value of the list * * float x_mean(); * returns the mean (average) x pixel value of the list * * float y_mean(); * returns the mean (average) y pixel value of the list * * PRIVATE DATA MEMBERS * TARG *targ[2000]; // array of 2000 pointers * static int counter; // counter of the target list * float ra_min,ra_max,dec_mima; // CCD window in arcsec * int hrcoll; // 1_hr_coll 0=sr_coll * float best_ra,best_dec,best_rot; // optimized position * float best_px,best_py; // optimized position * float as2pix; // CCD scale in arcsec to pixel * float PW_rot; // angle between pixel coords and * // wcs = CROTA2 * float best_rot; // angle between pixel coords and mask * int flag[20]; // contains the target list index * int pmos; // 0=mos 1=pmos * int kind; // 0=static 1=shift 2=rotat 3=all * int exposure_counter; // as in the tcl script fims.tcl * float alp0,del0; // for static method * float px0,py0; // for static method * float fits_x_scale, fits_y_scale;// coord orient * float focscale // focal filed scale in arcsec / mm * float slit_width // MOS common slit width, parsed to output * float slit_length // MXU common slit width, used in pixel_mxu_quality * int FORS // 1 or 2, parsed to output * * FILES * AP_list.C AP_list.h * * ENVIRONMENT * * COMMANDS * * RETURN VALUES * * CAUTIONS * * EXAMPLES * * SEE ALSO * * BUGS * *------------------------------------------------------------------------ */ #ifndef __cplusplus #error This is a C++ include file and cannot be used from plain C #endif #include #include #include #include #include #include #include #include "AP_targ.h" #include "AP_list.h" using namespace std; int TARG_set::counter=0; //------------------------------------------- TARG_set::TARG_set() //------------------------------------------- { best_rot=0; best_ra=0; best_dec=0; best_px=0; best_py=0; } //------------------------------------------- float TARG_set::a_mean() //------------------------------------------- { float a=0; for (int i=0; i< counter; i++) a += targ[i]->get_RA(); return a /= counter; } //------------------------------------------- float TARG_set::d_mean() //------------------------------------------- { float d=0; for (int i=0; i< counter; i++) d += targ[i]->get_DEC(); return d /= counter; } //------------------------------------------- float TARG_set::x_mean() //------------------------------------------- { float a=0; for (int i=0; i< counter; i++) a += targ[i]->get_px(); return a /= counter; } //------------------------------------------- float TARG_set::y_mean() //------------------------------------------- { float d=0; for (int i=0; i< counter; i++) d += targ[i]->get_py(); return d /= counter; } //------------------------------------------- void TARG_set::set_by_index(int i, char* ltn, float lta, float ltd, float ltx, float lty, float ltm, int ltp) //------------------------------------------- { if (strcmp(targ[i]->get_nam(),ltn)!=0) targ[i]->set_nam(ltn); if (lta!=0) targ[i]->set_RA(lta); if (ltd!=0) targ[i]->set_DEC(ltd); if (ltx!=0) targ[i]->set_px(ltx); if (lty!=0) targ[i]->set_py(lty); if (ltm!=0) targ[i]->set_MAG(ltm); if (ltp!=9) targ[i]->set_PRI(ltp); targ[i]->print(); } //------------------------------------------- void TARG_set::delete_by_index(int indax) //------------------------------------------- { for (int i=indax; i < counter-1; i++) targ[i]=targ[i+1]; counter--; } //------------------------------------------- void TARG_set::add (char* ltn, float lta, float ltd, float ltx, float lty, float ltm, int ltp) //------------------------------------------- { targ[counter++] = new TARG(ltn, lta, ltd, ltx, lty, ltm, ltp); } //------------------------------------------- void TARG_set::print() //------------------------------------------- { cout << "TargetList:" << endl; for (int i=0; i < counter; i++) targ[i]->print(); } //------------------------------------------- void TARG_set::delete_all () //------------------------------------------- { for (int i=0; i < counter; i++) { delete targ[i]; // cout << "object " << i << " deleted ..." << endl; } // cout << " deleted " << counter << " objects" << endl; counter=0; } //------------------------------------------- void TARG_set::load_from_file (char *filename) //------------------------------------------- { // load the .fims/LOG/APlist.asc file // this list contains in the the first line a number of control parameters // and the resting lines are the catalog (target list) int debug=0; int control_number,p=0,local_counter; char *n; // string pointer to target name float a,d,m=0,px,py,PA_rot; float aa,dummy; // ifstream fileId("/home/hummel/FITS/AutoPos.asc", ios::in); char dirname[260]; if (debug) cout << "Tach in AP.C load_from_file" << endl; sprintf(dirname,"%s%c%s",getenv("LOG_DIR"),'/',filename); if (debug) cout << " >" << dirname << "< >" << filename << "<" << endl; if (!strcmp(dirname,filename)) { cout << "Warning : environment LOG_DIR not found" << endl; } ifstream fileId(dirname, ios::in); // or use 'if (!fileId)' if (fileId.bad()) { cout << "File not found\n"; exit(8); } else { if (debug) cout << "File found\n"; delete_all(); fileId >> exposure_counter; if (debug) cout << "L: EXP_counter =" << exposure_counter << endl; fileId >> FORS; // read either 1 or 2 (FORS1 or FORS2) if (debug) cout << "L: FORS =" << FORS << endl; fileId >> control_number; // saved as mp_header, two flags, if table contains a column for // mag (magnitudes) and/or prio (priority) if (debug) cout << "L: MP_flag =" << control_number << endl; fileId >> kind; // saved as control // 0 static: take the current mask position and position angle // and only move slitlets // 1 shift: move the mask and optimize slitlets, do not rotate // 2 rotat: rotate the mask on the current position, do not shift // 3 all : rotate ans shift mask, optimze slit positions if (debug) cout << "L: kind =" << kind << endl; fileId >> focscale; // saved as focal scale in arcsec / mm if (debug) cout << "L: FocSca =" << focscale << endl; fileId >> pmos; // saved as 0=mos, 1=pmos 3=mxu if (debug) cout << "L: pmos =" << pmos << endl; fileId >> hrcoll; // int flag, if 1 -> high resolution modus // only the 8 inner slitlets can be used if (debug) cout << "L: hrcoll =" << hrcoll << endl; fileId >> as2pix; // scale in arcsec / pixel if (debug) cout << "L: scale =" << as2pix << endl; fileId >> PW_rot; // angle between pixel coords and WCS = CROTA1 if (debug) cout << "L: PW_rot =" << PW_rot << endl; fileId >> PA_rot; // position angle of the mask if (debug) cout << "L: PA_rot =" << PA_rot << endl; fileId >> alp0 >> del0; // RA and DEC of the current mask position if (debug) cout << "L: alp0 del0 =" << alp0 << " | " << del0 << endl; fileId >> px0 >> py0; // current mask position in pixels // in pixels if (debug) cout << "L: px0,py0 =" << px0 << " | " << py0 << endl; fileId >> dec_mima >> ra_min >> ra_max; // the target 'window', the area on the // CCD for which target are allowed to be selected if (debug) cout << "L: dec_mima=" << dec_mima << " arcsec, ra_min=" << ra_min << " arcsec, ra_max=" << ra_max << " arcsec " << endl; fileId >> fits_x_scale >> fits_y_scale; // scale x ><0, scale y ><0 fileId >> slit_width; // the common MOS slit width, parsed to output fileId >> slit_length; // the common MXU slit length in arcsec if (debug) cout << "L: MXU slitlength" << slit_length << endl; if (debug) cout << "L: orient = " << \ fits_x_scale << " | " << fits_y_scale << " | " << slit_width << endl; // no XOR available in C/C++ // if a mirror coord system occurs, Mask-pixel angle is different // aa=(fits_x_scale<0); bb=(fits_y_scale<0); // if ( !(aa && bb) && (aa || bb) ) { // MP_rot = PA_rot + PW_rot; // cout << " verkehrt " << MP_rot << endl; // } else { // MP_rot = -PW_rot - PA_rot; // cout << " normal " << MP_rot << endl; // } // special treatment for wrongly mounted intrument // PA_rot=PA_rot+mount_angle; aa=fits_x_scale*fits_y_scale; // equivalent to coord orient MP_rot = PW_rot - PA_rot*aa; // if (aa==-1.0) { cout << " verkehrt " << MP_rot << endl; } // if (aa==1.0) { cout << " normal " << MP_rot << endl; } // this was previously in pixel-quality // we adapt the wavelength borders according to the // current orientation of the coord system // x_flag is for the AP_mask (if MOS1 is in direction of // increasing or decerasing pixels if (fits_x_scale<0) { dummy=ra_min; ra_min = -ra_max; ra_max = -dummy; x_flag=1; if (debug) { cout << " -X Scale " << endl;} } else { if (debug) { cout << " +X Scale " << endl;} x_flag=0; } if (fits_y_scale<0) { dec_mima *= -1.0; if (debug) {cout << " -Y Scale " << endl;} } else { if (debug) {cout << " +Y Scale " << endl;} } n= new char[40]; // create target name string pointer during run-time fileId >> n; // for EOF problem // cout << n; local_counter=0; // just to be shure while (!fileId.eof() && local_counter<2000) { fileId >> a >> d >> px >> py; // cout << " " << a << " " << " " << d << " " << px << " " << py; switch (control_number) { case 0: break; case 1: { fileId >> p; /* cout << " " << p; */ } break; case 10: { fileId >> m; /* cout << " " << m; */ } break; case 11: { fileId >> m >> p; /* cout << " " << m << " " << p; */ } break; } // cout << endl; cout << setw(4) << setfill('_') << local_counter; cout << " (" << n << " | " << a << " | " << d << " | " << px << " | " << py << " | " << m << " | " << p << ")" << endl; // shift in wcs, rest in pixel add (n,a,d,px,py,m,p); local_counter++; fileId >> n; } } print(); if (debug) print(); if (debug) cout << "L: load ok" << endl; fileId.close(); if (debug) cout << dirname << endl; } //-------------------------------------- void TARG_set::SaveBest(char *filename, int slits[]) //-------------------------------------- { // filename : under which name to save the fims-input file // slits : index array slit[6]=126 means // MOS slit 6 is aasigned to catalog target 126 int i; int dd,mm; char sign; float l_ra,l_dec,work_dec; float ss; int a; // boolean stuff float PA_new; // position angle int debug=0; char dirname[260]; sprintf(dirname,"%s%c%s",getenv("LOG_DIR"),'/',filename); if (!strcmp(dirname,filename)) { cout << "Warning : environment LOG_DIR not found" << endl; } ofstream tpx(dirname, ios::out); // or use 'if (!fileId)' if (tpx.bad()) { cout << "File problems \n"; exit(8); } else { tpx.setf(ios::showpoint|ios::fixed); //fill with '0',, floating-point format tpx.precision(3); // here numbers after . if (debug) { cout.setf(ios::showpoint|ios::fixed); //fill with '0',, floating-point format cout.precision(3); // here numbers after . } tpx << "# created by AutoPos (TOS) V1.4.0" << endl; tpx << "# user date host" << endl; tpx << "# no further comments" << endl; tpx << "# --------------------------" << endl; // --RA--TEL-- // arithmitic deg.xxxx -> ddmmss.xx dd=(int)floor(best_ra/15.0); mm=(int)floor((best_ra/15.0-dd)*60.0); ss=(((best_ra/15-dd)*60.0)-floor((best_ra/15-dd)*60.0))*60.0; tpx << "TEL.TARG.ALPHA "; tpx.fill('0'); tpx.width(2); tpx << dd; tpx.fill('0'); tpx.width(2); tpx << mm; tpx.fill('0'); tpx.width(6); tpx << ss << endl; if (debug) { cout << "TEL.TARG.ALPHA "; cout.fill('0'); cout.width(2); cout << dd; cout.fill('0'); cout.width(2); cout << mm; cout.fill('0'); cout.width(6); cout << ss << endl; } // --DEC--TEL-- work_dec=fabs(best_dec); if (work_dec != best_dec) sign='-'; else sign=' '; dd=(int)floor(work_dec); mm=(int)floor((work_dec-dd)*60.0); ss=(((work_dec-dd)*60.0)-floor((work_dec-dd)*60.0))*60.0; tpx << "TEL.TARG.DELTA "; tpx.fill(' '); tpx.width(1); tpx << sign; tpx.fill('0'); tpx.width(2); tpx << dd; tpx.fill('0'); tpx.width(2); tpx << mm; tpx.fill('0'); tpx.width(6); tpx << ss << endl; if (debug) { cout << "TEL.TARG.DELTA "; cout.fill(' '); cout.width(1); cout << sign; cout.fill('0'); cout.width(2); cout << dd; cout.fill('0'); cout.width(2); cout << mm; cout.fill('0'); cout.width(6); cout << ss << endl; } // --PA--TEL-- // no XOR available in C/C++ // a=(fits_x_scale<0); b=(fits_y_scale<0); // if ( !(a && b) && (a || b) ) { // PA_new = -PW_rot + best_rot; // cout << " verkehrt " << PA_new << endl; // } else { // PA_new = -best_rot - PW_rot; // cout << " normal " << PA_new << endl; // } a=(int)(fits_x_scale*fits_y_scale); // equivalent to coord orient PA_new = a*(best_rot-PW_rot); // special treatment of wrongly mounted instruments // PA_new=PA_new-mount_angle; // if (a==-1.0) { cout << " verkehrt " << PA_new << endl; } // if (a==1.0) { cout << " normal " << PA_new << endl; } tpx << "TEL.ROT.OFFANGLE " << PA_new << endl; if (debug) { cout << "TEL.ROT.OFFANGLE " << PA_new << endl; } // we should also read Epoch from input file tpx << "TEL.TARG.EQUINOX " << 2000.0 << endl; tpx << "TEL.TARG.TYPE COORDINATE" << endl; if (debug) { cout << "TEL.TARG.EQUINOX " << 2000.0 << endl; cout << "TEL.TARG.TYPE COORDINATE" << endl; } // now pointing of mask ok, we begin slit positions, which // are different in (p)mos and mxu // instead of a new class and new inheritance we simply // use a switch switch (pmos) { case 0: case 1: { // M O S P M O S // --TARG-- // standard output will be read by tcl plugin for (i=1; i<=19; i++) { if (slits[i] > -1) { // --NAME-- : INS.TARG16.NAME HR_0815 tpx.precision(5); // here numbers after . tpx << "INS.TARG" << setw(1) << i << ".NAME " << targ[ slits[i] ]->get_nam() << endl; if (debug) { cout.precision(5); // here numbers after . cout << "INS.TARG" << setw(1) << i << ".NAME " << targ[ slits[i] ]->get_nam() << endl; } // --SLIT-WIDTH-TARG tpx.fill(' '); tpx.width(4); tpx.precision(1); tpx << "INS.TARG" << setw(1) << i << ".WID " << slit_width << endl; if (debug) { cout.fill(' '); cout.width(4); cout.precision(1); cout << "INS.TARG" << setw(1) << i << ".WID " << slit_width << endl; } // --RA--TARG-- : INS.TARG16.ALPHA 030503.123 cout.setf(ios::showpoint|ios::fixed); //fill with '0',, floating-point format cout.precision(3); // here numbers after . tpx.precision(3); // here numbers after . l_ra=targ[ slits[i] ]->get_RA(); dd=(int)floor(l_ra/15); mm=(int)floor((l_ra/15-dd)*60.0); ss=(((l_ra/15-dd)*60.0)-floor((l_ra/15-dd)*60.0))*60.0; tpx << "INS.TARG" << setw(1) << i << ".ALPHA "; tpx.fill('0'); tpx.width(2); tpx << dd; tpx.fill('0'); tpx.width(2); tpx << mm; tpx.fill('0'); tpx.width(6); tpx << ss << endl; if (debug) { cout << "INS.TARG" << setw(1) << i << ".ALPHA "; cout.fill('0'); cout.width(2); cout << dd; cout.fill('0'); cout.width(2); cout << mm; cout.fill('0'); cout.width(6); cout << ss << endl; } // Note : for mxu DEC must be the last of a slit // --DEC--TARG-- : INS.TARG16.DELTA -130525.473 tpx.precision(3); // here numbers after . l_dec=targ[ slits[i] ]->get_DEC(); work_dec=fabs(l_dec); if (work_dec != l_dec) sign='-'; else sign=' '; dd=(int)floor(work_dec); mm=(int)floor((work_dec-dd)*60.0); ss=(((work_dec-dd)*60.0)-floor((work_dec-dd)*60.0))*60.0; tpx << "INS.TARG" << setw(1) << i << ".DELTA "; tpx.fill(' '); tpx.width(1); tpx << sign; tpx.fill('0'); tpx.width(2); tpx << dd; tpx.fill('0'); tpx.width(2); tpx << mm; tpx.fill('0'); tpx.width(6); tpx << ss << endl; if (debug) { cout << "INS.TARG" << setw(1) << i << ".DELTA "; cout.fill(' '); cout.width(1); cout << sign; cout.fill('0'); cout.width(2); cout << dd; cout.fill('0'); cout.width(2); cout << mm; cout.fill('0'); cout.width(6); cout << ss << endl; } } } // slit-target assignment for (i=1; i<=19; i++) { if (slits[i] > -1) { tpx << "INS.MOS" << i << ".TARG " << i << endl; } } } break; case 2: { // --------------------- // M X U // --------------------- i=0; while (slits[i] != -1) i++; // --slit number-- : INS.SLIT.NUMBER 45 tpx << "INS.SLIT.NUMBER " << setw(1) << i << endl; if (debug) { cout << "INS.SLIT.NUMBER " << setw(1) << i << endl; } // here we begin the slit loop i=0; while (slits[i] != -1) { // --NAME-- : INS.TARG16.NAME HR_0815 tpx.precision(5); // here numbers after . tpx << "INS.TARG" << setw(1) << i+100 << ".NAME " << targ[ slits[i] ]->get_nam() << endl; if (debug) { cout.precision(5); // here numbers after . cout << "INS.TARG" << setw(1) << i+100 << ".NAME " << targ[ slits[i] ]->get_nam() << endl; } // --SHAPE-- : INS.TARG16.SHAPE STRAIGHT tpx << "INS.TARG" << setw(1) << i+100 << ".SHAPE STRAIGHT" << endl; if (debug) { cout << "INS.TARG" << setw(1) << i+100 << ".SHAPE STRAIGHT" << endl; } // --SLIT-WIDTH-TARG tpx.fill(' '); tpx.width(4); tpx.precision(1); tpx << "INS.TARG" << setw(1) << i+100 << ".WID " << slit_width << endl; if (debug) { cout.fill(' '); cout.width(4); cout.precision(1); cout << "INS.TARG" << setw(1) << i+100 << ".WID " << slit_width << endl; } // --SLIT-LENGTH-TARG tpx.fill(' '); tpx.width(4); tpx.precision(1); tpx << "INS.TARG" << setw(1) << i+100 << ".LEN " << slit_length << endl; if (debug) { cout.fill(' '); cout.width(4); cout.precision(1); cout << "INS.TARG" << setw(1) << i+100 << ".LEN " << slit_length << endl; } // --SLIT-ROT-TARG tpx.fill(' '); tpx.width(4); tpx.precision(1); tpx << "INS.TARG" << setw(1) << i+100 << ".ROT " << 0.000 << endl; if (debug) { cout.fill(' '); cout.width(4); cout.precision(1); cout << "INS.TARG" << setw(1) << i+100 << ".ROT " << 0.000 << endl; } // --RA--TARG-- : INS.TARG16.ALPHA 030503.123 cout.setf(ios::showpoint|ios::fixed); //fill with '0',, floating-point format cout.precision(3); // here numbers after . tpx.precision(3); // here numbers after . l_ra=targ[ slits[i] ]->get_RA(); dd=(int)floor(l_ra/15); mm=(int)floor((l_ra/15-dd)*60.0); ss=(((l_ra/15-dd)*60.0)-floor((l_ra/15-dd)*60.0))*60.0; tpx << "INS.TARG" << setw(1) << i+100 << ".ALPHA "; tpx.fill('0'); tpx.width(2); tpx << dd; tpx.fill('0'); tpx.width(2); tpx << mm; tpx.fill('0'); tpx.width(6); tpx << ss << endl; if (debug) { cout << "INS.TARG" << setw(1) << i+100 << ".ALPHA "; cout.fill('0'); cout.width(2); cout << dd; cout.fill('0'); cout.width(2); cout << mm; cout.fill('0'); cout.width(6); cout << ss << endl; } // Note : for mxu DEC must be the last of a slit // --DEC--TARG-- : INS.TARG16.DELTA -130525.473 tpx.precision(3); // here numbers after . l_dec=targ[ slits[i] ]->get_DEC(); work_dec=fabs(l_dec); if (work_dec != l_dec) sign='-'; else sign=' '; dd=(int)floor(work_dec); mm=(int)floor((work_dec-dd)*60.0); ss=(((work_dec-dd)*60.0)-floor((work_dec-dd)*60.0))*60.0; tpx << "INS.TARG" << setw(1) << i+100 << ".DELTA "; tpx.fill(' '); tpx.width(1); tpx << sign; tpx.fill('0'); tpx.width(2); tpx << dd; tpx.fill('0'); tpx.width(2); tpx << mm; tpx.fill('0'); tpx.width(6); tpx << ss << endl; if (debug) { cout << "INS.TARG" << setw(1) << i+100 << ".DELTA "; cout.fill(' '); cout.width(1); cout << sign; cout.fill('0'); cout.width(2); cout << dd; cout.fill('0'); cout.width(2); cout << mm; cout.fill('0'); cout.width(6); cout << ss << endl; } // cout << i << " " << slits[i] << endl; i++; } } break; default: break; } } tpx.close(); } //-------------------------------------- void TARG_set::fit_list() //-------------------------------------- { // selects one of the fitting methods // dependent on private variable kind // cout << "fit_list " << pmos << " " << kind << endl; if (pmos<2) { switch (kind) { case 0: Static_Method(); break; case 1: Shift_Method(); break; case 2: Rotat_Method(); break; // case 3: All_Method(); break; } } else { switch (kind) { case 0: Static_Method(); break; // case 1: Shift_Method(); break; // case 2: Rotat_Method(); break; // case 3: All_Method(); break; case 4: Simple_Method(); break; } } } //-------------------------------------- void TARG_set::Shift_Method() //-------------------------------------- { // try several mask positions alp0,del0 (not rotation) // and test the quality of the mask positions // the mask positions it self to be tested are // are simply the target positions int tn,i; int best_targ[20]; float a=0,best_quality=0; int debug=0; cout << endl; for (tn=0; tn < counter; tn++) { a=pixel_mos_quality(targ[tn]->get_px(),targ[tn]->get_py(),MP_rot); if (a>best_quality) { best_quality=a; best_ra=targ[tn]->get_RA(); best_dec=targ[tn]->get_DEC(); best_px=targ[tn]->get_px(); best_py=targ[tn]->get_py(); best_rot=0; // the target numbers for each of the 19 slits for (i=0; i<=19; i++) best_targ[i]=flag[i]; if (debug) { cout << tn << " " << a << " " << best_px << " " << best_py << " " << best_targ[0] << " " << best_targ[1] << " " << best_targ[2] << " " << best_targ[3] << " " << best_targ[4] << " " << best_targ[5] << " " << best_targ[6] << " " << best_targ[7] << " " << best_targ[8] << " " << best_targ[9] << " " << best_targ[10] << " " << best_targ[11] << " " << best_targ[12] << " " << best_targ[13] << " " << best_targ[14] << " " << best_targ[15] << " " << best_targ[16] << " " << best_targ[17] << " " << best_targ[18] << " " << best_targ[19] << endl; } } } cout << "shift_on_targets " << counter << " tests, quality=" << best_quality << " number of slits=" << best_targ[0] << endl; // simplex with two parameters RA, DEC /* a=wcs_quality(RA,DEC); if (a>best_quality) { best_quality=a; best_ra=targ[tn]->get_RA(); best_dec=targ[tn]->get_DEC(); best_rot=0; for (i=0; i<=19; i++) best_targ[i]=flag[i]; } */ best_rot=MP_rot; char name[30]; if (pmos) { sprintf(name,"%s%1d%s%1d%s","/pmos",FORS,".",exposure_counter,".auto_pos.fims"); } else { sprintf(name,"%s%1d%s%1d%s","/mos", FORS,".",exposure_counter,".auto_pos.fims"); } SaveBest(name, best_targ); } //-------------------------------------- void TARG_set::Shift_Method_ATF() //-------------------------------------- { // try several mask positions alp0,del0 (not rotation) // and test the quality of the mask positions // the mask positions it self to be tested are // are simply the target positions int tn,i; int best_targ[20]; float a=0,best_quality=0; for (tn=0; tn< counter; tn++) // for (tn=0; tn< 1 ; tn++) { a=wcs_quality(targ[tn]->get_RA(),targ[tn]->get_DEC()); // cout << tn << " qual= " << a << endl; if (a>best_quality) { best_quality=a; best_ra=targ[tn]->get_RA(); best_dec=targ[tn]->get_DEC(); best_px=targ[tn]->get_px(); best_py=targ[tn]->get_py(); best_rot=0; for (i=0; i<=19; i++) best_targ[i]=flag[i]; } } cout << "shift_on_targets " << counter << " tests, quality=" << best_quality << endl; // simplex with two parameters RA, DEC /* a=wcs_quality(RA,DEC); if (a>best_quality) { best_quality=a; best_ra=targ[tn]->get_RA(); best_dec=targ[tn]->get_DEC(); best_rot=0; for (i=0; i<=19; i++) best_targ[i]=flag[i]; } */ best_rot=MP_rot; char name[30]; if (pmos) { sprintf(name,"%s%1d%s%1d%s","/pmos",FORS,".",exposure_counter,".auto_pos.fims"); } else { sprintf(name,"%s%1d%s%1d%s","/mos",FORS,".", exposure_counter,".auto_pos.fims"); } SaveBest(name, best_targ); } //-------------------------------------- void TARG_set::NoNameFitMethod() //-------------------------------------- { } //-------------------------------------- void TARG_set::Static_Method() //-------------------------------------- { // if a mirror coord system occurs, Mask-pixel angle is different // px0 ,and py0 is the current telescope pointing in pixel coords // MP_rot is the current angle between mask and CCD pixel // if (pmos<2) { cout << " Static :" << pixel_mos_quality(px0,py0,MP_rot) << endl; } else { cout << " Static :" << pixel_mxu_quality(px0,py0,MP_rot, 0) << endl; } best_ra=alp0; best_px=0.0; best_dec=del0; best_py=0.0; best_rot=MP_rot; char name[30]; switch (pmos) { case 2: {sprintf(name,"%s%1d%s%1d%s","/mxu", FORS,".",exposure_counter,".auto_pos.fims"); } break; case 1: {sprintf(name,"%s%1d%s%1d%s","/pmos",FORS,".",exposure_counter,".auto_pos.fims"); } break; case 0: {sprintf(name,"%s%1d%s%1d%s","/mos", FORS,".",exposure_counter,".auto_pos.fims"); } break; } // cout << "filename=" << name << endl; // flag is the result field, filled in pixel_mos_quality SaveBest(name, flag); } //-------------------------------------- void TARG_set::Simple_Method() //-------------------------------------- { int i; char name[30]; //cout << " Static :" << pixel_mxu_quality(px0,py0,MP_rot,1) << endl; cout << " Static: with " << counter << endl; // AAG PPRS 10127 best_ra=alp0; best_px=0.0; best_dec=del0; best_py=0.0; best_rot=MP_rot; for (i=0; i ra_min && del_ra < ra_max) if (del_dec < dec_mima/2) { // cout << "Q: targ #" << tn << " inwindow (" << // targ[tn]->get_nam() << "," << // targ[tn]->get_RA() << "," << // targ[tn]->get_DEC() << ")" << endl; // this has to be optimized qual=0; iqual=0; l_q=0; // if nothing found for (slit_nr=1; slit_nr <= max_sl; slit_nr++) { // test only on the free slitlets if (flag[slit_nr] == -1) { // cout << " <" << slit_nr << "> "; // slitlet width=22'', slitlet width/2 = 11'' // distance between two slitlets 22.73'' // >actual slit-pos< >targpos< diff=(22.73*(slit_nr-10)-del_dec)/11.0; // cout << " #" << // slit_nr << " ''=" << 22.73*(slit_nr-10) << // " del_dec=" << del_dec << " diff=" << diff << endl; // 22 arcsec slitlet width // 1/e at 1/6 slitwidth // or just the distance, linear approach // if (fabs(22.0*(slit_nr-10)-targ[i]->get_DEC())<10.0) if (fabs(diff)<0.9) l_q=exp(-diff*diff); else l_q=0.0; // cout << " q=" << l_q << endl; if (l_q>qual) { qual=l_q; // new quality iqual=slit_nr; // the slitlet number ldist=diff; // distance from slitlet center } } } if (qual>0) { // the target could be associated to a slitlet flag[iqual]=tn; // slitlet occupied by tn alp[iqual]=targ[tn]->get_RA(); // save slitlet position del[iqual]=targ[tn]->get_DEC(); // save slitlet position dis[iqual]=ldist; // distance from centre found++; // cout << " slitlet #" << iqual << // " dis=" << dis[iqual] << // " qual=" << qual << // " name=" << targ[tn]->get_nam() << // " (" << alp[iqual] << "," << del[iqual] << ")" << endl; } qual_sum += qual; } } // cout << "Q: qual=" << qual_sum << // " (" << mask_ra << ", " << mask_dec << ") " << // " | " << found << " of 19 found" << endl; return qual_sum; } //------------------------------------------------------ float TARG_set::pixel_mos_quality(float mask_x, float mask_y, float mask_rot) //------------------------------------------------------ { // this function calculates the match quality of a // given mask position and PA (=position angle). // the mask position is decsribed by mask_x and mask_y (pixel // Each target in the target list is testet if it // matches the CCD window ( < full CCD) // Then we search for the closest slitlet // and set the occupy flag for the closest slitlet // ra,dec in deg (e.q. -72.3456) // del_ra,dec_mima in arcsec float x,y; float qual_sum=0.0; // the return value float l_q,diff,ldist=0,qual; int tn; // target number in the list int max_sl=19; // max of available slits =f(coll,lambda) int i,iqual,slit_nr; float alp[20],del[19]; // the positions of the slitlets float dis[20]; // distance of target from centeral slit line int found=0; // number of occupied slitlets int reoccup_slit=0; int debug=0; // create an instance of the mask // cout << "PQ arcsec: " << as2pix << " " << ra_min << " " // << ra_max << " " << dec_mima << endl; // x>0, ra>0 y>0,dec>0 AP_mask my(as2pix,focscale,ra_min,ra_max,dec_mima,x_flag); // the working mask // maintenance // my.test(); // position the mask // cout << "mask_rot = " << mask_rot << endl; my.reset_mask(mask_x,mask_y,mask_rot); // window array automatically created, // AP_mask::wx, AP_mask::wy if (debug) my.print_mask(); // init flag // -1 = free (and enabled), -2 = disabled, 0..counter = occupied flag[0]=-1; for (i=1; i<=19; i++) { flag[i]=-1; alp[i]=0; del[i]=0; dis[i]=99; } // these slits cannot be used in high resolution if (hrcoll) { for (i=1; i<=5; i++) flag[i]=-2; for (i=15; i<=19; i++) flag[i]=-2; } // even slits cannot be used in mode PMOS if (pmos) for (i=1; i<=19; i += 2) flag[i]=-2; // cout << "Q: " << hrcoll << " " << mask_ra << " " << mask_dec << endl; // for (i=0; i<=19; i++) // cout << " (" << i << "," << flag[i] << ")"; // cout << endl; // cout << "Q:begin ra_min=" << ra_min << " ra_max=" << ra_max << // " dec_mima/2=" << dec_mima/2 << endl; // Now the priority flag: // if there are two targets for one slit, take the // prio=1 target in stead of the target with the better quality // scan all target stars in the list; tn=target_number for (tn=0; tn< counter; tn++) { x=targ[tn]->get_px(); y=targ[tn]->get_py(); // copy target position in CCDwindow array // if the target is in the wavelength-confined CCD-window if (my.p_in_window(x,y)) { if (debug) { cout << tn << " (" << x << "," << y << ") in window" << endl; } qual=0; iqual=0; l_q=0; // if nothing found for (slit_nr=1; slit_nr <= max_sl; slit_nr++) { // test only on the free slitlets // if (flag[slit_nr] == -1) { // test on the free and on already occupied slit blades if (flag[slit_nr] >= -1) { // cout << " <" << slit_nr << "> "; // slit blade width=22'', slit blade width/2 = 11'' // distance between two slit blades 22.73'' // >actual slit-pos< >targpos< diff=my.pr_on_slitlet(x,y,slit_nr); // if (debug) { cout << " #" << slit_nr // << " diff=" << diff << endl; // } // 22 arcsec slit blade width // 1/e at 1/6 slitwidth // or just the distance, linear approach // if (fabs(22.0*(slit_nr-10)-targ[i]->get_py())<10.0) // l_q = local_quality simple Gaussian weight of the distance // l_q =exp(-diff*diff) or l_q=1.0-fabs(diff) if (fabs(diff)<0.9) l_q=1.0-fabs(diff); else l_q=0.0; // cout << " q=" << l_q << endl; if (l_q>qual) { qual=l_q; // new quality of the current target iqual=slit_nr; // the best slit blade number // for the current target ldist=diff; // distance from slit blade center } } // enabled slit blades } // slit blade loop if (qual>0) { // target 'tn' matches a slit blade if (debug) { cout << "Target " << tn << " best slit #" << iqual << " qual=" << qual << " dist=" << ldist << endl; } if (flag[iqual]== -1) { // the current slit is still free: insert anyway found++; // increment if the slit-blade was free qual_sum += ldist; reoccup_slit=1; if (debug) { cout << " occupy free slit" << endl;} } else { // the slit is already occupied by target flag[iqual] if (debug) { cout << "new/old prio " << targ[tn]->get_PRI() << targ[flag[iqual]]->get_PRI() << endl; } // now some priority logics // ------------------------ // new_prio>old_prio : replace // new_prio=dis[iqual] forget it if (targ[tn]->get_PRI() > targ[flag[iqual]]->get_PRI()) { // if new target priority larger than // occupying target priority // correct the sum, if already occupied qual_sum += ldist-dis[iqual]; reoccup_slit=1; if (debug) { cout << " better prio replaced" << endl;} } else { if (targ[tn]->get_PRI() < targ[flag[iqual]]->get_PRI()) { // worser, do not replace reoccup_slit=0; if (debug) { cout << "worser prio forget it" << endl;} } else { // equal priority, check on quality if (ldistget_RA(); // save slitlet position del[iqual]=targ[tn]->get_DEC(); // save slitlet position dis[iqual]=ldist; // distance from centre } /* cout << " slitlet #" << iqual << " dis=" << dis[iqual] << " qual=" << qual << " name=" << targ[tn]->get_nam() << " (" << alp[iqual] << "," << del[iqual] << ")" << endl; */ } // if in window else { if (debug) { cout << tn << " (" << x << "," << y << ") NOT in window" << endl; } } } // for all targets if (debug) { cout << "Q: qual=" << qual_sum << " (" << mask_x << ", " << mask_y << mask_rot << ") " << " | " << found << " of 19 occupied" << endl; } // the number of occupied slits flag[0]=found; return qual_sum; } //------------------------------------------------------ float TARG_set::pixel_mxu_quality(float mask_x, float mask_y, float mask_rot, int f) //------------------------------------------------------ { // this function calculates the match quality of a // given mask position and PA (=position angle). // the mask position is decsribed by mask_x and mask_y (pixel // Each target in the target list is testet if it // matches the CCD window ( < full CCD) // Then we search for the closest slitlet // and set the occupy flag for the closest slitlet // ra,dec in deg (e.q. -72.3456) // del_ra,dec_mima in arcsec float x,y; int tn; // target number in the list int max_sl=1000; // max of available slits =f(coll,lambda) int i; float xslit,yslit; // y-value of mxu slit int debug=0; int assigned; float sl; int s_pointer; // create an instance of the mask // cout << "PQ arcsec: " << as2pix << " " << ra_min << " " // << ra_max << " " << dec_mima << endl; // x>0, ra>0 y>0,dec>0 AP_mask my(as2pix,focscale,ra_min,ra_max,dec_mima,x_flag); // the working mask // position the mask // cout << "mask_rot = " << mask_rot << endl; // mask_x, maksy are the current position of the mask in CCD pixels // mask_rot is the position angle of the mask with respect to the // CCD pixel coords my.reset_mask(mask_x,mask_y,mask_rot); // window array automatically created, // AP_mask::wx, AP_mask::wy // if (debug) my.print_mask(); // default : not assigned to any catalog target position for (i=0; iget_py(); xslit = targ[flag[i]]->get_px(); // coord transformation of slit from frame to focfield pixels my.p2m(&xslit,&yslit); // y = ytarg = y in CCD frame pixel if ( y < (yslit+sl) && y > (yslit-sl) && f!=1) { if (debug) { cout << "target " << tn << " with y=" << y << " occupied by slit " << i << " with y=" << yslit << endl; } assigned=0; } } if (assigned) { if (debug) { cout << "Assign target " << tn << " with y=" << y << " to slit :" << s_pointer << endl; } flag[s_pointer]=tn; ++s_pointer; } } // if in window else { if (debug) { cout << tn << " (" << x << "," << y << ") NOT in window" << endl; } } } // for all targets if (debug) { for (i=0; i < s_pointer; i++) { cout << " slit= " << i << " target = " << flag[i] << endl; } } return s_pointer*1.0; }