""" A bunch of useful python scripts that: * return equatorial/Galactic coordinates of vertices of OGLE-IV fields * plot OGLE-IV fields/subfields in equatorial/Galactic coordinates P. Mroz, 26 May 2019 """ import numpy as np import matplotlib.pyplot as plt import matplotlib.colors as colors from matplotlib.patches import Polygon from matplotlib.collections import PatchCollection PIXEL_SIZE = 0.2585 # 1 pix = 0.2585 arcsec def o4_field (ra0, dec0): """ Returns the equatorial coordinates of vertices of an OGLE-IV field centered on (ra0, dec0). """ ra, dec = [], [] pix = PIXEL_SIZE / 3600.0 pixra = PIXEL_SIZE / 3600.0 / np.cos(dec0*np.pi/180.0) ra.append(ra0-7518.0*pixra) dec.append(dec0+4321.0*pix) ra.append(ra0-7518.0*pixra) dec.append(dec0+8610.0*pix) ra.append(ra0+7518.0*pixra) dec.append(dec0+8610.0*pix) ra.append(ra0+7518.0*pixra) dec.append(dec0+4321.0*pix) ra.append(ra0+9616.0*pixra) dec.append(dec0+4321.0*pix) ra.append(ra0+9616.0*pixra) dec.append(dec0-4321.0*pix) ra.append(ra0+7518.0*pixra) dec.append(dec0-4321.0*pix) ra.append(ra0+7518.0*pixra) dec.append(dec0-8610.0*pix) ra.append(ra0-7518.0*pixra) dec.append(dec0-8610.0*pix) ra.append(ra0-7518.0*pixra) dec.append(dec0-4321.0*pix) ra.append(ra0-9616.0*pixra) dec.append(dec0-4321.0*pix) ra.append(ra0-9616.0*pixra) dec.append(dec0+4321.0*pix) return np.array(ra), np.array(dec) def o4_subfield (ra0, dec0, chip): """ Returns the equatorial coordinates of vertices of an OGLE-IV subfield centered on (ra0, dec0). """ ra, dec = [], [] pix = PIXEL_SIZE / 3600.0 pixra = PIXEL_SIZE / 3600.0 / np.cos(dec0*np.pi/180.0) assert type(chip) is int, 'chip is not an integer: %r' % chip assert chip >=1 and chip <= 32, 'chip >=1 and chip <= 32' if chip >= 1 and chip <= 7: ra.append(ra0+7468.0*pixra-(chip-1.0)*2148.0*pixra) dec.append(dec0-8610.0*pix) ra.append(ra0+7468.0*pixra-(chip-1.0)*2148.0*pixra) dec.append(dec0-4508.0*pix) ra.append(ra0+5420.0*pixra-(chip-1.0)*2148.0*pixra) dec.append(dec0-4508.0*pix) ra.append(ra0+5420.0*pixra-(chip-1.0)*2148.0*pixra) dec.append(dec0-8610.0*pix) elif chip >= 8 and chip <= 16: ra.append(ra0+9616.0*pixra-(chip-8.0)*2148.0*pixra) dec.append(dec0-4134.0*pix) ra.append(ra0+9616.0*pixra-(chip-8.0)*2148.0*pixra) dec.append(dec0-32.0*pix) ra.append(ra0+7568.0*pixra-(chip-8.0)*2148.0*pixra) dec.append(dec0-32.0*pix) ra.append(ra0+7568.0*pixra-(chip-8.0)*2148.0*pixra) dec.append(dec0-4134.0*pix) elif chip >= 17 and chip <= 25: ra.append(ra0+9616.0*pixra-(chip-17.0)*2148.0*pixra) dec.append(dec0+4134.0*pix) ra.append(ra0+9616.0*pixra-(chip-17.0)*2148.0*pixra) dec.append(dec0+32.0*pix) ra.append(ra0+7568.0*pixra-(chip-17.0)*2148.0*pixra) dec.append(dec0+32.0*pix) ra.append(ra0+7568.0*pixra-(chip-17.0)*2148.0*pixra) dec.append(dec0+4134.0*pix) elif chip >= 26 and chip <= 32: ra.append(ra0+7468.0*pixra-(chip-26.0)*2148.0*pixra) dec.append(dec0+8610.0*pix) ra.append(ra0+7468.0*pixra-(chip-26.0)*2148.0*pixra) dec.append(dec0+4508.0*pix) ra.append(ra0+5420.0*pixra-(chip-26.0)*2148.0*pixra) dec.append(dec0+4508.0*pix) ra.append(ra0+5420.0*pixra-(chip-26.0)*2148.0*pixra) dec.append(dec0+8610.0*pix) return np.array(ra), np.array(dec) def equatorial_to_galactic (_ra, _dec): """ Transforming equatorial to Galactic coordinates""" deg = np.pi/180.0 ra_g = 192.85948*deg dec_g = 27.12825*deg l_NGP = 122.93192*deg ra = _ra*deg dec = _dec*deg sinb = np.cos(dec)*np.cos(dec_g)*np.cos(ra-ra_g)+np.sin(dec)*np.sin(dec_g) sinlcosb = np.cos(dec)*np.sin(ra-ra_g) coslcosb = np.sin(dec)*np.cos(dec_g)-np.cos(dec)*np.sin(dec_g)*np.cos(ra-ra_g) b = np.arcsin(sinb) l = np.arctan2(sinlcosb,coslcosb) l = l_NGP - l b /= deg l /= deg l[l>180] -= 360.0 return l, b #-----Reading coordinates of OGLE-IV fields from Table 6-----# field, _ra, _dec, nstars = np.loadtxt('table6.dat',unpack=True,dtype='S6,float,float,float',usecols=(0,1,2,5)) ra = dict(zip(field,_ra)) dec = dict(zip(field,_dec)) """ # EXAMPLE 1: Getting equatorial coordinates of vertices of field BLG660 ra_v, dec_v = o4_field(ra['BLG660'],dec['BLG660']) for (x,y) in zip(ra_v,dec_v): print '%8.4f %8.4f' % (x, y) plt.plot(ra_v,dec_v) plt.gca().invert_xaxis() plt.show() """ # EXAMPLE 2: Plotting OGLE-IV fields in the Galactic coordinates fig, ax = plt.subplots() patches = [] for f in field: x, y = o4_field(ra[f],dec[f]) xg, yg = equatorial_to_galactic(x,y) p = Polygon(np.column_stack((xg,yg)),closed=True) patches.append(p) ax.text(np.mean(xg),np.mean(yg),f.replace('BLG',''),va='center',ha='center') norm = colors.LogNorm(vmin=0.5, vmax=5.0) pc = PatchCollection(patches, cmap='viridis', norm=norm) pc.set_array(nstars) pc.set_edgecolor('none') ax.add_collection(pc) cbar = plt.colorbar(pc) cbar.set_label('Number of stars in OGLE-IV databases (million)') ax.set_xlim((13,-13)) ax.set_ylim((-18,8)) ax.set_xlabel('Galactic longitude') ax.set_ylabel('Galactic latitude') plt.show() """ # EXAMPLE 3: Plotting subfields of BLG500 in the equatorial coordinates fig, ax = plt.subplots() for chip in xrange(1,33): x, y = o4_subfield(ra['BLG500'],dec['BLG500'],chip) p = Polygon(np.column_stack((x,y)),closed=True) ax.add_patch(p) ax.text(np.mean(x),np.mean(y),str(chip),va='center',ha='center') ax.set_xlim((270,266)) ax.set_ylim((-30,-27)) ax.set_xlabel('Right Ascension') ax.set_ylabel('Declination') plt.show() """ """ # EXAMPLE 4: Plotting subfields of BLG500 in the Galactic coordinates fig, ax = plt.subplots() for chip in xrange(1,33): x, y = o4_subfield(ra['BLG500'],dec['BLG500'],chip) xg, yg = equatorial_to_galactic(x,y) p = Polygon(np.column_stack((xg,yg)),closed=True) ax.add_patch(p) ax.text(np.mean(xg),np.mean(yg),str(chip),va='center',ha='center') ax.set_xlim((2,0)) ax.set_ylim((-2,0)) ax.set_xlabel('Galactic longitude') ax.set_ylabel('Galactic latitude') plt.show() """