vircam_jmp_utils.c

00001 /* $Id: vircam_jmp_utils.c,v 1.33 2008/07/10 13:05:53 jim Exp $
00002  *
00003  * This file is part of the VIRCAM Pipeline
00004  * Copyright (C) 2005 Cambridge Astronomy Survey Unit
00005  *
00006  * This program is free software; you can redistribute it and/or modify
00007  * it under the terms of the GNU General Public License as published by
00008  * the Free Software Foundation; either version 2 of the License, or
00009  * (at your option) any later version.
00010  *
00011  * This program is distributed in the hope that it will be useful,
00012  * but WITHOUT ANY WARRANTY; without even the implied warranty of
00013  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00014  * GNU General Public License for more details.
00015  *
00016  * You should have received a copy of the GNU General Public License
00017  * along with this program; if not, write to the Free Software
00018  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
00019  */
00020 
00021 /*
00022  * $Author: jim $
00023  * $Date: 2008/07/10 13:05:53 $
00024  * $Revision: 1.33 $
00025  * $Name:  $
00026  */
00027 
00028 /* Includes */
00029 
00030 #ifdef HAVE_CONFIG_H
00031 #include <config.h>
00032 #endif
00033 
00034 #include <cpl.h>
00035 #include <math.h>
00036 #include <string.h>
00037 
00038 #include "vircam_utils.h"
00039 #include "vircam_pfits.h"
00040 #include "vircam_dfs.h"
00041 #include "vircam_mods.h"
00042 #include "vircam_stats.h"
00043 #include "vircam_fits.h"
00044 #include "vircam_tfits.h"
00045 #include "vircam_wcsutils.h"
00046 #include "vircam_paf.h"
00047 
00048 #include "vircam_jmp_utils.h"
00049 
00050 static char *vircam_jmp_outfile(const char *bname, int ind, int isfits);
00051 
00065 /*---------------------------------------------------------------------------*/
00084 /*---------------------------------------------------------------------------*/
00085 
00086 extern int vircam_jmp_save_simple(cpl_frameset *framelist,
00087                                   cpl_parameterlist *parlist) {
00088     cpl_propertylist *plist;
00089     int i,isdummy;
00090     cpl_frame *product_frame;
00091     char *fname;
00092     const char *base[] = {"","simple_jmp","simple_std","simple_mes"};
00093     const char *fctid = "vircam_jmp_save_simple";
00094 
00095     /* Initialise the product frameset */
00096 
00097     if (ps.product_frames_simple == NULL)
00098         ps.product_frames_simple = cpl_malloc(ps.nscience*sizeof(cpl_frame *));
00099 
00100     /* Loop for each of the frames in the ustep sequence */
00101 
00102     for (i = 0; i < ps.nscience; i++) {
00103         fname = vircam_jmp_outfile(base[recflag],i,1);
00104         isdummy = (vircam_fits_get_status(ps.sci_fits[i]) != VIR_OK);
00105 
00106         /* If we need to make a PHU then do that now based on the first frame
00107            in the input frame list */
00108 
00109         if (isfirst) {
00110 
00111             /* Create a new product frame object and define some tags */
00112 
00113             product_frame = cpl_frame_new();
00114             cpl_frame_set_filename(product_frame,fname);
00115             switch (recflag) {
00116             case RECSCI:
00117                 cpl_frame_set_tag(product_frame,VIRCAM_PRO_SIMPLE_SCI);
00118                 break;
00119             case RECSTD:
00120                 cpl_frame_set_tag(product_frame,VIRCAM_PRO_SIMPLE_STD);
00121                 break;
00122             default:
00123                 cpl_frame_set_tag(product_frame,VIRCAM_PRO_SIMPLE_SCI);
00124                 break;
00125             }
00126             cpl_frame_set_type(product_frame,CPL_FRAME_TYPE_IMAGE);
00127             cpl_frame_set_group(product_frame,CPL_FRAME_GROUP_PRODUCT);
00128             cpl_frame_set_level(product_frame,CPL_FRAME_LEVEL_FINAL);
00129 
00130             /* Set up the PHU header */
00131 
00132             plist = vircam_fits_get_phu(ps.sci_fits[i]);
00133             if (i == 0)
00134                 ps.phupaf = vircam_paf_phu_items(plist);;
00135             vircam_dfs_set_product_primary_header(plist,product_frame,
00136                                                   framelist,parlist,
00137                                                   vircam_recipename,
00138                                                   "PRO-1.15");
00139 
00140             /* 'Save' the PHU image */                   
00141 
00142             if (cpl_image_save(NULL,fname,CPL_BPP_8_UNSIGNED,plist,
00143                                CPL_IO_DEFAULT) != CPL_ERROR_NONE) {
00144                 cpl_msg_error(fctid,"Cannot save product PHU");
00145                 cpl_frame_delete(product_frame);
00146                 freespace(fname);
00147                 return(-1);
00148             }
00149             cpl_frameset_insert(framelist,product_frame);
00150             ps.product_frames_simple[i] = product_frame;
00151         }
00152 
00153         /* Get the extension property list */
00154 
00155         plist = vircam_fits_get_ehu(ps.sci_fits[i]);
00156         if (isdummy)
00157             vircam_dummy_property(plist);
00158 
00159         /* Fiddle with the header now */
00160 
00161         product_frame = ps.product_frames_simple[i];
00162         vircam_dfs_set_product_exten_header(plist,product_frame,framelist,
00163                                             parlist,vircam_recipename,
00164                                             "PRO-1.15");
00165         if (cpl_image_save(vircam_fits_get_image(ps.sci_fits[i]),fname,
00166                            CPL_BPP_IEEE_FLOAT,plist,
00167                            CPL_IO_EXTEND) != CPL_ERROR_NONE) {
00168             cpl_msg_error(fctid,"Cannot save product image extension");
00169             freespace(fname);
00170             return(-1);
00171         }
00172         
00173         /* Clear the allocated workspace */
00174 
00175         freespace(fname);
00176     }
00177 
00178     /* Get out of here */
00179 
00180     return(0);
00181 }
00182 
00183 /*---------------------------------------------------------------------------*/
00202 /*---------------------------------------------------------------------------*/
00203 
00204 extern int vircam_jmp_save_super(cpl_frameset *framelist,
00205                                  cpl_parameterlist *parlist) {
00206     cpl_propertylist *plist,*p;
00207     int i,isdummy,isdummyc;
00208     cpl_frame *product_frame;
00209     char *fname,*fnamec;
00210     const char *base[] = {"","super_jmp","super_std","super_mes"};
00211     const char *basec[] = {"","superc_jmp","superc_std","superc_mes"};
00212     vir_fits *ff,*ffc;
00213     cpl_image *fim,*fimc;
00214     const char *fctid = "vircam_jmp_save_super";
00215 
00216     /* Initialise the product frameset */
00217 
00218     if (ps.product_frames_super == NULL)
00219         ps.product_frames_super = cpl_malloc(ps.nustep_sets*sizeof(cpl_frame *));
00220     if (ps.product_frames_superc == NULL)
00221         ps.product_frames_superc = cpl_malloc(ps.nustep_sets*sizeof(cpl_frame *));
00222 
00223     /* Loop for each of the ustep sequences */
00224 
00225     for (i = 0; i < ps.nustep_sets; i++) {
00226 
00227         /* Work out which frame to base the output on. If this particular
00228            microstep sequence failed for whatever reason, there will be
00229            a dummy super frame. */
00230         
00231         ff = ps.ustep_sets[i].super;
00232         fim = vircam_fits_get_image(ff);
00233         ffc = ps.ustep_sets[i].superc;
00234         fimc = vircam_fits_get_image(ffc);
00235         isdummy = (vircam_fits_get_status(ff) != VIR_OK);
00236         isdummyc = (vircam_fits_get_status(ffc) != VIR_OK);
00237 
00238         /* Create some names */
00239 
00240         fname = vircam_jmp_outfile(base[recflag],i,1);
00241         fnamec = vircam_jmp_outfile(basec[recflag],i,1);
00242 
00243         /* If we need to make a PHU then do that now based on the first frame
00244            in the input frame list */
00245 
00246         if (isfirst) {
00247 
00248             /* Create a new product frame object and define some tags */
00249 
00250             product_frame = cpl_frame_new();
00251             cpl_frame_set_filename(product_frame,fname);
00252             switch (recflag) {
00253             case RECSCI:
00254                 cpl_frame_set_tag(product_frame,VIRCAM_PRO_INTER_SCI);
00255                 break;
00256             case RECSTD:
00257                 cpl_frame_set_tag(product_frame,VIRCAM_PRO_INTER_STD);
00258                 break;
00259             default:
00260                 cpl_frame_set_tag(product_frame,VIRCAM_PRO_INTER_SCI);
00261                 break;
00262             }
00263             cpl_frame_set_type(product_frame,CPL_FRAME_TYPE_IMAGE);
00264             cpl_frame_set_group(product_frame,CPL_FRAME_GROUP_PRODUCT);
00265             cpl_frame_set_level(product_frame,CPL_FRAME_LEVEL_FINAL);
00266 
00267             /* Set up the PHU header. */
00268 
00269             plist = vircam_fits_get_phu(ff);
00270             vircam_dfs_set_product_primary_header(plist,product_frame,
00271                                                   framelist,parlist,
00272                                                   vircam_recipename,
00273                                                   "PRO-1.15");
00274 
00275             /* 'Save' the PHU image */                   
00276 
00277             if (cpl_image_save(NULL,fname,CPL_BPP_8_UNSIGNED,plist,
00278                                CPL_IO_DEFAULT) != CPL_ERROR_NONE) {
00279                 cpl_msg_error(fctid,"Cannot save product PHU");
00280                 cpl_frame_delete(product_frame);
00281                 freespace(fname);
00282                 freespace(fnamec);
00283                 return(-1);
00284             }
00285             cpl_frameset_insert(framelist,product_frame);
00286             ps.product_frames_super[i] = product_frame;
00287 
00288             /* Now do exactly the same thing for the confidence map */
00289 
00290             product_frame = cpl_frame_new();
00291             cpl_frame_set_filename(product_frame,fnamec);
00292             switch (recflag) {
00293             case RECSCI:
00294                 cpl_frame_set_tag(product_frame,VIRCAM_PRO_CONF_SCI);
00295                 break;
00296             case RECSTD:
00297                 cpl_frame_set_tag(product_frame,VIRCAM_PRO_CONF_STD);
00298                 break;
00299             default:
00300                 cpl_frame_set_tag(product_frame,VIRCAM_PRO_CONF_SCI);
00301                 break;
00302             }
00303             cpl_frame_set_type(product_frame,CPL_FRAME_TYPE_IMAGE);
00304             cpl_frame_set_group(product_frame,CPL_FRAME_GROUP_PRODUCT);
00305             cpl_frame_set_level(product_frame,CPL_FRAME_LEVEL_FINAL);
00306 
00307             /* Set up the PHU header. */
00308 
00309             plist = vircam_fits_get_phu(ffc);
00310             vircam_dfs_set_product_primary_header(plist,product_frame,
00311                                                   framelist,parlist,
00312                                                   vircam_recipename,
00313                                                   "PRO-1.15");
00314 
00315             /* 'Save' the PHU image */                   
00316 
00317             if (cpl_image_save(NULL,fnamec,CPL_BPP_8_UNSIGNED,plist,
00318                                CPL_IO_DEFAULT) != CPL_ERROR_NONE) {
00319                 cpl_msg_error(fctid,"Cannot save product PHU");
00320                 cpl_frame_delete(product_frame);
00321                 freespace(fname);
00322                 freespace(fnamec);
00323                 return(-1);
00324             }
00325             cpl_frameset_insert(framelist,product_frame);
00326             ps.product_frames_superc[i] = product_frame;
00327         }
00328 
00329 
00330         /* Get the extension property list */
00331 
00332         p = cpl_propertylist_duplicate(vircam_fits_get_ehu(ff));
00333         if (isdummy)
00334             vircam_dummy_property(p);
00335 
00336         /* Fiddle with the header now */
00337 
00338         product_frame = ps.product_frames_super[i];
00339         vircam_dfs_set_product_exten_header(p,product_frame,framelist,parlist,
00340                                             vircam_recipename,"PRO-1.15");
00341         if (cpl_image_save(fim,fname,CPL_BPP_IEEE_FLOAT,p,CPL_IO_EXTEND) != 
00342             CPL_ERROR_NONE) {
00343             cpl_msg_error(fctid,"Cannot save product image extension");
00344             freespace(fname);
00345             freespace(fnamec);
00346             return(-1);
00347         }
00348         cpl_propertylist_delete(p);
00349         
00350         /* And now for the confidence map */
00351 
00352         p = cpl_propertylist_duplicate(vircam_fits_get_ehu(ffc));
00353         if (isdummyc)
00354             vircam_dummy_property(p);
00355 
00356         /* Fiddle with the header now */
00357 
00358         product_frame = ps.product_frames_superc[i];
00359         vircam_dfs_set_product_exten_header(p,product_frame,framelist,parlist,
00360                                             vircam_recipename,"PRO-1.15");
00361         if (cpl_image_save(fimc,fnamec,CPL_BPP_16_SIGNED,p,CPL_IO_EXTEND) != 
00362             CPL_ERROR_NONE) {
00363             cpl_msg_error(fctid,"Cannot save confidence map image extension");
00364             freespace(fname);
00365             freespace(fnamec);
00366             return(-1);
00367         }
00368         cpl_propertylist_delete(p);
00369         
00370         /* Clear the allocated workspace */
00371 
00372         freespace(fname);
00373         freespace(fnamec);
00374     }
00375 
00376     /* Get out of here */
00377 
00378     return(0);
00379 }
00380 
00381 /*---------------------------------------------------------------------------*/
00400 /*---------------------------------------------------------------------------*/
00401 
00402 extern int vircam_jmp_save_stack(cpl_frameset *framelist,
00403                                  cpl_parameterlist *parlist) {
00404     cpl_propertylist *plist,*p,*pafprop;
00405     int isdummy,isdummyc;
00406     vir_fits *ff,*ffc;
00407     cpl_image *fim,*fimc;
00408     cpl_frame *product_frame;
00409     char *fname,*fnamec,*fnamepaf;
00410     const char *base[] = {"","stack_jmp","stack_std","stack_mes"};
00411     const char *basec[] = {"","stackc_jmp","stackc_std","stackc_mes"};
00412     const char *fctid = "vircam_jmp_save_stack";
00413 
00414     /* Work out which frame to base the output on. If this particular
00415        jitter sequence failed for whatever reason, there will be
00416        a dummy stack frame. */
00417 
00418     ff = ps.stack_frame;
00419     fim = vircam_fits_get_image(ff);
00420     ffc = ps.stackc_frame;
00421     fimc = vircam_fits_get_image(ffc);
00422     isdummy = (vircam_fits_get_status(ff) != VIR_OK);
00423     isdummyc = (vircam_fits_get_status(ffc) != VIR_OK);
00424 
00425     /* Create some names */
00426 
00427     fname = vircam_jmp_outfile(base[recflag],0,1);
00428     fnamec = vircam_jmp_outfile(basec[recflag],0,1);
00429     fnamepaf = vircam_jmp_outfile(base[recflag],0,0);
00430 
00431     /* If we need to make a PHU then do that now based on the first frame
00432        in the input frame list */
00433 
00434     if (isfirst) {
00435 
00436         /* Create a new product frame object and define some tags */
00437 
00438         product_frame = cpl_frame_new();
00439         cpl_frame_set_filename(product_frame,fname);
00440         switch (recflag) {
00441         case RECSCI:
00442             cpl_frame_set_tag(product_frame,VIRCAM_PRO_JITTERED_SCI);
00443             break;
00444         case RECSTD:
00445             cpl_frame_set_tag(product_frame,VIRCAM_PRO_JITTERED_STD);
00446             break;
00447         default:
00448             cpl_frame_set_tag(product_frame,VIRCAM_PRO_JITTERED_SCI);
00449             break;
00450         }
00451         cpl_frame_set_type(product_frame,CPL_FRAME_TYPE_IMAGE);
00452         cpl_frame_set_group(product_frame,CPL_FRAME_GROUP_PRODUCT);
00453         cpl_frame_set_level(product_frame,CPL_FRAME_LEVEL_FINAL);
00454         ps.product_frame_stack = product_frame;
00455 
00456         /* Set up the PHU header. */
00457 
00458         plist = vircam_fits_get_phu(ff);
00459         vircam_dfs_set_product_primary_header(plist,product_frame,framelist,
00460                                               parlist,vircam_recipename,
00461                                               "PRO-1.15");
00462 
00463         /* 'Save' the PHU image */                       
00464 
00465         if (cpl_image_save(NULL,fname,CPL_BPP_8_UNSIGNED,plist,
00466                            CPL_IO_DEFAULT) != CPL_ERROR_NONE) {
00467             cpl_msg_error(fctid,"Cannot save product PHU");
00468             cpl_frame_delete(product_frame);
00469             freespace(fname);
00470             freespace(fnamec);
00471             return(-1);
00472         }
00473         cpl_frameset_insert(framelist,product_frame);
00474 
00475         /* Now do exactly the same thing for the confidence map */
00476 
00477         product_frame = cpl_frame_new();
00478         cpl_frame_set_filename(product_frame,fnamec);
00479         switch (recflag) {
00480         case RECSCI:
00481             cpl_frame_set_tag(product_frame,VIRCAM_PRO_CONF_SCI);
00482             break;
00483         case RECSTD:
00484             cpl_frame_set_tag(product_frame,VIRCAM_PRO_CONF_STD);
00485             break;
00486         default:
00487             cpl_frame_set_tag(product_frame,VIRCAM_PRO_CONF_SCI);
00488             break;
00489         }
00490         cpl_frame_set_type(product_frame,CPL_FRAME_TYPE_IMAGE);
00491         cpl_frame_set_group(product_frame,CPL_FRAME_GROUP_PRODUCT);
00492         cpl_frame_set_level(product_frame,CPL_FRAME_LEVEL_FINAL);
00493         ps.product_frame_stackc = product_frame;
00494 
00495         /* Set up the PHU header. */
00496 
00497         plist = vircam_fits_get_phu(ffc);
00498         vircam_dfs_set_product_primary_header(plist,product_frame,framelist,
00499                                               parlist,vircam_recipename,
00500                                               "PRO-1.15");
00501 
00502         /* 'Save' the PHU image */                       
00503 
00504         if (cpl_image_save(NULL,fnamec,CPL_BPP_8_UNSIGNED,plist,
00505                            CPL_IO_DEFAULT) != CPL_ERROR_NONE) {
00506             cpl_msg_error(fctid,"Cannot save product PHU");
00507             cpl_frame_delete(product_frame);
00508             freespace(fname);
00509             freespace(fnamec);
00510             return(-1);
00511         }
00512         cpl_frameset_insert(framelist,product_frame);
00513     }
00514 
00515     /* Get the extension property list */
00516 
00517     p = cpl_propertylist_duplicate(vircam_fits_get_ehu(ff));
00518     if (isdummy) {
00519         vircam_dummy_property(p);
00520         vircam_merge_propertylists(p,dummyqc);
00521     }
00522 
00523     /* Fiddle with the header now */
00524 
00525     product_frame = ps.product_frame_stack;
00526     vircam_dfs_set_product_exten_header(p,product_frame,framelist,parlist,
00527                                         vircam_recipename,"PRO-1.15");
00528     if (cpl_image_save(fim,fname,CPL_BPP_IEEE_FLOAT,p,CPL_IO_EXTEND) != 
00529         CPL_ERROR_NONE) {
00530         cpl_msg_error(fctid,"Cannot save product image extension");
00531         freespace(fname);
00532         freespace(fnamec);
00533         return(-1);
00534     }
00535 
00536     /* Write PAF */
00537 
00538     pafprop = vircam_paf_req_items(p);
00539     vircam_merge_propertylists(pafprop,ps.phupaf);
00540     vircam_paf_append(pafprop,vircam_fits_get_phu(ff),"ESO INS FILT1 NAME");
00541     vircam_paf_append(pafprop,p,"ESO DET NDIT");
00542     vircam_paf_append(pafprop,vircam_fits_get_phu(ff),"RA");
00543     vircam_paf_append(pafprop,vircam_fits_get_phu(ff),"DEC");
00544     vircam_paf_append(pafprop,vircam_fits_get_phu(ff),"ESO TEL AIRM START");
00545     vircam_paf_append(pafprop,vircam_fits_get_phu(ff),"ESO TEL GUID FWHM");
00546     vircam_paf_append(pafprop,vircam_fits_get_phu(ff),"ESO TEL AMBI RHUM");
00547     vircam_paf_append(pafprop,vircam_fits_get_phu(ff),"ESO OBS TARG NAME");
00548     if (vircam_paf_print(fnamepaf,vircam_recipepaf,"QC file",
00549                          pafprop) != VIR_OK)
00550         cpl_msg_warning(fctid,"Unable to save PAF for stack");
00551     cpl_propertylist_delete(pafprop);
00552 
00553     /* Quick tidy */
00554 
00555     cpl_propertylist_delete(p);
00556 
00557     /* And now for the confidence map */
00558 
00559     p = cpl_propertylist_duplicate(vircam_fits_get_ehu(ffc));
00560     if (isdummyc) 
00561         vircam_dummy_property(p);
00562 
00563     /* Fiddle with the header now */
00564 
00565     product_frame = ps.product_frame_stackc;
00566     vircam_dfs_set_product_exten_header(p,product_frame,framelist,parlist,
00567                                         vircam_recipename,"PRO-1.15");
00568     if (cpl_image_save(fimc,fnamec,CPL_BPP_16_SIGNED,p,CPL_IO_EXTEND) != 
00569         CPL_ERROR_NONE) {
00570         cpl_msg_error(fctid,"Cannot save product image extension");
00571         freespace(fname);
00572         freespace(fnamec);
00573         freepropertylist(p);
00574         return(-1);
00575     }
00576     cpl_propertylist_delete(p);
00577 
00578     /* Clear the allocated workspace */
00579 
00580     freespace(fname);
00581     freespace(fnamec);
00582     freespace(fnamepaf);
00583 
00584     /* Get out of here */
00585 
00586     return(0);
00587 }
00588 
00589 /*---------------------------------------------------------------------------*/
00608 /*---------------------------------------------------------------------------*/
00609 
00610 extern int vircam_jmp_save_catalogue(cpl_frameset *framelist,
00611                                      cpl_parameterlist *parlist) {
00612     cpl_frame *product_frame;
00613     int isdummy,i,status;
00614     cpl_table *ftab;
00615     cpl_propertylist *ehu,*phu,*ehu2,*pafprop;
00616     char *fname,*fnamepaf;
00617     const char *base[] = {"","catalogue_jmp","catalogue_std","catalogue_mes"};
00618     const char *fctid = "vircam_jmp_save_catalogue";
00619 
00620     /* Work out whether you have a table to save or not. If you don't
00621        then create a dummy table and save it */
00622 
00623     isdummy = 0;
00624     if (ps.outcat != NULL) {
00625         ftab = vircam_tfits_get_table(ps.outcat);
00626     } else {
00627         ftab = vircam_dummy_catalogue(2);
00628         isdummy = 1;
00629     }
00630 
00631     /* Make an attempt to get some header information */
00632 
00633     if (ps.outcat != NULL) {
00634         ehu = vircam_tfits_get_ehu(ps.outcat);
00635     } else if (ps.stack_frame != NULL) {
00636         ehu = vircam_fits_get_ehu(ps.stack_frame);
00637     } else {
00638         for (i = 0; i < ps.nscience; i++) {
00639             if ((ehu = vircam_fits_get_ehu(ps.sci_fits[i])) != NULL)
00640                 break;
00641         }
00642     }
00643     if (ps.outcat != NULL) {
00644         phu = vircam_tfits_get_phu(ps.outcat);
00645     } else if (ps.stack_frame != NULL) {
00646         phu = vircam_fits_get_phu(ps.stack_frame);
00647     } else {
00648         for (i = 0; i < ps.nscience; i++) {
00649             if ((phu = vircam_fits_get_phu(ps.sci_fits[i])) != NULL)
00650                 break;
00651         }
00652     }
00653 
00654     /* Create a name */
00655 
00656     fname = vircam_jmp_outfile(base[recflag],0,1);
00657     fnamepaf = vircam_jmp_outfile(base[recflag],0,0);
00658 
00659     /* If we need to make a PHU then do that now based on the first frame
00660        in the input frame list */
00661 
00662     if (isfirst) {
00663 
00664         /* Create a new product frame object and define some tags */
00665 
00666         product_frame = cpl_frame_new();
00667         cpl_frame_set_filename(product_frame,fname);
00668         switch (recflag) {
00669         case RECSCI:
00670             cpl_frame_set_tag(product_frame,VIRCAM_PRO_OBJCAT_SCI);
00671             break;
00672         case RECSTD:
00673             cpl_frame_set_tag(product_frame,VIRCAM_PRO_OBJCAT_STD);
00674             break;
00675         default:
00676             cpl_frame_set_tag(product_frame,VIRCAM_PRO_OBJCAT_SCI);
00677             break;
00678         }
00679         cpl_frame_set_type(product_frame,CPL_FRAME_TYPE_TABLE);
00680         cpl_frame_set_group(product_frame,CPL_FRAME_GROUP_PRODUCT);
00681         cpl_frame_set_level(product_frame,CPL_FRAME_LEVEL_FINAL);
00682         ps.product_frame_cat = product_frame;
00683 
00684         /* Set up the PHU header. */
00685 
00686         vircam_dfs_set_product_primary_header(phu,product_frame,framelist,
00687                                               parlist,vircam_recipename,
00688                                               "PRO-1.15");
00689 
00690         /* Copy the header and delete any WCS info */
00691 
00692         ehu2 = cpl_propertylist_duplicate(ehu);
00693         status = VIR_OK;
00694         if (isdummy) {
00695             vircam_dummy_property(ehu2);
00696             vircam_merge_propertylists(ehu2,dummyqc);
00697             (void)vircam_removewcs(ehu2,&status);
00698         }
00699 
00700         /* Setup the extension header */
00701 
00702         vircam_dfs_set_product_exten_header(ehu2,product_frame,framelist,
00703                                             parlist,vircam_recipename,
00704                                             "PRO-1.15");
00705 
00706         /* 'Save' the PHU image */                       
00707 
00708         if (cpl_table_save(ftab,phu,ehu2,fname,CPL_IO_DEFAULT) != 
00709             CPL_ERROR_NONE) {
00710             cpl_msg_error(fctid,"Cannot save product PHU");
00711             cpl_frame_delete(product_frame);
00712             freepropertylist(ehu2);
00713             freespace(fname);
00714             return(-1);
00715         }
00716 
00717         /* Write PAF */
00718 
00719         pafprop = vircam_paf_req_items(ehu2);
00720         vircam_merge_propertylists(pafprop,ps.phupaf);
00721         vircam_paf_append(pafprop,phu,"ESO INS FILT1 NAME");
00722         vircam_paf_append(pafprop,ehu,"ESO DET NDIT");
00723         vircam_paf_append(pafprop,phu,"RA");
00724         vircam_paf_append(pafprop,phu,"DEC");
00725         vircam_paf_append(pafprop,phu,"ESO TEL AIRM START");
00726         vircam_paf_append(pafprop,phu,"ESO TEL GUID FWHM");
00727         vircam_paf_append(pafprop,phu,"ESO TEL AMBI RHUM");
00728         vircam_paf_append(pafprop,phu,"ESO OBS TARG NAME");
00729         if (vircam_paf_print(fnamepaf,vircam_recipepaf,"QC file",
00730                              pafprop) != VIR_OK)
00731             cpl_msg_warning(fctid,"Unable to save PAF for catalogue");
00732         cpl_propertylist_delete(pafprop);
00733 
00734         /* Quick tidy */
00735 
00736         freepropertylist(ehu2);
00737         if (isdummy)
00738             cpl_table_delete(ftab);
00739         cpl_frameset_insert(framelist,product_frame);
00740 
00741     } else {
00742 
00743 
00744         /* Fiddle with the header now */
00745 
00746         product_frame = ps.product_frame_cat;
00747         vircam_dfs_set_product_exten_header(ehu,product_frame,framelist,
00748                                             parlist,vircam_recipename,
00749                                             "PRO-1.15");
00750         
00751         /* Copy the header and delete any WCS info */
00752 
00753         ehu2 = cpl_propertylist_duplicate(ehu);
00754         status = VIR_OK;
00755         if (isdummy) {
00756             vircam_dummy_property(ehu2);
00757             vircam_merge_propertylists(ehu2,dummyqc);
00758             (void)vircam_removewcs(ehu2,&status);
00759         }
00760         
00761         /* Save the table */
00762 
00763         if (cpl_table_save(ftab,NULL,ehu2,fname,CPL_IO_EXTEND) != 
00764             CPL_ERROR_NONE) {
00765             cpl_msg_error(fctid,"Cannot save product image extension");
00766             freepropertylist(ehu2);
00767             freespace(fname);
00768             return(-1);
00769         }
00770 
00771         /* Write PAF */
00772 
00773         pafprop = vircam_paf_req_items(ehu2);
00774         vircam_merge_propertylists(pafprop,ps.phupaf);
00775         vircam_paf_append(pafprop,phu,"ESO INS FILT1 NAME");
00776         vircam_paf_append(pafprop,ehu,"ESO DET NDIT");
00777         vircam_paf_append(pafprop,phu,"RA");
00778         vircam_paf_append(pafprop,phu,"DEC");
00779         vircam_paf_append(pafprop,phu,"ESO TEL AIRM START");
00780         vircam_paf_append(pafprop,phu,"ESO TEL GUID FWHM");
00781         vircam_paf_append(pafprop,phu,"ESO TEL AMBI RHUM");
00782         vircam_paf_append(pafprop,phu,"ESO OBS TARG NAME");
00783         if (vircam_paf_print(fnamepaf,vircam_recipepaf,"QC file",
00784                              pafprop) != VIR_OK)
00785             cpl_msg_warning(fctid,"Unable to save PAF for catalogue");
00786         cpl_propertylist_delete(pafprop);
00787 
00788         /* Quick tidy */
00789 
00790         freepropertylist(ehu2);
00791         if (isdummy)
00792             cpl_table_delete(ftab);
00793     }
00794     freespace(fname);
00795     freespace(fnamepaf);
00796 
00797     /* Get out of here */
00798 
00799     return(0);
00800 }
00801 
00802 /*---------------------------------------------------------------------------*/
00821 /*---------------------------------------------------------------------------*/
00822 
00823 extern int vircam_jmp_save_illum(cpl_frameset *framelist,
00824                                  cpl_parameterlist *parlist) {
00825     cpl_frame *product_frame;
00826     int isdummy,i,status;
00827     cpl_table *ftab;
00828     cpl_propertylist *ehu,*phu,*ehu2,*pafprop;
00829     const char *fname = "illum.fits";
00830     const char *fnamepaf = "illum";
00831     const char *fctid = "vircam_jmp_save_illum";
00832 
00833     /* Work out whether you have a table to save or not. If you don't
00834        then create a dummy table and save it */
00835 
00836     isdummy = 0;
00837     if (ps.illcor != NULL) {
00838         ftab = vircam_tfits_get_table(ps.illcor);
00839     } else {
00840         ftab = vircam_illcor_newtab(1);
00841         isdummy = 1;
00842     }
00843 
00844     /* Make an attempt to get some header information */
00845 
00846     if (ps.illcor != NULL) {
00847         ehu = vircam_tfits_get_ehu(ps.illcor);
00848     } else {
00849         for (i = 0; i < ps.nscience; i++) {
00850             if ((ehu = vircam_fits_get_ehu(ps.sci_fits[i])) != NULL)
00851                 break;
00852         }
00853     }
00854     if (ps.illcor != NULL) {
00855         phu = vircam_tfits_get_phu(ps.illcor);
00856     } else {
00857         for (i = 0; i < ps.nscience; i++) {
00858             if ((phu = vircam_fits_get_phu(ps.sci_fits[i])) != NULL)
00859                 break;
00860         }
00861     }
00862 
00863     /* If we need to make a PHU then do that now based on the first frame
00864        in the input frame list */
00865 
00866     if (isfirst) {
00867 
00868         /* Create a new product frame object and define some tags */
00869 
00870         product_frame = cpl_frame_new();
00871         cpl_frame_set_filename(product_frame,fname);
00872         switch (recflag) {
00873         case RECSTD:
00874             cpl_frame_set_tag(product_frame,VIRCAM_PRO_ILLCOR_STD);
00875             break;
00876         case RECMES:
00877             cpl_frame_set_tag(product_frame,VIRCAM_PRO_ILLCOR_MES);
00878             break;
00879         default:
00880             cpl_frame_set_tag(product_frame,VIRCAM_PRO_ILLCOR_STD);
00881             break;
00882         }
00883         cpl_frame_set_type(product_frame,CPL_FRAME_TYPE_TABLE);
00884         cpl_frame_set_group(product_frame,CPL_FRAME_GROUP_PRODUCT);
00885         cpl_frame_set_level(product_frame,CPL_FRAME_LEVEL_FINAL);
00886         ps.product_frame_illcor = product_frame;
00887 
00888         /* Set up the PHU header. */
00889 
00890         if (ps.phupaf != NULL) 
00891             ps.phupaf = vircam_paf_phu_items(phu);
00892         vircam_dfs_set_product_primary_header(phu,product_frame,framelist,
00893                                               parlist,vircam_recipename,
00894                                               "PRO-1.15");
00895 
00896         /* Copy the header and delete any WCS info */
00897 
00898         ehu2 = cpl_propertylist_duplicate(ehu);
00899         status = VIR_OK;
00900         (void)vircam_removewcs(ehu2,&status);
00901         if (isdummy) {
00902             vircam_dummy_property(ehu2);
00903             vircam_merge_propertylists(ehu2,dummyqc);
00904         }
00905             
00906 
00907         /* Setup the extension header */
00908 
00909         vircam_dfs_set_product_exten_header(ehu2,product_frame,framelist,
00910                                             parlist,vircam_recipename,
00911                                             "PRO-1.15");
00912 
00913         /* 'Save' the PHU image */                       
00914 
00915         if (cpl_table_save(ftab,phu,ehu2,fname,CPL_IO_DEFAULT) != 
00916             CPL_ERROR_NONE) {
00917             cpl_msg_error(fctid,"Cannot save product PHU");
00918             cpl_frame_delete(product_frame);
00919             freepropertylist(ehu2);
00920             return(-1);
00921         }
00922 
00923         /* Write the PAF */
00924 
00925         pafprop = vircam_paf_req_items(ehu2);
00926         vircam_merge_propertylists(pafprop,ps.phupaf);
00927         vircam_paf_append(pafprop,phu,"ESO INS FILT1 NAME");
00928         vircam_paf_append(pafprop,ehu,"ESO DET NDIT");
00929         vircam_paf_append(pafprop,phu,"RA");
00930         vircam_paf_append(pafprop,phu,"DEC");
00931         vircam_paf_append(pafprop,phu,"ESO TEL AIRM START");
00932         vircam_paf_append(pafprop,phu,"ESO TEL GUID FWHM");
00933         vircam_paf_append(pafprop,phu,"ESO TEL AMBI RHUM");
00934         vircam_paf_append(pafprop,phu,"ESO OBS TARG NAME");
00935         if (vircam_paf_print((char *)fnamepaf,vircam_recipepaf,"QC file",
00936                              pafprop) != VIR_OK)
00937             cpl_msg_warning(fctid,"Unable to save PAF for illcor table");
00938         cpl_propertylist_delete(pafprop);
00939         
00940         /* Quick tidy */
00941 
00942         if (isdummy)
00943             cpl_table_delete(ftab);
00944         freepropertylist(ehu2);
00945         cpl_frameset_insert(framelist,product_frame);
00946 
00947     } else {
00948 
00949 
00950         /* Fiddle with the header now */
00951 
00952         product_frame = ps.product_frame_illcor;
00953         vircam_dfs_set_product_exten_header(ehu,product_frame,framelist,
00954                                             parlist,vircam_recipename,
00955                                             "PRO-1.15");
00956         
00957         /* Copy the header and delete any WCS info */
00958 
00959         ehu2 = cpl_propertylist_duplicate(ehu);
00960         status = VIR_OK;
00961         (void)vircam_removewcs(ehu2,&status);
00962         if (isdummy) {
00963             vircam_dummy_property(ehu2);
00964             vircam_merge_propertylists(ehu2,dummyqc);
00965         }
00966         
00967         /* Save the table */
00968 
00969         if (cpl_table_save(ftab,NULL,ehu2,fname,CPL_IO_EXTEND) != 
00970             CPL_ERROR_NONE) {
00971             cpl_msg_error(fctid,"Cannot save product image extension");
00972             freepropertylist(ehu2);
00973             return(-1);
00974         }
00975 
00976         /* Write the PAF */
00977 
00978         pafprop = vircam_paf_req_items(ehu2);
00979         vircam_merge_propertylists(pafprop,ps.phupaf);
00980         vircam_paf_append(pafprop,phu,"ESO INS FILT1 NAME");
00981         vircam_paf_append(pafprop,ehu,"ESO DET NDIT");
00982         vircam_paf_append(pafprop,phu,"RA");
00983         vircam_paf_append(pafprop,phu,"DEC");
00984         vircam_paf_append(pafprop,phu,"ESO TEL AIRM START");
00985         vircam_paf_append(pafprop,phu,"ESO TEL GUID FWHM");
00986         vircam_paf_append(pafprop,phu,"ESO TEL AMBI RHUM");
00987         vircam_paf_append(pafprop,phu,"ESO OBS TARG NAME");
00988         if (vircam_paf_print((char *)fnamepaf,vircam_recipepaf,"QC file",
00989                              pafprop) != VIR_OK)
00990             cpl_msg_warning(fctid,"Unable to save PAF for illcor table");
00991         cpl_propertylist_delete(pafprop);
00992         
00993         /* Quick tidy */
00994 
00995         freepropertylist(ehu2);
00996         if (isdummy)
00997             cpl_table_delete(ftab);
00998     }
00999 
01000     /* Get out of here */
01001 
01002     return(0);
01003 }
01004 
01005 /*---------------------------------------------------------------------------*/
01021 /*---------------------------------------------------------------------------*/
01022 
01023 extern void vircam_jmp_ustep_seq(void) {
01024     int nalloc,i,match,j,ustepnum,nustep;
01025     vir_fits *ff;
01026     cpl_propertylist *plist;
01027     const char *fctid = "vircam_jmp_ustep_seq";
01028 
01029     /* Allocate an initial amount of workspace for the microstep sequence
01030        image sets */
01031 
01032     nalloc = INITALLOC;
01033     ps.ustep_sets = cpl_malloc(nalloc*sizeof(ustep_set));
01034     ps.nustep_sets = 0;
01035 
01036     /* Loop for each frame and get the microstep sequence number from the
01037        primary header */
01038 
01039     for (i = 0; i < ps.nscience; i++) {
01040         ff = ps.sci_fits[i];
01041         plist = vircam_fits_get_phu(ff);
01042         if (vircam_pfits_get_ustepnum(plist,&ustepnum) != VIR_OK) {
01043             cpl_msg_error(fctid,"No microstep number in %s",
01044                           vircam_fits_get_filename(ff));
01045             vircam_fits_set_error(ff,VIR_FATAL);
01046             continue;
01047         }
01048 
01049         /* See if this sequence number matches any of the others we've
01050            already defined. If it does, then simply add this frame into
01051            that sequences frameset */
01052         
01053         match = 0;
01054         for (j = 0; j < ps.nustep_sets; j++) {
01055             if (ustepnum == ps.ustep_sets[j].ustep_number) {
01056                 match = 1;
01057                 ps.ustep_sets[j].f[ps.ustep_sets[j].nframes] = ff;
01058                 ps.ustep_sets[j].nframes += 1;
01059                 if (vircam_fits_get_status(ff) != VIR_FATAL) 
01060                     ps.ustep_sets[j].ngood += 1;
01061                 break;
01062             }
01063         }
01064 
01065         /* If it doesn't match increment the number of sets and check to 
01066            make sure that we haven't overrun our allocation for ustep sets. */
01067 
01068         if (! match) {
01069             if (ps.nustep_sets+1 == nalloc) {
01070                 nalloc += INITALLOC;
01071                 ps.ustep_sets = cpl_realloc(ps.ustep_sets,nalloc*sizeof(ustep_set));
01072             }
01073 
01074             /* Now define this ustep set */
01075             
01076             (void)vircam_pfits_get_nusteps(plist,&nustep);
01077             ps.ustep_sets[ps.nustep_sets].f = cpl_malloc(nustep*sizeof(vir_fits *));
01078             ps.ustep_sets[ps.nustep_sets].ustep_number = ustepnum;
01079             ps.ustep_sets[ps.nustep_sets].nustep = nustep;
01080             ps.ustep_sets[ps.nustep_sets].status = VIR_OK;
01081             ps.ustep_sets[ps.nustep_sets].super = NULL;
01082             ps.ustep_sets[ps.nustep_sets].superc = NULL;
01083             ps.ustep_sets[ps.nustep_sets].nframes = 0;
01084             ps.ustep_sets[ps.nustep_sets].f[0] = ff;
01085             ps.ustep_sets[ps.nustep_sets].nframes = 1;
01086             ps.ustep_sets[ps.nustep_sets].ngood = 0;
01087             if (vircam_fits_get_status(ff) != VIR_FATAL)
01088                 ps.ustep_sets[ps.nustep_sets].ngood += 1;
01089             ps.nustep_sets++;
01090         }    
01091     }
01092 
01093     /* Fix the allocation to what we need and throw the rest away */
01094 
01095     ps.ustep_sets = cpl_realloc(ps.ustep_sets,
01096                                 ps.nustep_sets*sizeof(ustep_set));
01097 
01098     /* Loop through each of the defined sets and see if each is complete */
01099 
01100     for (i = 0; i < ps.nustep_sets; i++) {
01101         if (ps.ustep_sets[i].ngood == 0) {
01102             cpl_msg_warning(fctid,"Microstep sequence %d has no input",
01103                             ps.ustep_sets[i].ustep_number);
01104             ps.ustep_sets[i].status = VIR_FATAL;
01105         } else if (ps.ustep_sets[i].ngood != ps.ustep_sets[i].nustep) {
01106             cpl_msg_warning(fctid,"Microstep sequence %d incomplete",
01107                             ps.ustep_sets[i].ustep_number);
01108             ps.ustep_sets[i].status = VIR_WARN;
01109         }
01110     }
01111 }
01112 
01113 
01114 /*---------------------------------------------------------------------------*/
01130 /*---------------------------------------------------------------------------*/
01131 
01132 extern void vircam_jmp_interleave(void) {
01133     int i,refset,k,nstep,status,nk,*d;
01134     long npts;
01135     float val;
01136     double refx,refy,refra,refdec,x,y;
01137     cpl_image *fi,*outimage,*outconf;
01138     cpl_propertylist *plist;
01139     vir_fits *ff,**tmp;
01140     cpl_wcs *wcs;
01141     cpl_array *dims;
01142     const char *fctid = "vircam_jmp_interleave";
01143 
01144     /* Work out the microstep sequences */
01145 
01146     vircam_jmp_ustep_seq(); 
01147    
01148     /* Get workspace to contain the output superframes */
01149 
01150     ps.dith_input = cpl_malloc(ps.nustep_sets*sizeof(vir_fits*));
01151     ps.dithc_input = cpl_malloc(ps.nustep_sets*sizeof(vir_fits*));
01152             
01153     /* Loop for each of the ustep sets */
01154 
01155     ps.ndith = 0;
01156     for (i = 0; i < ps.nustep_sets; i++) {
01157         if (ps.ustep_sets[i].status == VIR_FATAL) {
01158             outimage = vircam_dummy_image(ps.ustep_sets[i].f[0]);
01159             outconf = vircam_dummy_image(ps.fconf);
01160             ff = vircam_fits_wrap(outimage,ps.ustep_sets[i].f[0],NULL,NULL);
01161             ps.dith_input[ps.ndith] = ff;
01162             vircam_fits_set_error(ff,VIR_FATAL);
01163             ps.ustep_sets[i].super = ff;
01164             ff = vircam_fits_wrap(outconf,ps.ustep_sets[i].f[0],NULL,NULL);
01165             ps.dithc_input[ps.ndith++] = ff;
01166             vircam_fits_set_error(ff,VIR_FATAL);
01167             ps.ustep_sets[i].superc = ff;
01168             continue;
01169         }
01170 
01171         /* Work out the offsets in the sequence from the WCS in the header.
01172            Fail the ustep sequence if any of it's components have an
01173            unreadable WCS */
01174 
01175         refset = 0;
01176         for (k = 0; k < ps.ustep_sets[i].nframes; k++) {
01177             ff = ps.ustep_sets[i].f[k];
01178             if (vircam_fits_get_status(ff) == VIR_FATAL)
01179                 continue;
01180             wcs = cpl_wcs_new_from_propertylist(vircam_fits_get_ehu(ff));
01181             if (wcs == NULL) {
01182                 cpl_msg_error(fctid,"Unable to open WCS structure %s",
01183                               vircam_fits_get_fullname(ff));
01184                 vircam_fits_set_error(ff,VIR_FATAL);
01185                 continue;
01186             }
01187     
01188             /* Get the background value for this image */
01189 
01190             fi = vircam_fits_get_image(ff);
01191             npts = vircam_getnpts(fi);
01192             val = vircam_med(cpl_image_get_data_float(fi),NULL,npts);
01193             cpl_propertylist_update_float(vircam_fits_get_ehu(ff),
01194                                           "ESO DRS BACKMED",val);
01195 
01196             /* If this is the first frame, then set up the reference coords. */
01197     
01198             if (refset == 0) {
01199                 refset = 1;
01200                 dims = cpl_wcs_get_image_dims(wcs);
01201                 d = cpl_array_get_data_int(dims);
01202                 refx = 0.5*(double)d[0];
01203                 refy = 0.5*(double)d[1];
01204                 vircam_xytoradec(wcs,refx,refy,&refra,&refdec);
01205                 cpl_propertylist_update_double(vircam_fits_get_ehu(ff),
01206                                                "ESO DRS XOFFMICRO",0.0);
01207                 cpl_propertylist_update_double(vircam_fits_get_ehu(ff),
01208                                                "ESO DRS YOFFMICRO",0.0);
01209                 cpl_wcs_delete(wcs);
01210                 continue;
01211             }
01212     
01213             /* Take the reference ra and dec and see where that occurs on
01214                the program image in x,y space */
01215     
01216             vircam_radectoxy(wcs,refra,refdec,&x,&y);
01217             x = refx - x;
01218             y = refy - y;
01219             cpl_propertylist_update_double(vircam_fits_get_ehu(ff),
01220                                            "ESO DRS XOFFMICRO",x);
01221             cpl_propertylist_update_double(vircam_fits_get_ehu(ff),
01222                                            "ESO DRS YOFFMICRO",y);
01223             cpl_wcs_delete(wcs);
01224         }
01225 
01226         /* Get a temporary workspace to hold all the good files in the list */
01227 
01228         tmp = cpl_malloc(ps.ustep_sets[i].nframes*sizeof(vir_fits *));
01229         nk = 0;
01230         for (k = 0; k < ps.ustep_sets[i].nframes; k++) {
01231             ff = ps.ustep_sets[i].f[k];
01232             if (vircam_fits_get_status(ff) != VIR_FATAL)
01233                 tmp[nk++] = ff;
01234         }
01235         if (nk < ps.ustep_sets[i].nframes) {
01236             cpl_msg_error(fctid,"A frame in this ustep sequence failed");
01237             ps.ustep_sets[i].status = VIR_WARN;
01238         }
01239 
01240         /* Otherwise interleave them */
01241 
01242         status = VIR_OK;
01243         nstep = (int)sqrt((double)(ps.ustep_sets[i].nustep));
01244         (void)vircam_interleave(tmp,nk,&(ps.fconf),1,nstep,&plist,&outimage,
01245                                 &outconf,&status);
01246         freespace(tmp);
01247         if (status != VIR_OK) {
01248             cpl_msg_error(fctid,"Interleaving failed for ugroup %d extn %d",
01249                           ps.ustep_sets[i].ustep_number,
01250                           vircam_fits_get_nexten(ps.ustep_sets[i].f[0]));
01251             freepropertylist(plist);
01252             freeimage(outimage);
01253             freeimage(outconf);
01254             outimage = vircam_dummy_image(ps.ustep_sets[i].f[0]);
01255             outconf = vircam_dummy_image(ps.fconf);
01256             ff = vircam_fits_wrap(outimage,ps.ustep_sets[i].f[0],NULL,NULL);
01257             vircam_fits_set_error(ff,VIR_FATAL);
01258             ps.ustep_sets[i].super = ff;
01259             ff = vircam_fits_wrap(outconf,ps.ustep_sets[i].f[0],NULL,NULL);
01260             vircam_fits_set_error(ff,VIR_FATAL);
01261             ps.ustep_sets[i].superc = ff;
01262         }
01263 
01264         /* Wrap the output results and store them away into the list of
01265            frames that will be dithered later on. */
01266 
01267         ps.dith_input[ps.ndith] = vircam_fits_wrap(outimage,
01268                                                    ps.ustep_sets[i].f[0],NULL,
01269                                                    plist);
01270         ps.dithc_input[ps.ndith++] = vircam_fits_wrap(outconf,
01271                                                       ps.ustep_sets[i].f[0],
01272                                                       NULL,plist);
01273         ps.ndithc = ps.ndith;
01274         ps.ustep_sets[i].super = ps.dith_input[ps.ndith - 1];
01275         ps.ustep_sets[i].superc = ps.dithc_input[ps.ndithc - 1];
01276         freepropertylist(plist);
01277     }
01278 }
01279 
01280 /*---------------------------------------------------------------------------*/
01300 /*---------------------------------------------------------------------------*/
01301 
01302 extern void vircam_jmp_dither_offsets(void) {
01303     int status,i,ustepnum,nmatch,refset;
01304     float *xoffs,*yoffs,xoff,yoff,filtfwhm;
01305     cpl_wcs *wcsref,*wcs;
01306     vir_fits *ff,*ffc;
01307     vir_tfits *catref,*outcat;
01308     cpl_table *cr,*oc;
01309     const char *fctid = "vircam_jmp_dither_offsets";
01310 
01311     /* Is there anything to dither? If not, then get out of here. NB: We
01312        don't have to check the status of the input dither files as the
01313        calling routine won't have included any files with bad status in
01314        the input list for this routine. */
01315 
01316     if (ps.ndith == 0)
01317         return;
01318 
01319     /* If there is only 1 then set the offsets to zero and get out of here */
01320 
01321     if (ps.ndith == 1) {
01322         cpl_propertylist_update_double(vircam_fits_get_ehu(ps.dith_input[0]),
01323                                        "ESO DRS XOFFDITHER",(double)0.0);
01324         cpl_propertylist_update_double(vircam_fits_get_ehu(ps.dith_input[0]),
01325                                        "ESO DRS YOFFDITHER",(double)0.0);
01326         return;
01327     }
01328 
01329     /* Initialise the status variable and get some workspace for the offsets */
01330 
01331     status = VIR_OK;
01332     xoffs = cpl_malloc(ps.ndith*sizeof(float));
01333     yoffs = cpl_malloc(ps.ndith*sizeof(float));
01334 
01335     /* Loop for all the input files and get the FITS WCS information. */
01336 
01337     refset = 0;
01338     wcsref = NULL;
01339     for (i = 0; i < ps.ndith; i++) {
01340         if (vircam_fits_get_status(ps.dith_input[i]) == VIR_FATAL)
01341             continue;   
01342         wcs = cpl_wcs_new_from_propertylist(vircam_fits_get_ehu(ps.dith_input[i]));
01343         (void)vircam_pfits_get_ustepnum(vircam_fits_get_phu(ps.dith_input[i]),
01344                                         &ustepnum);
01345 
01346         /* If we can't get a WCS for this image, then signal that with
01347            a warning */
01348 
01349         if (wcs == NULL) {
01350             cpl_msg_warning(fctid,"Unable to get WCS for ustep %d",ustepnum);
01351             xoffs[i] = 0.0;
01352             yoffs[i] = 0.0;
01353             vircam_fits_set_error(ps.dith_input[i],VIR_WARN);
01354             continue;
01355         }
01356 
01357         /* If we don't have a reference WCS yet, then make the current frame
01358            WCS the reference and move on */
01359 
01360         if (! refset) {
01361             xoffs[i] = 0.0;
01362             yoffs[i] = 0.0;
01363             refset = 1;
01364             wcsref = wcs;
01365             continue;
01366         }
01367 
01368         /* Right, assuming we're here, then we need to work out the xy 
01369            differences */
01370 
01371         (void)vircam_diffxywcs(wcs,wcsref,&xoff,&yoff,&status);
01372 
01373         /* Did it work? If not the set a warning status for this file */
01374 
01375         if (status != VIR_OK) {
01376             xoffs[i] = 0.0;
01377             yoffs[i] = 0.0;
01378             cpl_msg_warning(fctid,"Unable to WCS difference for %d",ustepnum);
01379         } else {
01380             xoffs[i] = xoff;
01381             yoffs[i] = yoff;
01382         }
01383         cpl_wcs_delete(wcs);
01384     }
01385     if (wcsref != NULL)
01386         cpl_wcs_delete(wcsref);
01387 
01388     /* Now generate a catalogue for each of the input images */
01389 
01390     catref = NULL;
01391     filtfwhm = (interlv ? 3.5 : 2.0);
01392     for (i = 0; i < ps.ndith; i++) {
01393         if (vircam_fits_get_status(ps.dith_input[i]) == VIR_FATAL)
01394             continue;
01395         status = VIR_OK;
01396         ff = ps.dith_input[i];
01397         if (ps.ndithc != 1)
01398             ffc = ps.dithc_input[i];
01399         else
01400             ffc = ps.dithc_input[0];
01401         (void)vircam_pfits_get_ustepnum(vircam_fits_get_phu(ff),&ustepnum);
01402         outcat = NULL;
01403         (void)vircam_imcore(ff,ffc,25,5.0,0,3.5,64,1,filtfwhm,&outcat,
01404                             &status);
01405 
01406         /* If we get a bad status from imcore and this frame has already
01407            failed in the WCS stage, then mark it with VIR_FATAL status
01408            and move to the next one */
01409 
01410         if (status != VIR_OK && vircam_fits_get_status(ff) != VIR_OK) {
01411             vircam_fits_set_error(ff,VIR_FATAL);
01412             freetfits(outcat);
01413             cpl_msg_error(fctid,"Unable to get offsets for %d",ustepnum);
01414             continue;
01415 
01416         /* If we get bad status, but this file has a perfectly good offset
01417            from the WCS then just go with it that and issue a warning */
01418 
01419         } else if (status != VIR_OK && vircam_fits_get_status(ff) == VIR_OK) {
01420             vircam_fits_set_error(ff,VIR_WARN);
01421             freetfits(outcat);
01422             cpl_msg_error(fctid,"Unable to get object offset for %d. Going with WCS offset",
01423                           ustepnum);
01424             continue;
01425 
01426         /*  OK, we got a good status from imcore. See if the reference 
01427             catalogue has already been defined. If it hasn't then this 
01428             catalogue becomes the reference */
01429 
01430         } else {
01431             oc = vircam_tfits_get_table(outcat);
01432             cpl_table_add_scalar(oc,"X_coordinate",(double)xoffs[i]);
01433             cpl_table_add_scalar(oc,"Y_coordinate",(double)yoffs[i]);
01434             if (catref == NULL) {
01435                 catref = outcat;
01436                 cr = oc;
01437                 continue;
01438             }
01439 
01440             /* If this isn't the reference file, then do the cross match */
01441 
01442 
01443             (void)vircam_matchxy(oc,cr,10.0,&xoff,&yoff,&nmatch,&status);
01444             freetfits(outcat);
01445     
01446             /* If we got a bad result with the cross match and this current
01447                file is already got a warning status, then flag it as fatal
01448                status */
01449     
01450             if ((nmatch == 0 || status == VIR_FATAL) && 
01451                 vircam_fits_get_status(ff) != VIR_OK) {
01452                 xoff = 0.0;
01453                 yoff = 0.0;
01454                 vircam_fits_set_error(ff,VIR_FATAL);
01455                 cpl_msg_error(fctid,"Unable to match stars for %d",ustepnum);
01456     
01457             /* If we got a bad result with the cross match and this current
01458                file has a good status, then just go with the WCS offset
01459                we already have */
01460     
01461             } else if ((nmatch == 0 || status == VIR_FATAL) && 
01462                        vircam_fits_get_status(ff) == VIR_OK) {
01463                 vircam_fits_set_error(ff,VIR_WARN);
01464                 xoff = 0.0;
01465                 yoff = 0.0;
01466                 cpl_msg_warning(fctid,"Unable to match stars for %d. Going with WCS offsets",
01467                                 ustepnum);
01468             }
01469     
01470             /* Add the current offsets to what we already have */
01471     
01472             xoffs[i] += xoff;
01473             yoffs[i] += yoff;
01474         }
01475     }
01476 
01477     /* Write the results to the headers */
01478 
01479     for (i = 0; i < ps.ndith; i++) {
01480         if (vircam_fits_get_status(ps.dith_input[i]) == VIR_FATAL)
01481             continue;
01482         ff = ps.dith_input[i];
01483         cpl_propertylist_update_double(vircam_fits_get_ehu(ff),
01484                                        "ESO DRS XOFFDITHER",(double)xoffs[i]);
01485         cpl_propertylist_update_double(vircam_fits_get_ehu(ff),
01486                                        "ESO DRS YOFFDITHER",(double)yoffs[i]);
01487     }
01488 
01489     /* Get rid of some workspace */
01490 
01491     freespace(xoffs);
01492     freespace(yoffs);
01493     freetfits(catref);
01494 }               
01495 
01496 /*---------------------------------------------------------------------------*/
01512 /*---------------------------------------------------------------------------*/
01513 
01514 extern void vircam_jmp_dither_images(void) {
01515     int status,ngood,i,n;
01516     vir_fits **d,**dc,*ff;
01517     cpl_propertylist *dither_ehu,*dither_phu;
01518     cpl_image *outdither,*outditherc;
01519     const char *fctid = "vircam_jmp_dither_images";
01520 
01521     /* Count how many good input images there are. */
01522 
01523     ngood = 0;
01524     for (i = 0; i < ps.ndith; i++)
01525         if (vircam_fits_get_status(ps.dith_input[i]) != VIR_FATAL)
01526             ngood++;
01527     
01528     /* If there are none, then get out of here. The output jittered image
01529        and it's confidence map should still be set to NULL */
01530 
01531     if (ngood == 0) {
01532         cpl_msg_error(fctid,"No good input images for jittering");
01533 
01534 
01535         outdither = vircam_dummy_image(ps.sci_fits[0]);
01536         outditherc = vircam_dummy_image(ps.dithc_input[0]);
01537         ff = vircam_fits_wrap(outdither,ps.sci_fits[0],NULL,NULL);
01538         ps.stack_frame = ff;
01539         vircam_fits_set_error(ff,VIR_FATAL);
01540         ff = vircam_fits_wrap(outditherc,ps.dithc_input[0],NULL,NULL);
01541         vircam_fits_set_error(ff,VIR_FATAL);
01542         ps.stackc_frame = ff;
01543         return;
01544     }
01545      
01546     /* Now transfer the input images over to new arrays... */
01547 
01548     d = cpl_malloc(ngood*sizeof(vir_fits *));
01549     dc = cpl_malloc(ngood*sizeof(vir_fits *));
01550     n = 0;
01551     for (i = 0; i < ps.ndith; i++) {
01552         if (vircam_fits_get_status(ps.dith_input[i]) != VIR_FATAL) {
01553             d[n] = ps.dith_input[i];
01554             if (ps.ndithc != 1) 
01555                 dc[n++] = ps.dithc_input[i];
01556             else 
01557                 dc[n++] = ps.dithc_input[0];
01558         }
01559     }
01560 
01561     /* Dither the remaining images */
01562 
01563     status = VIR_OK;
01564     (void)vircam_imdither(d,dc,ngood,ngood,5.0,5.0,&dither_ehu,&outdither,
01565                           &outditherc,&status);
01566     dither_phu = cpl_propertylist_duplicate(vircam_fits_get_phu(d[0]));
01567     if (status != VIR_OK) {
01568         freeimage(outdither);
01569         freeimage(outditherc);
01570         cpl_msg_error(fctid,"Error jittering to output");
01571     } else {
01572         ps.stack_frame = vircam_fits_wrap(outdither,d[0],dither_phu,
01573                                           dither_ehu);
01574         ps.stackc_frame = vircam_fits_wrap(outditherc,dc[0],dither_phu,
01575                                            dither_ehu);
01576     }
01577 
01578     /* Tidy and exit */
01579 
01580     freepropertylist(dither_phu);
01581     freepropertylist(dither_ehu);
01582     freespace(d);
01583     freespace(dc);
01584 }
01585 
01586 /*---------------------------------------------------------------------------*/
01601 /*---------------------------------------------------------------------------*/
01602 
01603 extern void vircam_jmp_catalogue(void) {
01604     const char *fctid = "vircam_jmp_catalogue";
01605     int status;
01606     vir_tfits *outtab;
01607     float filtfwhm;
01608 
01609     /* Check to see if there has been a jitter frame defined */
01610 
01611     if (ps.stack_frame == NULL || vircam_fits_get_status(ps.stack_frame) == VIR_FATAL) {
01612         cpl_msg_error(fctid,"No stack image available. No catalogue generated");
01613         return;
01614     }
01615 
01616     /* Generate the catalogue */
01617 
01618     status = VIR_OK;
01619     filtfwhm = (interlv ? 3.5 : 2);
01620     (void)vircam_imcore(ps.stack_frame,ps.stackc_frame,
01621                         vircam_jmp_config.ipix,
01622                         vircam_jmp_config.threshold,
01623                         vircam_jmp_config.icrowd,
01624                         vircam_jmp_config.rcore,
01625                         vircam_jmp_config.nbsize,2,filtfwhm,
01626                         &outtab,&status);
01627 
01628     /* If it failed, then get rid of workspaces and get out of here with
01629        an error message. If it went well, then wrap the result in a
01630        vir_tfits object */
01631 
01632     if (status != VIR_OK) {
01633         cpl_msg_error(fctid,"Error generating catalogue");
01634         freetfits(outtab);
01635         vircam_fits_set_error(ps.stack_frame,VIR_FATAL);
01636     } else {
01637         ps.outcat = outtab;
01638     }
01639                           
01640     return;
01641 }
01642 
01643 /*---------------------------------------------------------------------------*/
01659 /*---------------------------------------------------------------------------*/
01660     
01661 extern void vircam_jmp_matched_stds(void) {
01662     int status;
01663     const char *fctid = "vircam_jmp_matched_stds";
01664     cpl_table *stdscat;
01665 
01666     /* Initialise status */
01667 
01668     status = VIR_OK;
01669 
01670     /* Check that we have a catalogue */
01671     
01672     if (ps.outcat == NULL) {
01673         cpl_msg_error(fctid,"No input catalogue found");
01674         return;
01675     }
01676 
01677     /* Get some standard stars */
01678 
01679     (void)vircam_getstds(vircam_fits_get_ehu(ps.stack_frame),1,
01680                          current_catpath,current_cat,&stdscat,&status);
01681     if (status != VIR_OK) {
01682         freetable(stdscat);
01683         cpl_msg_error(fctid,"Failed to find any standards");
01684         return;
01685     }
01686 
01687     /* Now match this against the catalogue */
01688 
01689     (void)vircam_matchstds(vircam_tfits_get_table(ps.outcat),stdscat,
01690                            10.0,&(ps.matchstds),&status);
01691     freetable(stdscat);
01692     if (status != VIR_OK) {
01693         freetable(ps.matchstds);
01694         cpl_msg_error(fctid,"Failed to match standards to catalogue");
01695         return;
01696     }
01697 }
01698 
01699 /*---------------------------------------------------------------------------*/
01716 /*---------------------------------------------------------------------------*/
01717     
01718 extern void vircam_jmp_wcsfit(void) {
01719     int status,n,i;
01720     const char *fctid = "vircam_jmp_wcsfit";
01721     float *ra,*dec,*x,*y;
01722     double r,d;
01723     cpl_table *cat;
01724     cpl_wcs *wcs;
01725 
01726     /* Initialise status */
01727 
01728     status = VIR_OK;
01729 
01730     /* Check that we have a catalogue */
01731     
01732     if (ps.matchstds == NULL) {
01733         cpl_msg_error(fctid,"No input matched standards catalogue found");
01734         return;
01735     }
01736 
01737     /* Fit the plate solution */
01738 
01739     (void)vircam_platesol(vircam_fits_get_ehu(ps.stack_frame),
01740                           vircam_tfits_get_ehu(ps.outcat),ps.matchstds,
01741                           6,1,&status);
01742     if (status != VIR_OK) {
01743         cpl_msg_error(fctid,"Failed to fit WCS");
01744         return;
01745     }
01746 
01747     /* Update the RA and DEC of the objects in the object catalogue */
01748 
01749     cat = vircam_tfits_get_table(ps.outcat);
01750     n = cpl_table_get_nrow(cat);
01751     wcs = cpl_wcs_new_from_propertylist(vircam_fits_get_ehu(ps.stack_frame));
01752     if (wcs == NULL) {
01753         cpl_msg_error(fctid,"Failed to fill RA and Dec in catalogue");
01754         return;
01755     }
01756     x = cpl_table_get_data_float(cat,"X_coordinate"); 
01757     y = cpl_table_get_data_float(cat,"Y_coordinate");
01758     ra = cpl_table_get_data_float(cat,"RA");
01759     dec = cpl_table_get_data_float(cat,"DEC");
01760     for (i = 0; i < n; i++) {
01761         vircam_xytoradec(wcs,x[i],y[i],&r,&d);
01762         ra[i] = (float)r;
01763         dec[i] = (float)d;
01764     }
01765     cpl_wcs_delete(wcs);
01766 }
01767 
01768 /*---------------------------------------------------------------------------*/
01784 /*---------------------------------------------------------------------------*/
01785 
01786 extern void vircam_jmp_photcal(void) {
01787     int status;
01788     const char *fctid = "vircam_jmp_photcal";
01789     char filt[32];
01790     cpl_propertylist *pl;
01791 
01792     /* Initialise status */
01793 
01794     status = VIR_OK;
01795 
01796     /* Check that we have a catalogue */
01797     
01798     if (ps.matchstds == NULL || cpl_table_get_nrow(ps.matchstds) == 0) {
01799         cpl_msg_error(fctid,"No input matched standards catalogue found");
01800         return;
01801     }
01802 
01803     /* What filter is this? */
01804 
01805     if (vircam_pfits_get_filter(vircam_fits_get_phu(ps.stack_frame),filt) != VIR_OK) {
01806         cpl_msg_error(fctid,"No filter name in stack header");
01807         return;
01808     }
01809 
01810     /* Fit the photometric calibration */
01811     
01812     pl = vircam_tfits_get_ehu(ps.outcat);
01813     (void)vircam_photcal(&(ps.stack_frame),&(ps.matchstds),&pl,1,filt,
01814                          ps.tphottab,&status);
01815     if (status != VIR_OK) {
01816         cpl_msg_error(fctid,"Failed to fit photometric zeropoint");
01817         return;
01818     }
01819 }
01820 
01821 /*---------------------------------------------------------------------------*/
01838 /*---------------------------------------------------------------------------*/
01839 
01840 extern void vircam_jmp_bpm2conf(void) {
01841     cpl_image *im;
01842     int i,n,*data;
01843 
01844     /* Get the image */
01845 
01846     im = vircam_fits_get_image(ps.fconf);
01847     n = cpl_image_get_size_x(im)*cpl_image_get_size_y(im);
01848     data = cpl_image_get_data_int(im);
01849 
01850     /* Convert it now */ 
01851 
01852     for (i = 0; i < n; i++)
01853         data[i] = (data[i] == 1 ? 0 : 100);
01854         
01855 }
01856 
01857 /*---------------------------------------------------------------------------*/
01874 /*---------------------------------------------------------------------------*/
01875 
01876 extern void vircam_jmp_skycor(void) {
01877     int i,ngood,status;
01878     long npts;
01879     vir_fits **ftmp,*ff;
01880     const char *fctid = "vircam_jmp_skycor";
01881     unsigned char *rejmask,*rejplus;
01882     cpl_propertylist *drs;
01883     cpl_image *skyimg,*fim;
01884     float *data,med,sig;
01885 
01886     /* Sort out all images with good status */
01887 
01888     ftmp = cpl_malloc(ps.nscience*sizeof(vir_fits *));
01889     ngood = 0;
01890     for (i = 0; i < ps.nscience; i++) {
01891         ff = ps.sci_fits[i];
01892         if (vircam_fits_get_status(ff) != VIR_FATAL)
01893             ftmp[ngood++] = ff;
01894     }
01895 
01896     /* If there aren't any good images, then get out of here now */
01897 
01898     if (ngood == 0) {
01899         freespace(ftmp);
01900         cpl_msg_error(fctid,"Sky correction impossible. No good science frames available");
01901         return;
01902     }
01903 
01904     /* Combine all the good science images */
01905 
01906     status = VIR_OK;
01907     (void)vircam_imcombine(ftmp,ngood,2,1,1,3.0,&skyimg,&rejmask,
01908                            &rejplus,&drs,&status);
01909     freespace(rejmask);
01910     freespace(rejplus);
01911     freepropertylist(drs);
01912     freespace(ftmp);
01913 
01914     /* Normalise the sky frame to zero median */
01915 
01916     data = cpl_image_get_data_float(skyimg);
01917     npts = cpl_image_get_size_x(skyimg)*cpl_image_get_size_y(skyimg);
01918     vircam_qmedsig(data,NULL,npts,5.0,1,-1000.0,65535.0,&med,&sig);
01919     for (i = 0; i < npts; i++)
01920         data[i] -= med;
01921 
01922     /* Subtract the normalised sky frame from the science images */
01923 
01924     for (i = 0; i < ps.nscience; i++) {
01925         ff = ps.sci_fits[i];
01926         fim = vircam_fits_get_image(ff);
01927         if (vircam_fits_get_status(ff) != VIR_FATAL) {
01928             drs = vircam_fits_get_ehu(ff);
01929             cpl_image_subtract(fim,skyimg);         
01930             cpl_propertylist_update_bool(drs,"ESO DRS SKYCOR",TRUE);
01931             cpl_propertylist_set_comment(drs,"ESO DRS SKYCOR",
01932                                          "Image has been sky corrected");
01933         }
01934     }
01935 
01936     /* Clean up and get out of here */
01937 
01938     cpl_image_delete(skyimg);
01939     return;
01940 }
01941 
01942 /*---------------------------------------------------------------------------*/
01958 /*---------------------------------------------------------------------------*/
01959 
01960 extern void vircam_jmp_get_readnoise_gain(int jext, float *readnoise, 
01961                                           float *gain) {
01962     cpl_propertylist *p_rg;
01963     const char *fctid = "vircam_jmp_get_readnoise_gain";
01964 
01965     /* Load the propertylist */
01966 
01967     p_rg = cpl_propertylist_load(cpl_frame_get_filename(ps.readgain_file),
01968                                  jext);
01969 
01970     /* Check the readnoise property type and read it */
01971 
01972     switch (cpl_propertylist_get_type(p_rg,"ESO QC READNOISE")) {
01973     case CPL_TYPE_FLOAT:
01974         *readnoise = cpl_propertylist_get_float(p_rg,"ESO QC READNOISE");
01975         break;
01976     case CPL_TYPE_DOUBLE:
01977         *readnoise = (float)cpl_propertylist_get_double(p_rg,
01978                                                         "ESO QC READNOISE");
01979         break;
01980     default:
01981         cpl_error_reset();
01982         *readnoise = 25.0;
01983         cpl_msg_error(fctid,"Unable to get READNOISE estimate, guessing %g\n",
01984                       *readnoise);
01985     }
01986 
01987     /* Now the gain */
01988 
01989     switch (cpl_propertylist_get_type(p_rg,"ESO QC GAIN")) {
01990     case CPL_TYPE_FLOAT:
01991         *gain = cpl_propertylist_get_float(p_rg,"ESO QC GAIN");
01992         break;
01993     case CPL_TYPE_DOUBLE:
01994         *gain = (float)cpl_propertylist_get_double(p_rg,"ESO QC GAIN");
01995         break;
01996     default:
01997         cpl_error_reset();
01998         *gain = 1.0;
01999         cpl_msg_error(fctid,"Unable to get GAIN estimate, guessing %g\n",
02000                       *gain);
02001     }
02002     cpl_propertylist_delete(p_rg);
02003 }
02004 
02005 /*---------------------------------------------------------------------------*/
02020 /*---------------------------------------------------------------------------*/
02021 
02022 extern void vircam_jmp_illum(void) {
02023     int ngood,i,status,ii;
02024     float illcor_rms;
02025     vir_fits **ftmp,*ff,*ffc;
02026     char filt[32];
02027     cpl_table **mstds,*stdscat,*ms,*illcor,*ot;
02028     vir_tfits *outtab;
02029     cpl_propertylist **pl,*phu,*ehu;
02030     const char *fctid = "vircam_jmp_illum";
02031 
02032     /* Set some default values */
02033 
02034     ps.illcor = NULL;
02035     
02036     /* Sort out all images with good status */
02037 
02038     ftmp = cpl_malloc(ps.nscience*sizeof(vir_fits *));
02039     ngood = 0;
02040     for (i = 0; i < ps.nscience; i++) {
02041         ff = ps.sci_fits[i];
02042         if (vircam_fits_get_status(ff) != VIR_FATAL)
02043             ftmp[ngood++] = ff;
02044     }
02045 
02046     /* If there aren't any good images, then get out of here now */
02047 
02048     if (ngood == 0) {
02049         freespace(ftmp);
02050         cpl_msg_error(fctid,"Illumination correction impossible. No good science frames available");
02051         return;
02052     }
02053 
02054     /* What filter is this? */
02055 
02056     if (vircam_pfits_get_filter(vircam_fits_get_phu(ftmp[0]),filt) != VIR_OK) {
02057         cpl_msg_error(fctid,"No filter name in stack header");
02058         freespace(ftmp);
02059         return;
02060     }
02061 
02062     /* Get some workspace for the various arrays you need for the 
02063        illumination correction routine */
02064 
02065     mstds = cpl_malloc(ngood*sizeof(cpl_table *));
02066     for (i = 0; i < ngood; i++)
02067         mstds[i] = NULL;
02068     pl = cpl_malloc(ngood*sizeof(cpl_propertylist *));
02069     for (i = 0; i < ngood; i++)
02070         pl[i] = NULL;
02071 
02072     /* For each of the good input frames, do a catalogue generation and
02073        get some matched standards */
02074 
02075     ffc = ps.fconf;
02076     for (i = 0; i < ngood; i++) {
02077         status = VIR_OK;
02078         ff = ftmp[i];
02079         (void)vircam_imcore(ff,ffc,vircam_jmp_config.ipix,
02080                             1.5*vircam_jmp_config.threshold,0,
02081                             vircam_jmp_config.rcore,vircam_jmp_config.nbsize,2,
02082                             3.5,&outtab,&status);
02083         pl[i] = cpl_propertylist_duplicate(vircam_tfits_get_ehu(outtab));
02084 
02085         /* Get some standard stars */
02086 
02087         (void)vircam_getstds(vircam_fits_get_ehu(ff),1,current_catpath,
02088                              current_cat,&stdscat,&status);
02089         if (status == VIR_FATAL) {
02090             freetfits(outtab);
02091             freespace(ftmp);
02092             for (ii = 0; ii < ngood; ii++) {
02093                 freetable(mstds[ii]);
02094                 freepropertylist(pl[ii]);
02095             }
02096             freespace(mstds);
02097             freespace(pl);
02098             cpl_msg_error(fctid,"Illumination correction fails");
02099             return;
02100         } else if (status == VIR_WARN) {
02101             freetfits(outtab);
02102             continue;
02103         }
02104 
02105         /* Now match this against the catalogue */
02106 
02107         ot = vircam_tfits_get_table(outtab);
02108         (void)vircam_matchstds(ot,stdscat,10.0,&ms,&status);
02109         if (status == VIR_FATAL) {
02110             freetable(stdscat);
02111             freetfits(outtab);
02112             freespace(ftmp);
02113             for (ii = 0; ii < ngood; ii++) {
02114                 freetable(mstds[ii]);
02115                 freepropertylist(pl[ii]);
02116             }
02117             freespace(mstds);
02118             freespace(pl);
02119             cpl_msg_error(fctid,"%s",cpl_error_get_message());
02120             return;
02121         }
02122         mstds[i] = ms;
02123         freetfits(outtab);
02124         freetable(stdscat);
02125     }
02126     
02127     /* Call the illumination routine */
02128 
02129     status = VIR_OK;
02130     (void)vircam_illum(ftmp,mstds,pl,ngood,filt,ps.tphottab,128,&illcor,
02131                        &illcor_rms,&status);
02132     
02133     /* Wrap the result */
02134     
02135     phu = cpl_propertylist_duplicate(vircam_fits_get_phu(ftmp[0]));
02136     ehu = cpl_propertylist_duplicate(vircam_fits_get_ehu(ftmp[0]));
02137     ps.illcor = vircam_tfits_wrap(illcor,NULL,phu,ehu);
02138     cpl_propertylist_update_float(ehu,"ESO QC ILLUMCOR_RMS",illcor_rms);
02139     cpl_propertylist_set_comment(ehu,"ESO QC ILLUMCOR_RMS",
02140                                  "RMS of illumination correction map");
02141     
02142     /* Tidy up */
02143 
02144     for (i = 0; i < ngood; i++) {
02145         freetable(mstds[i]);
02146         freepropertylist(pl[i]);
02147     }
02148     freespace(mstds);
02149     freespace(pl);
02150     freespace(ftmp);
02151 }
02152 
02153 /*---------------------------------------------------------------------------*/
02175 /*---------------------------------------------------------------------------*/
02176 
02177 static char *vircam_jmp_outfile(const char *bname, int ind, int isfits) {
02178     int nf;
02179     char *fname;
02180 
02181     /* Count up how much space you need for the output string. The 'isfits'
02182        block accounts for ".fits", an underscore and an EOS*/
02183 
02184     nf = strlen(bname);
02185     if (ind == 0) 
02186         nf++;
02187     else 
02188         nf += ((int)log10((double)ind)+1);
02189     if (isfits) 
02190         nf += 7;
02191     else
02192         nf += 2;
02193 
02194     /* Get the space for the filename */
02195 
02196     fname = cpl_malloc(nf);
02197 
02198     /* Now write the name */
02199 
02200     if (isfits) 
02201         (void)snprintf(fname,nf,"%s_%d.fits",bname,ind);
02202     else
02203         (void)snprintf(fname,nf,"%s_%d",bname,ind);
02204     return(fname);
02205 }
02206 
02207 /*---------------------------------------------------------------------------*/
02222 /*---------------------------------------------------------------------------*/
02223 
02224 extern void vircam_jmp_init(void) {
02225 
02226     /* Level 0 stuff */
02227 
02228     ps.labels = NULL;
02229     ps.master_dark = NULL;
02230     ps.master_twilight_flat = NULL;
02231     ps.master_conf = NULL;
02232     ps.mask = NULL;
02233     ps.chantab = NULL;
02234     ps.phottab = NULL;
02235     ps.tphottab = NULL;
02236     ps.readgain_file = NULL;
02237     ps.science_frames = NULL;
02238     ps.product_frames_simple = NULL;
02239     ps.product_frames_super = NULL;
02240     ps.product_frames_superc = NULL;
02241     ps.product_frame_stack = NULL;
02242     ps.product_frame_stackc = NULL;
02243     ps.product_frame_cat = NULL;
02244     ps.product_frame_illcor = NULL;
02245     ps.phupaf = NULL;
02246     ps.gaincors = NULL;
02247     ps.catpath = NULL;
02248     ps.catname = NULL;
02249     ps.catpath2 = NULL;
02250     ps.catname2 = NULL;
02251 
02252     /* Level 1 stuff */
02253 
02254     ps.fdark = NULL;
02255     ps.fflat = NULL;
02256     ps.fconf = NULL;
02257     ps.fchantab = NULL;
02258     ps.nscience = 0;
02259     ps.sci_fits = NULL;
02260     ps.nustep_sets = 0;
02261     ps.ustep_sets = NULL;
02262     ps.ndith = 0;
02263     ps.ndithc = 0;
02264     ps.dith_input = NULL;
02265     ps.dithc_input = NULL;
02266     ps.stack_frame = NULL;
02267     ps.stackc_frame = NULL;
02268     ps.outcat = NULL;
02269     ps.illcor = NULL;
02270 }
02271 
02272 /*---------------------------------------------------------------------------*/
02293 /*---------------------------------------------------------------------------*/
02294 
02295 extern void vircam_jmp_tidy(int level) {
02296     int i;
02297 
02298     /* Level 1 stuff */
02299 
02300     freefits(ps.fdark); 
02301     freefits(ps.fflat); 
02302     freefits(ps.fconf);
02303     freetfits(ps.fchantab);
02304     freefitslist(ps.sci_fits,ps.nscience);
02305     ps.nscience = 0;
02306     for (i = 0; i < ps.nustep_sets; i++) {
02307         freespace(ps.ustep_sets[i].f);
02308         freefits(ps.ustep_sets[i].super);
02309         freefits(ps.ustep_sets[i].superc);
02310     }
02311     freespace(ps.ustep_sets);
02312     ps.nustep_sets = 0;
02313     freespace(ps.dith_input);
02314     ps.ndith = 0;
02315     freespace(ps.dithc_input);
02316     ps.ndithc = 0;
02317 
02318     freefits(ps.stack_frame);
02319     freefits(ps.stackc_frame);
02320     freetfits(ps.outcat);
02321     freetable(ps.matchstds);
02322     freetfits(ps.illcor);
02323 
02324     if (level == 1)
02325         return;
02326     
02327     /* Level 0 stuff */
02328 
02329     freespace(ps.labels);
02330     freeframe(ps.master_dark);
02331     freeframe(ps.master_twilight_flat);
02332     freeframe(ps.master_conf);
02333     freemask(ps.mask);
02334     freeframe(ps.chantab);
02335     freeframe(ps.phottab);
02336     freeframe(ps.readgain_file);
02337     freetable(ps.tphottab);
02338     freeframeset(ps.science_frames);
02339     freepropertylist(ps.phupaf);
02340     freespace(ps.product_frames_simple); /* NB: We only have to delete the */
02341     freespace(ps.product_frames_super);  /* arrays and not the frames */
02342     freespace(ps.product_frames_superc); /* as these get passed back to esorex */
02343     freespace(ps.gaincors);
02344     freespace(ps.catpath);
02345     freespace(ps.catname);
02346     freespace(ps.catpath2);
02347     freespace(ps.catname2);
02348 }
02349 
02352 /*
02353 
02354 $Log: vircam_jmp_utils.c,v $
02355 Revision 1.33  2008/07/10 13:05:53  jim
02356 Modified to use v4.2 version of cpl_wcs
02357 
02358 Revision 1.32  2008/06/20 11:13:35  jim
02359 Fixed dodgy call to cpl_wcs_get_image_dims
02360 
02361 Revision 1.31  2008/05/06 08:40:10  jim
02362 Modified to use cpl_wcs interface
02363 
02364 Revision 1.30  2007/11/22 12:34:08  jim
02365 Added line to vircam_jmp_save_illum to make sure that phupaf is defined
02366 when _save_simple isn't called
02367 
02368 Revision 1.29  2007/10/25 17:34:00  jim
02369 Modified to remove lint warnings
02370 
02371 Revision 1.28  2007/10/19 09:25:10  jim
02372 Fixed problems with missing includes
02373 
02374 Revision 1.27  2007/10/19 06:55:06  jim
02375 Modifications made to use new method for directing the recipes to the
02376 standard catalogues using the sof
02377 
02378 Revision 1.26  2007/10/15 12:50:28  jim
02379 Modified for compatibility with cpl_4.0
02380 
02381 Revision 1.25  2007/06/13 08:10:49  jim
02382 Modified to allow for different output file names depending upon the calling
02383 recipe
02384 
02385 Revision 1.24  2007/05/15 08:54:07  jim
02386 Fixed small bug in stack_save and dither_offsets. Also modified dither_offsets
02387 to just send back zeros if there is only one dither frame
02388 
02389 Revision 1.23  2007/05/08 10:41:49  jim
02390 Added gaincor variables
02391 
02392 Revision 1.22  2007/05/02 09:15:37  jim
02393 Modified to use new api for vircam_imcore and vircam_platesol
02394 
02395 Revision 1.21  2007/04/30 09:40:01  jim
02396 Added vircam_paf_append
02397 
02398 Revision 1.20  2007/04/13 12:29:11  jim
02399 Fixed bug in save_simple which mean that ps.pafphu was allocated for every
02400 file rather than just the first one
02401 
02402 Revision 1.19  2007/04/04 16:05:59  jim
02403 Modified to make paf information a bit more correct
02404 
02405 Revision 1.18  2007/04/04 10:34:55  jim
02406 Modified to use new dfs tags
02407 
02408 Revision 1.17  2007/03/14 22:08:54  jim
02409 Fixed typo
02410 
02411 Revision 1.16  2007/03/14 14:49:13  jim
02412 Fixed problem with missing paf files in jmp recipes if detlive = F. Also
02413 fixed problem where extra dummy products were being created
02414 
02415 Revision 1.15  2007/03/13 09:50:53  jim
02416 Fixed bug in vircam_jmp_save_illum where PAF wasn't being saved for the
02417 first extension
02418 
02419 Revision 1.14  2007/03/06 14:06:24  jim
02420 Fixed missing paf write
02421 
02422 Revision 1.13  2007/03/06 13:49:48  jim
02423 Created static routine vircam_jmp_outfile to create a 'predictable'
02424 output filename for each product
02425 
02426 Revision 1.12  2007/03/01 12:42:42  jim
02427 Modified slightly after code checking
02428 
02429 Revision 1.11  2007/02/25 06:34:20  jim
02430 Plugged memory leak
02431 
02432 Revision 1.10  2007/02/15 06:59:38  jim
02433 Added ability to write QC paf files
02434 
02435 Revision 1.9  2007/02/06 13:11:12  jim
02436 Fixed entry for PRO dictionary in cpl_dfs_set_product_header
02437 
02438 Revision 1.8  2007/01/17 23:54:00  jim
02439 Plugged some memory leaks
02440 
02441 Revision 1.7  2006/12/19 13:30:01  jim
02442 Fixed vircam_jmp_illum to initialise pointers
02443 
02444 Revision 1.6  2006/12/18 12:51:20  jim
02445 Tightened up some of the error reporting
02446 
02447 Revision 1.5  2006/12/15 09:58:28  jim
02448 read noise and gain keywords were wrong...
02449 
02450 Revision 1.4  2006/11/29 12:28:45  jim
02451 Modified so that the correct recipe names would appear in the headers of
02452 data products
02453 
02454 Revision 1.3  2006/11/28 22:10:22  jim
02455 Added illcor_rms to calling sequence of vircam_illum
02456 
02457 Revision 1.2  2006/11/28 20:56:31  jim
02458 Added vircam_jmp_illum and vircam_jmp_save_illum
02459 
02460 Revision 1.1  2006/11/27 11:54:37  jim
02461 Initial entry
02462 
02463 
02464 */
02465 

Generated on Wed Apr 10 04:01:56 2013 for VIRCAM Pipeline by  doxygen 1.5.1